qudarap/tests/fake_propagate.py

0001 #!/usr/bin/env python
0002
0003 import os, numpy as np, textwrap
0004 import logging
0005 log = logging.getLogger(__name__)
0006 logging.basicConfig(level=logging.INFO)
0007 # non-standard logging config location to catch logging from CSGFoundry
0008 # OK because this is a "main" script
0009
0010 from opticks.ana.fold import Fold
0011 from opticks.ana.pvplt import *
0012 from opticks.ana.p import *
0013 from opticks.ana.r import *
0014
0015 from opticks.sysrap.sevt import SEvt
0016
0017
0018
0019 if __name__ == '__main__':
0020     fold = os.path.expandvars("$TMP/GEOM/$GEOM/QSimTest/ALL0_fake_propagate/A000")
0021     t = SEvt.Load(fold, symbol="t", NEVT=0)
0022     print(repr(t))
0023
0024     PIDX = int(os.environ.get("PIDX","-1"))
0025     FONT_SIZE= int(os.environ.get("FONT_SIZE","36"))
0026     FACTOR = float(os.environ.get("FACTOR","60"))
0027
0028
0029
0030 if 0:
0031     p = t.f.photon
0032     r = t.f.record
0033     seq = t.f.seq
0034     prd = t.f.prd
0035     h = t.f.hit if hasattr(t.f,"hit") else None
0036     h_meta = t.f.hit_meta
0037
0038
0039     ## HMM:makes more sense to put this meta data on the domain, for when no hits
0040     #pyhit = hit__(p, hitmask)  # hits selected in python
0041     #assert np.all( pyhit == h )
0042
0043     if not r is None:
0044         r_post = r[:,:,0]
0045         r_pos = r[:,:,0,:3].reshape(-1,3)
0046         r_mom = r[:,:,1,:3].reshape(-1,3)
0047         r_pol = r[:,:,2,:3].reshape(-1,3)
0048         r_flag = flag__(r.reshape(-1,4,4))
0049
0050         r_flag_label = hm.label( r_flag )
0051
0052         r_cells = make_record_cells( r )
0053         # indices into flattened records
0054
0055
0056         ## try heuristic to arrange viewpoint and scale
0057         ## HMM: how to automate picking the view ?
0058         ## tuple: camera location, focus point, viewup vector
0059
0060         r_pos_max = np.max( r_pos, axis=0 )
0061         r_pos_min = np.min( r_pos, axis=0 )
0062         r_pos_rng = r_pos_max - r_pos_min
0063         r_pos_avg = np.average(r_pos, axis=0 )
0064         r_pos_vec = np.array( [-1,0,0], dtype=np.float32  )
0065         r_pos_scale = r_pos_rng.max()
0066         camloc = r_pos_avg + r_pos_scale*r_pos_vec
0067         camera_position = [camloc, r_pos_avg, (0,0,1)]
0068         os.environ["ZOOM"] = "%g" % (1./r_pos_scale )
0069
0070
0071         PLOT = "PLOT" in os.environ
0072         if PLOT:
0073             r_poly = pv.PolyData()
0074             r_poly.points = r_pos
0075             r_poly.lines = r_cells.ravel()  ## line connectivity array
0076
0077             r_poly["flag_label"] = r_flag_label
0078             r_tube = r_poly.tube(radius=1)
0079
0080             pl = pvplt_plotter()
0081             pl.add_mesh( r_tube )
0082             pvplt_polarized( pl, r_pos, r_mom, r_pol, factor=FACTOR )
0083
0084             pl.add_point_labels(r_poly, "flag_label", point_size=20, font_size=FONT_SIZE )
0085
0086             pl.camera_position = camera_position
0087
0088
0089             pl.show()
0090         pass
0091     pass
0092
0093
0094     s = str(p[:,:3])
0095     a = np.array( s.split("\n") + [""] ).reshape(-1,4)
0096
0097     for i in range(len(a)):
0098         if not (PIDX == -1 or PIDX == i): continue
0099         if PIDX > -1: print("PIDX %d " % PIDX)
0100
0101         if not r is None:
0102             print("r[%(i)d,:,:3]" % locals())
0103             print(r[i,:,:3])
0104             print("\n\nbflagdesc_(r[i,j])")
0105             for j in range(len(r[i])):
0106                 expr = "bflagdesc_(r[%(i)d,%(j)d])" % locals()
0107                 print(expr)
0108                 print(eval(expr))
0109             pass
0110         pass
0111
0112         print("\n")
0113         print("p")
0114         print("\n".join(a[i]))
0115
0116         expr = "bflagdesc_(p[%(i)d]) " % locals()
0117         print(expr)
0118         print(eval(expr))
0119         print("\n")
0120
0121         if not seq is None:
0122             print(seqhis_(seq[i,0]))
0123             print("\n")
0124         pass
0125         print("\n\n")
0126     pass
0127
0128
0129     EXPR = list(filter(None,textwrap.dedent(r"""
0130     np.c_[cf.sim.stree.standard.bnd_names]
0131     np.all( cf.sim.stree.standard.bnd_names == cf.sim.extra.GGeo.bnd_names )
0132     np.c_[cf.sim.stree.suname[np.char.startswith(cf.sim.stree.suname, "Hama")]]
0133     np.c_[cf.sim.stree.suname[np.char.startswith(cf.sim.stree.suname, "NNVT")]]
0134     """).split("\n")))
0135
0136     for expr in EXPR:
0137         print(expr)
0138         print(eval(expr))
0139     pass
0140
0141
0142
0143