eic-opticks/ana/gsplt.py

0001 #!/usr/bin/env python
0002 """
0003
0004 TODO:
0005
0006 * genstep deltaPosition arrows
0007
0008 """
0009
0010 from opticks.ana.gs import GS
0011 from opticks.analytic.gdml import GDML
0012 import os, numpy as np
0013 import matplotlib.pyplot as plt
0014 import matplotlib.lines as mlines
0015 from collections import OrderedDict as odict
0016
0017 X,Y,Z,T = 0,1,2,3
0018
0019
0020 if __name__ == '__main__':
0021
0022     g = GDML.parse("$OPTICKS_PREFIX/tds_ngt_pcnk.gdml")
0023
0024     lvs = odict()
0025     lvs[0] = g.find_one_volume("lPoolLining")
0026     lvs[1] = g.find_one_volume("lInnerWater")
0027     lvs[2] = g.find_one_volume("lAcrylic")
0028     lvs[3] = g.find_one_volume("lOuterWaterPool")
0029
0030     plt.close()
0031     plt.ion()
0032     fig = plt.figure(figsize=(6,6))
0033     fig.suptitle("suptitle")
0034
0035     plt.title("plt.title")
0036
0037     ax = fig.add_subplot(111)
0038     ax.set_xlim(-26000,26000)
0039     ax.set_ylim(-26000,26000)
0040
0041     for lv in lvs.values():
0042         s = lv.solid
0043         kwa = {}
0044         sh = s.as_shape(**kwa)
0045         for pt in sh.patches():
0046             ax.add_patch(pt)
0047         pass
0048     pass
0049
0050     ysli = 1000   # slice around y=0
0051     hpos = lvs[1].physvol_xyz('pLPMT_Hamamatsu','position')
0052     npos = lvs[1].physvol_xyz('pLPMT_NNVT','position')
0053     wpos = lvs[3].physvol_xyz("mask_PMT_20inch_vetolMaskVirtual_phys", "position")
0054
0055     hpos_ = hpos[np.abs(hpos[:,1]) < ysli]
0056     npos_ = npos[np.abs(npos[:,1]) < ysli]
0057     wpos_ = wpos[np.abs(wpos[:,1]) < ysli]
0058
0059     ax.scatter( hpos_[:,0], hpos_[:,2], 5.0 )
0060     ax.scatter( npos_[:,0], npos_[:,2], 5.0 )
0061     ax.scatter( wpos_[:,0], wpos_[:,2], 5.0 )
0062
0063
0064     path = "$TMP/evt/g4live/natural/1/gs.npy"
0065     args = GS.parse_args(__doc__, path=path)
0066     gs = GS(args.path)
0067
0068     xyzt = gs.xyzt
0069     ijku = gs.deltaPositionLength
0070     pdgc = gs.pdgCode
0071
0072     sli = np.abs(xyzt[:,Y]) < ysli
0073     nsli = np.count_nonzero(sli)
0074
0075     xyzt_ = xyzt[sli]
0076     ijku_ = ijku[sli]
0077
0078     #ax.scatter(xyzt_[:,0],xyzt_[:,2], 0.1)
0079
0080     kwa = dict(linewidth=0.1)
0081
0082     for i in range(nsli):
0083         l = mlines.Line2D(
0084                    [xyzt_[i,X], xyzt_[i,X] + ijku_[i,X]],
0085                    [xyzt_[i,Z], xyzt_[i,Z] + ijku_[i,Z]], **kwa )
0086         ax.add_line(l)
0087     pass
0088
0089
0090
0091
0092