CSG/tests/CSGSimtraceTest.py

0001 #!/usr/bin/env python
0002 """
0003 CSG/tests/CSGSimtraceTest.py
0004 ===============================
0005
0006 HMM: almost identical to extg4/tests/X4SimtraceTest.py
0007
0008 TODO: more duplication avoidance but keep distinct mains
0009
0010 """
0011 import os, logging, builtins, numpy as np
0012 log = logging.getLogger(__name__)
0013 from opticks.ana.fold import Fold
0014
0015 import matplotlib.pyplot as mp
0016 from opticks.ana.fold import Fold
0017 from opticks.sysrap.sframe import sframe , X, Y, Z
0018 from opticks.CSG.Values import Values
0019 from opticks.ana.pvplt import mpplt_simtrace_selection_line, mpplt_hist, mpplt_parallel_lines_auto, mpplt_add_shapes
0020 from opticks.ana.eget import efloatarray_
0021
0022
0023 if __name__ == '__main__':
0024
0025     #fmt = '[%(asctime)s] p%(process)s {%(pathname)s:%(lineno)d} %(levelname)s - %(message)s'
0026     fmt = '{%(pathname)s:%(lineno)d} %(levelname)s - %(message)s'
0027     logging.basicConfig(level=logging.INFO, format=fmt)
0028
0029     SYMBOLS = os.environ.get("SYMBOLS", None)
0030     FOLD = os.environ.get("FOLD", None)
0031
0032     if SYMBOLS is None and not FOLD is None:
0033         s = Fold.Load(symbol="s")
0034         t = None
0035         fr = s.sframe
0036         s_label = os.environ["GEOM"]
0037         SYMBOLS = "S"
0038         sv = Values.Find("$FOLD", symbol="sv") if not s_label is None else None
0039     elif not SYMBOLS is None:
0040         log.info("SYMBOLS defined %s proceed with Fold.MultiLoad", SYMBOLS)
0041         ff = Fold.MultiLoad()
0042         log.info(" ff %s " % str(ff))
0043         assert len(ff) > 0
0044         frs = list(filter(None, map(lambda f:f.sframe, ff)))
0045         assert len(frs) > 0
0046         fr = frs[0]       ## HMM: picking first frame, maybe need to form composite bbox from all frames ?
0047         sv = None
0048     else:
0049         assert 0
0050     pass
0051
0052     fig, ax = fr.mp_subplots(mp)
0053
0054     if not s is None:
0055         s_hit = s.simtrace[:,0,3]>0
0056         s_pos = s.simtrace[s_hit][:,1,:3]
0057
0058         ## investigate unexpected top cap intersects : FIXED z2cap TYPO
0059         ##
0060         ##                                                                    intersect pos x > 120.  t > 0
0061         ## e = np.logical_and( s.simtrace[:,2,0] > 100., np.logical_and( s.simtrace[:,1,0] > 120. , s.simtrace[:,0,3]>0 ))
0062         ##
0063         ## e_ori = s.simtrace[e][:100,2,:3]
0064         ## e_dir = s.simtrace[e][:100,3,:3]
0065         ## fr.mp_scatter(e_ori, label="e_ori", s=2 )
0066         ## fr.mp_arrow(  e_ori, 10*e_dir, label="e_ori,e_dir", s=2 )
0067     pass
0068
0069     if s_label.startswith("nmskSolidMaskVirtual") and not sv is None:
0070         r1=sv.get("SolidMask.SolidMaskVirtual.rOuter2.mask_radiu_virtual")
0071         r2=sv.get("SolidMask.SolidMaskVirtual.rOuter3.mask_radiu_virtual/2")
0072         z1=sv.get("SolidMask.SolidMaskVirtual.zPlane2.htop_out/2")
0073         z2=sv.get("SolidMask.SolidMaskVirtual.zPlane3.htop_out+MAGIC_virtual_thickness")
0074         z0=(z2*r1-z1*r2)/(r1-r2)
0075         ax.set_ylim( -240, z0+50 )
0076     pass
0077
0078     mpplt_add_shapes(ax)  # ELLIPSE0 ELLIPSE1 RECTANGLE0 RECTANGLE1
0079
0080
0081     if not s is None and "UNEXPECTED" in os.environ:
0082
0083         if s_label.startswith("nmskSolidMaskVirtual"):
0084             w_label, w = "apex glancers",  np.logical_and( np.abs(s.simtrace[:,1,0]) < 220, np.abs(s.simtrace[:,1,2]-98) < 1 )
0085             #w_label, w = "quadratic precision loss", np.logical_and( np.abs(s.simtrace[:,1,0] - (-214)) < 5, np.abs(s.simtrace[:,1,2] - (115)) < 5 )
0086
0087             w_simtrace = s.simtrace[w]
0088             w_path = "/tmp/simtrace_sample.npy"
0089             np.save(w_path, w_simtrace)
0090         else:
0091             w = np.logical_and( np.abs(s.simtrace[:,1,2]) > 0.20 , s.simtrace[:,0,3]>0 )
0092             w_simtrace = s.simtrace[w][::10]
0093         pass
0094         log.info("UNEXPECTED w_simtrace : %s " % str(w_simtrace.shape))
0095         w_ori = w_simtrace[:,2,:3]
0096         fr.mp_scatter(w_ori, label="w_ori", s=1 )
0097         mpplt_simtrace_selection_line(ax, w_simtrace, axes=fr.axes, linewidths=2)
0098     pass
0099
0100     if not s is None:
0101         mpplt_parallel_lines_auto( ax, fr.bbox.T, fr.axes, linestyle="dashed" )
0102     pass
0103     if not s is None and hasattr(s,"genstep") and "GS" in os.environ:
0104         s_gs  = s.genstep[:,5,:3]
0105         fr.mp_scatter(s_gs, label="s_gs", s=1 )
0106     pass
0107
0108     if not s is None and hasattr(s,"simtrace_selection"):
0109        sts = s.simtrace_selection
0110     elif not s is None and "SELECTION" in os.environ:
0111         #w = np.logical_and( s.simtrace[:,0,3]>0, np.logical_and( s.simtrace[:,1,Z] > -38.9, s.simtrace[:,1,Z] < -20. ))
0112         w = s.simtrace[:,1,X] > 264.5
0113         wi = np.where(w)[0]
0114         sts = s.simtrace[w][:50]
0115     else:
0116         sts = None
0117     pass
0118
0119     if not sts is None:
0120         mpplt_simtrace_selection_line(ax, sts, axes=fr.axes, linewidths=2)
0121     pass
0122
0123     d = getattr(s, "CSGDebug_Cylinder", None) if not s is None else None
0124
0125     log.info("SYMBOLS %s " % str(SYMBOLS))
0126     if not SYMBOLS is None:
0127         for A in list(SYMBOLS):
0128             a = A.lower()
0129             if hasattr(builtins, a):
0130                 fold = getattr(builtins, a)
0131                 label = getattr(builtins, "%s_label" % a )
0132
0133             elif a in globals():
0134                 fold = globals()[a]
0135                 label = globals()["%s_label" % a]
0136             else:
0137                 log.fatal(" a %s not in builtins or globals " % a )
0138                 assert 0
0139             pass
0140             a_offset = efloatarray_("%s_OFFSET" % A, "0,0,0")
0141             log.info("A %s a %s label %s" % ( A, a, label ) )
0142             a_hit = fold.simtrace[:,0,3]>0
0143             a_pos = a_offset + fold.simtrace[a_hit][:,1,:3]
0144             label = label.split("__")[0] if "__" in label else label
0145             fr.mp_scatter(a_pos, label="%s:%s" % (A,label), s=1 )
0146             pass
0147         pass
0148     pass
0149     fr.mp_legend()
0150     fig.show()
0151 pass
0152
0153