eic-opticks/ana/ucf.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 ucf.py : Comparing Opticks and Geant4 random consumption
0023 =============================================================
0024
0025 Compares the randoms consumed by Geant4 and Opticks and the
0026 positions where the consumption happens.  Relies on mask
0027 running of single photons
0028
0029 ::
0030
0031     ucf.py 1230
0032     UCF_PINDEX_DEV=-1 ucf.py 1230
0033     UCF_PINDEX_DEV=1230 ucf.py 1230
0034
0035
0036 Parse the kernel print log in a u_rng centric fashion::
0037
0038      tboolean-;tboolean-box --okg4 --align --mask 1230 --pindex 0 --pindexlog -DD
0039
0040           ## write kernel pindexlog for photon 1230
0041
0042     ucf.py 1230
0043
0044           ## parse the log and compare with expected rng sequence
0045
0046
0047 NB this is invoked from CRandomEngine::preTrack when using --mask option,
0048 as printindexlog stdout redirection borks stdout for subprocesses
0049 this script must not write to stdout.
0050
0051 cfg4/CRandomEngine.cc::
0052
0053     367         const char* cmd = BStr::concat<unsigned>("ucf.py ", mask_index, NULL );
0054     368         LOG(info) << "CRandomEngine::preTrack : START cmd \"" << cmd << "\"";
0055     369         int rc = SSys::run(cmd) ;  // NB must not write to stdout, stderr is ok though
0056     370         assert( rc == 0 );
0057     371         LOG(info) << "CRandomEngine::preTrack : DONE cmd \"" << cmd << "\"";
0058
0059
0060 """
0061 from __future__ import print_function
0062 import os, sys, re
0063
0064
0065 def print_(s):
0066     stream = sys.stderr
0067     print(s, file=stream)
0068
0069
0070 class U(object):
0071
0072     XRNG = []   # set by UCF.__init__
0073
0074     @classmethod
0075     def find(cls, u, tolerance=1e-6):
0076         """
0077         :return: index of u within XRNG sequence, or None if not found
0078         """
0079         idxf = None
0080         for i in range(len(cls.XRNG)):
0081             if abs(u-cls.XRNG[i]) < tolerance:
0082                 idxf = i
0083                 break
0084             pass
0085         pass
0086         return idxf
0087
0088     def __init__(self, idx, lab, val, lines):
0089         """
0090         :param idx: index into XRNG sequence
0091         :param lab:
0092         :param val:
0093         :param lines:
0094         """
0095         cls = self.__class__
0096         self.idx = idx
0097         self.lab = lab
0098         self.val = val
0099         self.fval = float(val)
0100
0101         assert idx < len(cls.XRNG)
0102         xval = cls.XRNG[idx]
0103
0104         idxf = cls.find(self.fval)
0105         idxd = idxf - idx
0106
0107         self.idxf = idxf
0108         self.idxd = idxd
0109
0110         self.xval = xval
0111         self.lines = lines[:]      ## must copy, not reference
0112         self.tail = []
0113
0114     def _get_hdr(self):
0115         fval = "%.9f" % self.fval
0116         xval = "%.9f" % self.xval
0117         mrk = "    " if fval == xval else "%+3d*" % self.idxd
0118         return " [%3d|%3d] %50s : %2s : %s : %s : %d " % ( self.idx, self.idxf, self.lab, mrk, fval, xval , len(self.lines) )
0119     hdr = property(_get_hdr)
0120
0121     def __str__(self):
0122         return "\n".join(self.lines + [self.hdr,""] + self.tail)
0123
0124     def __repr__(self):
0125         return self.hdr
0126
0127
0128 class UCF(list):
0129     @classmethod
0130     def rngpath(cls):
0131         return os.path.expandvars("$TMP/TRngBufTest_0.npy" )
0132     @classmethod
0133     def rngpathtxt(cls, pindex):
0134         return os.path.expandvars("$TMP/TRngBufTest_%s.txt" % pindex )
0135     @classmethod
0136     def printlogpath(cls, pindex):
0137         """
0138         :return: path to single photon log, obtained by redirection of OptiX output stream
0139         """
0140         return os.path.expandvars("$TMP/ox_%s.log" % pindex )
0141
0142     @classmethod
0143     def loadrngtxt(cls, pindex):
0144         """
0145         workaround lldb python failing to import numpy
0146         """
0147         trng_ = cls.rngpathtxt(pindex)
0148         trng = os.path.expandvars(trng_)
0149         assert os.path.exists(trng), (trng, trng_, "non-existing-trng")
0150         return map(float, file(trng).readlines())
0151
0152     def __init__(self, pindex):
0153         """
0154         :param pindex: photon record index
0155
0156         1. collects random consumption reported in OptiX single photon log into this list
0157         2. loads the expected sequence of randoms for this photon
0158
0159         """
0160         list.__init__(self)
0161
0162         evar = "UCF_PINDEX_DEV"
0163         upindex = int(os.environ.get(evar, pindex))
0164         if upindex != pindex:
0165             print_("WARNING evar active %s " % evar )
0166         pass
0167         path = self.printlogpath(upindex)
0168         xrng = self.loadrngtxt(pindex)
0169
0170         print_("path %s " % path )
0171
0172
0173         U.XRNG = xrng
0174
0175         self.pindex = pindex
0176         self.path = path
0177         self.parse(path)
0178
0179
0180     PTN = re.compile("u_(\S*):\s*(\S*)\s*")
0181
0182     def parse(self, path):
0183         """
0184         Parses the single photon log, collecting consumption lines
0185         into this list as U instances.
0186         """
0187         self.lines = map(lambda line:line.rstrip(),file(path).readlines())
0188         curr = []
0189         for i, line in enumerate(self.lines):
0190             m = self.PTN.search(line)
0191             #print "%2d : %s" % ( i, line)
0192             curr.append(line)
0193             if m is None or line[0] == "#": continue
0194
0195             sname = m.group(1)
0196             srng = m.group(2)
0197
0198             idx = len(self)
0199             u = U(idx, sname, srng, curr )
0200             curr[:] = []  ## clear collected lines
0201
0202             self.append(u)
0203         pass
0204         self[-1].tail = curr[:]
0205
0206
0207     def _get_hdr(self):
0208         return " %7d : %s  " % ( self.pindex, self.path )
0209     hdr = property(_get_hdr)
0210
0211     def __str__(self):
0212         return "\n".join([self.hdr]+map(str, self))
0213
0214     def __repr__(self):
0215         return "\n".join([self.hdr]+map(repr, self))
0216
0217
0218 if __name__ == '__main__':
0219
0220
0221     pindex = int(sys.argv[1]) if len(sys.argv) > 1 else 1230
0222
0223     import numpy as np
0224
0225     rng = np.load(UCF.rngpath())
0226     xrng = rng[pindex].ravel()
0227     #print_(str(xrng))
0228
0229     ## workaround lack of numpy in system python used by lldb
0230     ## by writing rng for the pindex to txt file
0231     trng = UCF.rngpathtxt(pindex)
0232     np.savetxt(trng, xrng, delimiter=",")
0233     xrng2 = np.array(UCF.loadrngtxt(pindex), dtype=np.float64)
0234     #print_(xrng2)
0235     assert np.all(xrng2 == xrng)
0236
0237     ucf = UCF( pindex )
0238
0239     #print_(str(ucf))
0240     print_(repr(ucf))
0241
0242