eic-opticks/ana/truncation.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 truncation.py
0023 ===============
0024
0025 Pseudo duplication of oxrap/cu/generate.cu for easy thinking about truncation,
0026 and implementing seqhis steered propagations.
0027
0028 record_max
0029      has a hard limit at 16, from 4bits*16 = 64 bits
0030      hopefully before this becomes a problem GPU/CUDA/OptiX
0031      will natively support 128bit uints
0032
0033      There are workarounds using uint4 and inline PTX asm assembly however
0034
0035      * :google:`cuda 128 bit int`
0036      * http://stackoverflow.com/questions/6162140/128-bit-integer-on-cuda
0037
0038 bounce_max
0039      is more flexible, but for clean truncation using bounce_max = record_max - 1
0040      is highly favored (ie maximally 15)
0041
0042 bounce
0043      not really the right good word as SC,AB,RE dont involve a "bounce",
0044      trace is more appropriate as a raytrace is always done for every turn of the loop
0045
0046
0047 Why use bounce at all, why not just use slot in the loop ?
0048    as step-by-step recording is just a temporary measure whilst debugging,
0049    the real thing in the photon that gets recorded after the bounce loop
0050
0051
0052 """
0053 import logging
0054 log = logging.getLogger(__name__)
0055 index = 0
0056
0057 BREAK = "BREAK"
0058 CONTINUE = "CONTINUE"
0059 PASS = "PASS"
0060
0061 #SEQHIS = "TO BT BT SC BT BT SA"
0062 #SEQHIS = "TO BT BT BT BT DR BT BT BT BT BT BT BT BT SA"
0063 SEQHIS = "TO BT BT BT BT DR BT BT BT BT SC BT BT BT BT SA"
0064 SQ = SEQHIS.split()
0065 tru = 0
0066
0067 def truth():
0068     global tru
0069     try:
0070        sq = SQ[tru]
0071     except IndexError:
0072        sq = None
0073     pass
0074     tru += 1
0075     return sq
0076
0077
0078 def RSAVE(msg, slot, slot_offset, bounce):
0079     global index
0080     log.info(" RSAVE  [%2d] (%10s) slot %2d slot_offset %2d bounce %2d " % (index, msg, slot, slot_offset, bounce))
0081     index += 1
0082
0083
0084 def propagate_to_boundary(sq):
0085     if sq == "AB":
0086         rc = BREAK
0087     elif sq == "RE" or sq == "SC":
0088         rc = CONTINUE
0089     else:
0090         rc = PASS
0091     pass
0092     log.info(" propagate_to_boundary sq %s rc %s " % (sq, rc) )
0093     return rc
0094
0095
0096 def propagate_at_surface(sq):
0097     if sq in "SD SA".split():
0098         rc = BREAK
0099     elif sq in "DR SR".split():
0100         rc = CONTINUE
0101     else:
0102         assert False, sq
0103
0104     log.info(" propagate_at_surface sq %s rc %s " % (sq, rc) )
0105     return rc
0106
0107 def propagate_at_boundary(sq):
0108     assert sq in "BR BT".split()
0109     rc = CONTINUE
0110     log.info(" propagate_at_boundary sq %s rc %s " % (sq, rc) )
0111     return rc
0112
0113
0114 def generate():
0115     record_max = 16
0116     bounce_max = 15
0117     MAXREC = record_max
0118
0119     photon_id = 1000
0120     slot_min = photon_id*MAXREC
0121     slot_max = slot_min + MAXREC - 1
0122     log.info(" SEQHIS %s   photon_id %7d slot_min %2d slot_max %2d " % (SEQHIS, photon_id, slot_min, slot_max))
0123
0124     bounce = 0
0125     slot = 0
0126
0127     flag = truth()
0128     assert flag == 'TO'
0129
0130
0131     while bounce < bounce_max:
0132
0133         print "\n"
0134         # bounce here up to bounce_max - 1
0135         bounce += 1
0136         # bounce here up to bounce_max
0137
0138         # raytrace happens here
0139         if slot < MAXREC:
0140             slot_offset = slot_min + slot
0141         else:
0142             slot_offset = slot_max ;
0143         pass
0144         RSAVE("inside", slot, slot_offset, bounce)
0145         slot += 1
0146
0147         flag = truth()
0148
0149         command = propagate_to_boundary(flag)
0150         if command == BREAK: break
0151         if command == CONTINUE: continue
0152         if command == PASS: pass
0153
0154         if flag in "SD SA DR SR".split():
0155             command = propagate_at_surface(flag)
0156             if command == BREAK: break
0157             if command == CONTINUE: continue
0158         elif flag in "BR BT".split():
0159             propagate_at_boundary(flag)
0160         else:
0161             assert 0, flag
0162     pass
0163
0164
0165     if bounce == bounce_max:
0166         log.warning(" exited while with: bounce == bounce_max (%d) ... about to post-loop write into topslot (0-based) (%d) " % (bounce, slot))
0167
0168
0169     if slot < MAXREC:
0170         slot_offset = slot_min + slot
0171     else:
0172         slot_offset = slot_max ;
0173     pass
0174     RSAVE("after", slot, slot_offset, bounce)
0175     pass
0176
0177
0178
0179
0180 if __name__ == '__main__':
0181     logging.basicConfig(level=logging.INFO)
0182     generate()
0183