eic-opticks/ana/GScintillatorLib.py

0001 #!/usr/bin/env python
0002 #
0003 # Copyright (c) 2019 Opticks Team. All Rights Reserved.
0004 #
0005 # This file is part of Opticks
0006 # (see https://bitbucket.org/simoncblyth/opticks).
0007 #
0008 # Licensed under the Apache License, Version 2.0 (the "License");
0009 # you may not use this file except in compliance with the License.
0010 # You may obtain a copy of the License at
0011 #
0012 #   http://www.apache.org/licenses/LICENSE-2.0
0013 #
0014 # Unless required by applicable law or agreed to in writing, software
0015 # distributed under the License is distributed on an "AS IS" BASIS,
0016 # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
0017 # See the License for the specific language governing permissions and
0018 # limitations under the License.
0019 #
0020 """
0021 ::
0022
0023    an ; ipython -i $(which GScintillatorLib.py)
0024    mkdir -p ~/simoncblyth.bitbucket.io/env/presentation/ana/GScintillatorLib
0025    cp /tmp/ana/GScintillatorLib/*.png ~/simoncblyth.bitbucket.io/env/presentation/ana/GScintillatorLib/
0026
0027 """
0028 import os, numpy as np
0029 from opticks.ana.main import opticks_main
0030 from opticks.ana.nload import np_load
0031 from opticks.ana.key import keydir
0032 import matplotlib.pyplot as plt
0033
0034 if __name__ == '__main__':
0035
0036     ok = opticks_main()
0037     kd = keydir(os.environ["OPTICKS_KEY"])
0038     aa,aa_paths = np_load(os.path.join(kd,"GScintillatorLib/GScintillatorLib.npy"))
0039     fc,fc_paths = np_load(os.path.join(kd,"GScintillatorLib/LS/FASTCOMPONENT.npy"))
0040     sc,sc_paths = np_load(os.path.join(kd,"GScintillatorLib/LS/SLOWCOMPONENT.npy"))
0041     a0 = aa[0,:,0]
0042     a1 = aa[1,:,0]
0043     a2 = aa[2,:,0]
0044     b = np.linspace(0,1,len(a0))
0045
0046     print("aa:%s" % str(aa.shape))
0047     print("a0:%s" % str(a0.shape))
0048     print("a1:%s" % str(a1.shape))
0049     print("a2:%s" % str(a2.shape))
0050     print(" b:%s" % str(b.shape))
0051     print("fc:%s" % str(fc.shape))
0052     print("sc:%s" % str(sc.shape))
0053
0054     assert aa.shape == (1, 4096, 1) or  aa.shape == (3, 4096, 1)
0055     assert np.all( fc == sc )
0056     assert a0.shape == (4096,)
0057     assert a1.shape == (4096,)
0058     assert a2.shape == (4096,)
0059
0060     fig = plt.figure(figsize=ok.figsize)
0061     plt.title("Inverted Cumulative Distribution Function for Scintillator Wavelength Generation " )
0062     ax = fig.add_subplot(1,1,1)
0063     ax.plot( b, a0, label="a0" )
0064     #ax.plot( b, a1, label="a1" )
0065     #ax.plot( b, a2, label="a2" )
0066
0067     ax.set_ylabel("Wavelength (nm)")
0068     ax.set_xlabel("Probability")
0069
0070     xlim = ax.get_xlim()
0071     ylim = ax.get_ylim()
0072
0073     plhs = 0.05
0074     prhs = 0.95
0075
0076     wlhs = np.interp( plhs, b, a0 )
0077     wrhs = np.interp( prhs, b, a0 )
0078
0079     ax.plot( [plhs,plhs], ylim, linestyle="dashed" )
0080     ax.plot( xlim       , [wlhs,wlhs], linestyle="dashed" )
0081
0082
0083     ax.plot( [prhs,prhs], ylim, linestyle="dashed" )
0084     ax.plot( xlim       , [wrhs,wrhs], linestyle="dashed" )
0085
0086     fig.show()
0087
0088     path="/tmp/ana/GScintillatorLib/icdf.png"
0089     fold=os.path.dirname(path)
0090     if not os.path.isdir(fold):
0091        os.makedirs(fold)
0092     pass
0093     log.info("save to %s " % path)
0094     fig.savefig(path)
0095
0096
0097