sysrap/tests/sreport_ab.py

0001 #!/usr/bin/env python
0002
0003 import os, numpy as np
0004 from opticks.ana.fold import Fold
0005 from opticks.ana.npmeta import NPMeta
0006 from opticks.sysrap.tests.sreport import Substamp
0007 from opticks.sysrap.tests.sreport import RUN_META
0008
0009 MODE = 2
0010 PLOT =  os.environ.get("PLOT","?plot")
0011 COMMANDLINE = os.environ.get("COMMANDLINE", "")
0012
0013
0014 if MODE != 0:
0015     from opticks.ana.pvplt import *
0016 pass
0017
0018
0019 class AB_Substamp_ALL_Etime_vs_Photon(object):
0020     def __init__(self, a, b):
0021         title = RUN_META.AB_Title(a,b)
0022
0023         af = a.substamp.a
0024         at = Substamp.ETime(af)
0025         an = a.submeta_NumPhotonCollected.a
0026         ah = Substamp.Subcount(af, "hit")
0027
0028         self.at = at
0029         self.an = an
0030         self.ah = ah
0031
0032
0033         bf = b.substamp.a
0034         bt = Substamp.ETime(bf)
0035         bn = b.submeta_NumPhotonCollected.a
0036         bh = Substamp.Subcount(bf, "hit")
0037
0038         self.bt = bt
0039         self.bn = bn
0040         self.bh = bh
0041
0042         assert np.all( an == bn )
0043         photon = an[:,0]/1e6
0044
0045
0046         fontsize = 20
0047         YSCALE_ = ["log", "linear" ]
0048         YSCALE = os.environ.get("YSCALE", YSCALE_[1])
0049         assert YSCALE in YSCALE_
0050
0051         SLICE = list(map(int, os.environ.get("SLICE", "3:12").split(":")))
0052
0053         YMIN = float(os.environ.get("YMIN", "1e-2"))
0054         YMAX = float(os.environ.get("YMAX", "45"))
0055
0056         if MODE == 2:
0057             fig, axs = mpplt_plotter(nrows=1, ncols=1, label=title, equal=False)
0058             ax = axs[0]
0059
0060             sli = slice(*SLICE)
0061             #sli = slice(3,12)
0062             #sli = slice(None)
0063
0064             deg = 1  # 1:lin 2:par
0065             a_fit = np.poly1d(np.polyfit(photon[sli], at[sli], deg))
0066             a_fit_label = "%s : linefit( slope %10.3f  intercept %10.3f )" % (a.symbol, a_fit.coef[0], a_fit.coef[1])
0067             self.a_fit = a_fit
0068
0069             b_fit = np.poly1d(np.polyfit(photon[sli], bt[sli], deg))
0070             b_fit_label = "%s : linefit( slope %10.3f  intercept %10.3f )" % (b.symbol, b_fit.coef[0], b_fit.coef[1])
0071             self.b_fit = b_fit
0072
0073             a_label = RUN_META.GPUMeta(a)
0074             b_label = RUN_META.GPUMeta(b)
0075
0076             ax.scatter( photon, at, label=a_label )
0077             ax.plot( photon[sli], a_fit(photon[sli]), label=a_fit_label )
0078
0079             ax.scatter( photon, bt, label=b_label )
0080             ax.plot( photon[sli], b_fit(photon[sli]), label=b_fit_label )
0081
0082             pass
0083             ax.set_xlim( -5, 105 );
0084             ax.set_ylim( YMIN, YMAX );  # 50*200 = 1e4
0085             ax.set_yscale(YSCALE)
0086             ax.set_ylabel("Event time (seconds)", fontsize=fontsize )
0087             ax.set_xlabel("Number of Photons (Millions)", fontsize=fontsize )
0088             ax.legend()
0089             ax.legend()
0090             fig.show()
0091         pass
0092
0093
0094
0095
0096
0097 if __name__ == '__main__':
0098
0099     print("[sreport_ab.py:PLOT[%s]" % PLOT )
0100
0101     asym = os.path.expandvars("$A")
0102     print("[sreport_ab.py:fold = Fold.Load A_SREPORT_FOLD [%s]" % asym )
0103     a = Fold.Load("$A_SREPORT_FOLD", symbol=asym )
0104     print("]sreport_ab.py:a = Fold.Load" )
0105
0106     print("[sreport_ab.py:repr(a)" )
0107     print(repr(a))
0108     print("]sreport_ab.py:repr(a)" )
0109
0110
0111     bsym = os.path.expandvars("$B")
0112     print("[sreport_ab.py:fold = Fold.Load B_SREPORT_FOLD [%s]" % bsym )
0113     b = Fold.Load("$B_SREPORT_FOLD", symbol=bsym )
0114     print("]sreport_ab.py:b = Fold.Load" )
0115
0116     print("[sreport_ab.py:repr(b)" )
0117     print(repr(b))
0118     print("]sreport_ab.py:repr(b)" )
0119
0120     print("]sreport_ab.py:PLOT[%s]" % PLOT )
0121
0122
0123     if PLOT.startswith("AB_Substamp_ALL_Etime_vs_Photon"):
0124         ab = AB_Substamp_ALL_Etime_vs_Photon(a,b)
0125     pass
0126
0127
0128
0129
0130
0131
0132
0133