qudarap/tests/QScintTest0.py

0001 #!/usr/bin/env python
0002 """
0003
0004
0005 * histo plotting of domain compressed wavelength ?
0006
0007 """
0008
0009 import os, sys, logging, numpy as np
0010 import matplotlib.pyplot as plt
0011 from opticks.ana.nload import np_load
0012 from opticks.ana.key import keydir
0013
0014 log = logging.getLogger(__name__)
0015
0016 if __name__ == '__main__':
0017
0018     ok = os.environ["OPTICKS_KEY"]
0019     kd = keydir(ok)
0020     aa = np.load(os.path.join(kd,"GScintillatorLib/GScintillatorLib.npy"))
0021     fc = np.load(os.path.join(kd,"GScintillatorLib/LS/FASTCOMPONENT.npy"))
0022     sc = np.load(os.path.join(kd,"GScintillatorLib/LS/SLOWCOMPONENT.npy"))
0023
0024     print("aa:%s" % str(aa.shape))
0025     print("fc:%s" % str(fc.shape))
0026     print("sc:%s" % str(sc.shape))
0027     assert aa.shape == (1, 4096, 1)
0028     assert np.all( fc == sc )
0029
0030     fold = "/tmp/QScintTest"
0031     w = np.load(os.path.join(fold, "wavelength.npy"))
0032
0033
0034     p = np.load(os.path.join(fold, "photon.npy"))
0035     dir = p[:,1,:3]
0036     pol = p[:,2,:3]
0037     dirpol = (dir*pol).sum(axis=1)   # check transverse pol
0038     dirpol_max = np.abs(dirpol).max()
0039     assert dirpol_max < 1e-6, dirpol_max
0040
0041
0042     plt.ion()
0043
0044     wd = np.linspace(60,820,256) - 1.
0045     # reduce bin edges by 1nm to avoid aliasing artifact in the histogram
0046
0047     mid = (wd[:-1]+wd[1:])/2.     # bin middle
0048
0049     counts, edges = np.histogram(w, bins=wd )
0050     fcounts = counts.astype(np.float32)
0051     fcounts  /= fcounts.sum()
0052     fcounts  /= 2.98   # bin width nm
0053
0054     plt.close()
0055
0056     plt.plot( edges[:-1], fcounts, drawstyle="steps-mid")
0057
0058
0059     _fc = fc[:,1]
0060     _fc /= _fc.sum()
0061     _fc /= 20.       # bin width nm
0062
0063     plt.plot( fc[:,0], _fc )
0064
0065     #plt.plot( mid,  pl )
0066
0067     plt.axis( [w.min() - 100, w.max() + 100, 0, fcounts.max()*1.1 ])
0068
0069
0070     #plt.hist(w, bins=256)   # 256 is number of unique wavelengths (from record compression)
0071
0072