File indexing completed on 2025-08-06 08:10:58
0001
0002
0003 import os
0004 import argparse
0005 import time
0006 from datetime import datetime
0007
0008 from orion.client import build_experiment
0009 from orion.storage.base import get_storage
0010
0011 import acts
0012 from acts import (
0013 SurfaceMaterialMapper,
0014 VolumeMaterialMapper,
0015 Navigator,
0016 Propagator,
0017 StraightLineStepper,
0018 MaterialMapJsonConverter,
0019 )
0020 from acts.examples import (
0021 Sequencer,
0022 WhiteBoard,
0023 AlgorithmContext,
0024 ProcessCode,
0025 RootMaterialTrackReader,
0026 MaterialMapping,
0027 JsonMaterialWriter,
0028 JsonFormat,
0029 )
0030 from acts.examples.odd import getOpenDataDetector
0031
0032
0033 def runMaterialMappingNoTrack(
0034 trackingGeometry,
0035 decorators,
0036 outputDir,
0037 inputDir,
0038 mapName="material-map",
0039 mapSurface=True,
0040 mapVolume=True,
0041 format=JsonFormat.Json,
0042 readCachedSurfaceInformation=False,
0043 s=None,
0044 ):
0045 """
0046 Implementation of the material mapping that doesn't write the material tracks.
0047 Used to create the material map that will then be used to compute the material variance.
0048
0049 trackingGeometry : The tracking geometry
0050 decorators : The decorators for the tracking geometry
0051 outputDir : Output directory for the material map
0052 inputDir : Input directory containing the Material track
0053 mapName : Name of the material map
0054 mapSurface : Is material being mapped onto surfaces ?
0055 mapVolume : Is material being mapped onto volumes ?
0056 format : Json format used to write the material map (json, cbor, ...)
0057 readCachedSurfaceInformation : If set to true it will be assumed that the surface has already been associated with each material interaction in the input file.
0058 """
0059
0060 s = s or Sequencer(numThreads=1)
0061 for decorator in decorators:
0062 s.addContextDecorator(decorator)
0063 wb = WhiteBoard(acts.logging.INFO)
0064 context = AlgorithmContext(0, 0, wb)
0065 for decorator in decorators:
0066 assert decorator.decorate(context) == ProcessCode.SUCCESS
0067
0068
0069 s.addReader(
0070 RootMaterialTrackReader(
0071 level=acts.logging.INFO,
0072 outputMaterialTracks="material-tracks",
0073 fileList=[
0074 os.path.join(
0075 inputDir,
0076 "optimised-material-map_tracks.root"
0077 if readCachedSurfaceInformation
0078 else "geant4_material_tracks.root",
0079 )
0080 ],
0081 readCachedSurfaceInformation=readCachedSurfaceInformation,
0082 )
0083 )
0084
0085 stepper = StraightLineStepper()
0086 mmAlgCfg = MaterialMapping.Config(context.geoContext, context.magFieldContext)
0087 mmAlgCfg.trackingGeometry = trackingGeometry
0088 mmAlgCfg.inputMaterialTracks = "material-tracks"
0089
0090 if mapSurface:
0091 navigator = Navigator(
0092 trackingGeometry=trackingGeometry,
0093 resolveSensitive=True,
0094 resolveMaterial=True,
0095 resolvePassive=True,
0096 )
0097 propagator = Propagator(stepper, navigator)
0098 mapper = SurfaceMaterialMapper(level=acts.logging.INFO, propagator=propagator)
0099 mmAlgCfg.materialSurfaceMapper = mapper
0100
0101 if mapVolume:
0102 navigator = Navigator(
0103 trackingGeometry=trackingGeometry,
0104 )
0105 propagator = Propagator(stepper, navigator)
0106 mapper = VolumeMaterialMapper(
0107 level=acts.logging.INFO, propagator=propagator, mappingStep=999
0108 )
0109 mmAlgCfg.materialVolumeMapper = mapper
0110
0111 jmConverterCfg = MaterialMapJsonConverter.Config(
0112 processSensitives=True,
0113 processApproaches=True,
0114 processRepresenting=True,
0115 processBoundaries=True,
0116 processVolumes=True,
0117 context=context.geoContext,
0118 )
0119
0120 jmw = JsonMaterialWriter(
0121 level=acts.logging.VERBOSE,
0122 converterCfg=jmConverterCfg,
0123 fileName=os.path.join(outputDir, mapName),
0124 writeFormat=format,
0125 )
0126
0127 mmAlgCfg.materialWriters = [jmw]
0128
0129 s.addAlgorithm(MaterialMapping(level=acts.logging.INFO, config=mmAlgCfg))
0130
0131 return s
0132
0133
0134 def runMaterialMappingVariance(
0135 binMap,
0136 events,
0137 job,
0138 inputPath,
0139 pathExp,
0140 pipeResult,
0141 readCachedSurfaceInformation=False,
0142 ):
0143 """
0144 Run the material mapping and compute the variance for each bin of each surface
0145 Return a dict with the GeometryId value of the surface as a key that stores
0146 a list of pairs corresponding to the variance and number of tracks associated with each bin of the surface
0147
0148 binMap : Map containing the binning for each surface
0149 events : Number of event to use in the mapping
0150 job : ID of the job
0151 inputPath : Directory containing the input geantino track and the json geometry
0152 pathExp : Material mapping optimisation path
0153 pipeResult : Pipe to send back the score to the main python instance
0154 readCachedSurfaceInformation : Are surface information stored in the material track. Switch to true if the mapping was already performed to improve the speed.
0155 """
0156 print(
0157 datetime.now().strftime("%H:%M:%S") + " Start mapping for job " + str(job),
0158 flush=True,
0159 )
0160 mapName = "material-map-" + str(job)
0161 mapSurface = True
0162 mapVolume = False
0163
0164
0165 matDeco = acts.IMaterialDecorator.fromFile(
0166 str(os.path.join(inputPath, "geometry-map.json"))
0167 )
0168 detectorTemp, trackingGeometryTemp, decoratorsTemp = getOpenDataDetector(matDeco)
0169 matMapDeco = acts.MappingMaterialDecorator(
0170 tGeometry=trackingGeometryTemp, level=acts.logging.ERROR
0171 )
0172
0173 matMapDeco.setBinningMap(binMap)
0174
0175 del detectorTemp
0176 del trackingGeometryTemp
0177 del decoratorsTemp
0178
0179
0180 detector, trackingGeometry, decorators = getOpenDataDetector(matMapDeco)
0181
0182
0183 sMap = acts.examples.Sequencer(
0184 events=events, numThreads=1, logLevel=acts.logging.INFO
0185 )
0186
0187
0188 runMaterialMappingNoTrack(
0189 trackingGeometry,
0190 decorators,
0191 outputDir=pathExp,
0192 inputDir=inputPath,
0193 mapName=mapName,
0194 format=JsonFormat.Cbor,
0195 mapVolume=mapVolume,
0196 readCachedSurfaceInformation=readCachedSurfaceInformation,
0197 s=sMap,
0198 )
0199 sMap.run()
0200 del sMap
0201 del detector
0202 del trackingGeometry
0203 del decorators
0204
0205
0206 print(
0207 datetime.now().strftime("%H:%M:%S")
0208 + " Job "
0209 + str(job)
0210 + ": second pass to compute the variance",
0211 flush=True,
0212 )
0213
0214 cborMap = os.path.join(pathExp, (mapName + ".cbor"))
0215 matDecoVar = acts.IMaterialDecorator.fromFile(cborMap)
0216 detectorVar, trackingGeometryVar, decoratorsVar = getOpenDataDetector(matDecoVar)
0217 s = acts.examples.Sequencer(events=events, numThreads=1, logLevel=acts.logging.INFO)
0218 for decorator in decoratorsVar:
0219 s.addContextDecorator(decorator)
0220 wb = WhiteBoard(acts.logging.ERROR)
0221 context = AlgorithmContext(0, 0, wb)
0222 for decorator in decoratorsVar:
0223 assert decorator.decorate(context) == ProcessCode.SUCCESS
0224
0225
0226 reader = RootMaterialTrackReader(
0227 level=acts.logging.ERROR,
0228 outputMaterialTracks="material-tracks",
0229 fileList=[
0230 os.path.join(
0231 inputPath,
0232 "optimised-material-map_tracks.root"
0233 if readCachedSurfaceInformation
0234 else "geant4_material_tracks.root",
0235 )
0236 ],
0237 readCachedSurfaceInformation=readCachedSurfaceInformation,
0238 )
0239 s.addReader(reader)
0240 stepper = StraightLineStepper()
0241 mmAlgCfg = MaterialMapping.Config(context.geoContext, context.magFieldContext)
0242 mmAlgCfg.trackingGeometry = trackingGeometryVar
0243 mmAlgCfg.inputMaterialTracks = "material-tracks"
0244
0245 if mapSurface:
0246 navigator = Navigator(
0247 trackingGeometry=trackingGeometryVar,
0248 resolveSensitive=True,
0249 resolveMaterial=True,
0250 resolvePassive=True,
0251 )
0252 propagator = Propagator(stepper, navigator)
0253 surfaceCfg = SurfaceMaterialMapper.Config(
0254 computeVariance=True
0255 )
0256 mapper = SurfaceMaterialMapper(
0257 config=surfaceCfg, level=acts.logging.ERROR, propagator=propagator
0258 )
0259 mmAlgCfg.materialSurfaceMapper = mapper
0260
0261 if mapVolume:
0262 navigator = Navigator(
0263 trackingGeometry=trackingGeometryVar,
0264 )
0265 propagator = Propagator(stepper, navigator)
0266 mapper = VolumeMaterialMapper(
0267 level=acts.logging.ERROR, propagator=propagator, mappingStep=999
0268 )
0269 mmAlgCfg.materialVolumeMapper = mapper
0270
0271 mapping = MaterialMapping(level=acts.logging.ERROR, config=mmAlgCfg)
0272 s.addAlgorithm(mapping)
0273 s.run()
0274
0275
0276 score = dict()
0277 for key in binMap:
0278 objective = 0
0279 nonZero = 0
0280 binParameters = mapping.scoringParameters(key)
0281
0282
0283 for parameters in binParameters:
0284 if parameters[1] != 0:
0285 objective += parameters[0] * (1.0 + 1.0 / parameters[1])
0286 nonZero += 1
0287 if nonZero != 0:
0288 objective = objective / nonZero
0289 score[key] = [dict(name="surface_score", type="objective", value=objective)]
0290 print(
0291 datetime.now().strftime("%H:%M:%S")
0292 + " Mapping over for job "
0293 + str(job)
0294 + " : now sending score",
0295 flush=True,
0296 )
0297 pipeResult.send(score)
0298
0299 del mapping
0300 del s
0301 del detectorVar
0302 del trackingGeometryVar
0303 del decoratorsVar
0304 os.remove(cborMap)
0305
0306
0307 def surfaceExperiment(key, nbJobs, pathDB, pathResult, pipeBin, pipeResult, doPloting):
0308 """
0309 This function create an experiment for a given single surface
0310 Due to how Orion is implemented only one DB can exist per job, this thus need to be call using pythons multiprocessing to circumvent the issue.
0311
0312 key : Id of the surface corresponding to this experiment
0313 nbJobs : Total number of jobs to be executed simultaneously
0314 pathDB : Path to the databases
0315 pathResult : Path to write the result of the optimisation
0316 pipeBin : Pipe use to send the experiment binning to the main python instance
0317 pipeResult : Pipe to receive the result of the optimisation
0318 doPloting : true if we want to plot the result of the optimisation and obtain the optimal material map
0319 """
0320
0321 storage = {
0322 "database": {
0323 "name": "database_" + str(key),
0324 "type": "pickleddb",
0325 "host": os.path.join(pathDB, "database_" + str(key) + ".pkl"),
0326 "timeout": 43200,
0327 },
0328 }
0329
0330
0331
0332 space = {
0333 "x": "uniform(1, 240, discrete=True)",
0334 "y": "uniform(1, 240, discrete=True)",
0335 }
0336
0337 experiments = build_experiment(
0338 "s_" + str(key),
0339 version="1",
0340 space=space,
0341 algorithms="random",
0342 storage=storage,
0343 max_idle_time=43200,
0344 )
0345
0346 store = get_storage()
0347 store.delete_trials(uid=experiments.id, where={"status": {"$ne": "completed"}})
0348 store.release_algorithm_lock(uid=experiments.id)
0349
0350
0351 trials = dict()
0352 binMap = dict()
0353 for job in range(nbJobs):
0354 trials[job] = experiments.suggest()
0355 binMap[job] = (trials[job].params["x"], trials[job].params["y"])
0356 print(
0357 datetime.now().strftime("%H:%M:%S")
0358 + " Binning for job "
0359 + str(job)
0360 + " and surface "
0361 + str(key)
0362 + " has been selected",
0363 flush=True,
0364 )
0365 pipeBin.send(binMap[job])
0366 print(
0367 datetime.now().strftime("%H:%M:%S")
0368 + " All binning for surface "
0369 + str(key)
0370 + " has been sent",
0371 flush=True,
0372 )
0373
0374 for job in range(nbJobs):
0375 score = pipeResult.recv()
0376 print(
0377 datetime.now().strftime("%H:%M:%S")
0378 + " Received score for job "
0379 + str(job)
0380 + " and surface "
0381 + str(key),
0382 flush=True,
0383 )
0384 experiments.observe(trials[job], score)
0385 print(
0386 datetime.now().strftime("%H:%M:%S")
0387 + " Score for job "
0388 + str(job)
0389 + " and surface "
0390 + str(key)
0391 + " has been written",
0392 flush=True,
0393 )
0394
0395
0396 if doPloting:
0397 print(
0398 datetime.now().strftime("%H:%M:%S")
0399 + " All the jobs are over. Now creating the optimisation plots",
0400 flush=True,
0401 )
0402
0403 pathExpSurface = os.path.join(pathResult, "b_" + str(key))
0404 if not os.path.isdir(pathExpSurface):
0405 os.makedirs(pathExpSurface)
0406
0407 regret = experiments.plot.regret()
0408 regret.write_html(pathExpSurface + "/regret.html")
0409
0410 parallel_coordinates = experiments.plot.parallel_coordinates()
0411 parallel_coordinates.write_html(pathExpSurface + "/parallel_coordinates.html")
0412
0413 lpi = experiments.plot.lpi()
0414 lpi.write_html(pathExpSurface + "/lpi.html")
0415
0416 partial_dependencies = experiments.plot.partial_dependencies()
0417 partial_dependencies.write_html(pathExpSurface + "/partial_dependencies.html")
0418
0419
0420 df = experiments.to_pandas()
0421 best = df.iloc[df.objective.idxmin()]
0422 print(
0423 datetime.now().strftime("%H:%M:%S")
0424 + " Best score for surface "
0425 + str(key)
0426 + " : "
0427 + str(best),
0428 flush=True,
0429 )
0430 resultBinMap = (best.x, best.y)
0431 pipeBin.send(resultBinMap)
0432
0433
0434 if "__main__" == __name__:
0435 print(datetime.now().strftime("%H:%M:%S") + " Starting")
0436
0437 parser = argparse.ArgumentParser()
0438 parser.add_argument(
0439 "--numberOfJobs", nargs="?", default=2, type=int
0440 )
0441 parser.add_argument(
0442 "--topNumberOfEvents", nargs="?", default=100, type=int
0443 )
0444 parser.add_argument(
0445 "--inputPath", nargs="?", default=os.getcwd(), type=str
0446 )
0447 parser.add_argument(
0448 "--outputPath", nargs="?", default="", type=str
0449 )
0450 parser.add_argument(
0451 "--doPloting", action="store_true"
0452 )
0453 parser.add_argument(
0454 "--readCachedSurfaceInformation", action="store_true"
0455 )
0456 parser.set_defaults(doPloting=False)
0457 parser.set_defaults(readCachedSurfaceInformation=False)
0458 args = parser.parse_args()
0459
0460
0461 pathExp = os.path.join(args.outputPath, "Mapping")
0462 pathDB = os.path.join(pathExp, "Database")
0463 pathResult = os.path.join(pathExp, "Result")
0464 if not os.path.isdir(pathExp):
0465 os.makedirs(pathExp)
0466 if not os.path.isdir(pathDB):
0467 os.makedirs(pathDB)
0468 if not os.path.isdir(pathResult):
0469 os.makedirs(pathResult)
0470
0471
0472 matDeco = acts.IMaterialDecorator.fromFile(
0473 str(os.path.join(args.inputPath, "geometry-map.json"))
0474 )
0475 detector, trackingGeometry, decorators = getOpenDataDetector(matDeco)
0476
0477
0478 matMapDeco = acts.MappingMaterialDecorator(
0479 tGeometry=trackingGeometry, level=acts.logging.WARNING
0480 )
0481 binDict = matMapDeco.binningMap()
0482 del detector, decorators
0483
0484
0485 from multiprocessing import Process, Pipe
0486
0487 binPipes_child = dict()
0488 resultPipes_child = dict()
0489 scorePipes_child = dict()
0490
0491 binPipes_parent = dict()
0492 resultPipes_parent = dict()
0493 scorePipes_parent = dict()
0494
0495 expJob = dict()
0496 OptiJob = dict()
0497
0498
0499
0500
0501 for key in binDict:
0502 binPipes_parent[key], binPipes_child[key] = Pipe()
0503 scorePipes_parent[key], scorePipes_child[key] = Pipe()
0504 expJob[key] = Process(
0505 target=surfaceExperiment,
0506 args=(
0507 key,
0508 args.numberOfJobs,
0509 pathDB,
0510 pathResult,
0511 binPipes_child[key],
0512 scorePipes_child[key],
0513 args.doPloting,
0514 ),
0515 )
0516 expJob[key].start()
0517
0518
0519 for job in range(args.numberOfJobs):
0520 resultPipes_parent[job], resultPipes_child[job] = Pipe()
0521 binMap = dict()
0522
0523 for key in binDict:
0524 binMap[key] = binPipes_parent[key].recv()
0525 print(
0526 datetime.now().strftime("%H:%M:%S")
0527 + " Binning for job "
0528 + str(job)
0529 + " have been selected, now running the mapping",
0530 flush=True,
0531 )
0532
0533 OptiJob[job] = Process(
0534 target=runMaterialMappingVariance,
0535 args=(
0536 binMap,
0537 args.topNumberOfEvents,
0538 job,
0539 args.inputPath,
0540 pathExp,
0541 resultPipes_child[job],
0542 args.readCachedSurfaceInformation,
0543 ),
0544 )
0545 OptiJob[job].start()
0546
0547
0548 for job in range(args.numberOfJobs):
0549 scores = resultPipes_parent[job].recv()
0550 print(
0551 datetime.now().strftime("%H:%M:%S")
0552 + " Retried score for job "
0553 + str(job),
0554 flush=True,
0555 )
0556 for key in binDict:
0557 score = scores[key]
0558 scorePipes_parent[key].send(score)
0559 print(
0560 datetime.now().strftime("%H:%M:%S")
0561 + " Job number "
0562 + str(job)
0563 + " is over",
0564 flush=True,
0565 )
0566
0567 if args.doPloting:
0568
0569
0570 print(
0571 datetime.now().strftime("%H:%M:%S")
0572 + " Running the material mapping to obtain the optimised material map",
0573 flush=True,
0574 )
0575 resultBinMap = dict()
0576 for key in binDict:
0577 resultBinMap[key] = binPipes_parent[key].recv()
0578 matMapDeco.setBinningMap(resultBinMap)
0579
0580
0581 resultDetector, resultTrackingGeometry, resultDecorators = getOpenDataDetector(
0582 matMapDeco
0583 )
0584
0585
0586 rMap = acts.examples.Sequencer(
0587 events=args.topNumberOfEvents, numThreads=1, logLevel=acts.logging.INFO
0588 )
0589
0590
0591 from material_mapping import runMaterialMapping
0592
0593 runMaterialMapping(
0594 resultTrackingGeometry,
0595 resultDecorators,
0596 outputDir=args.outputPath,
0597 inputDir=args.inputPath,
0598 mapName="optimised-material-map",
0599 format=JsonFormat.Json,
0600 mapVolume=False,
0601 s=rMap,
0602 )
0603 rMap.run()
0604 del rMap
0605 del resultDetector, resultTrackingGeometry, resultDecorators
0606 print(
0607 datetime.now().strftime("%H:%M:%S")
0608 + " Waiting for all the score to have been stored",
0609 flush=True,
0610 )
0611
0612
0613 deadline = time.time() + 600
0614
0615 for key in binDict:
0616 timeout = max(deadline - time.time(), 0)
0617 expJob[key].join(timeout=timeout)
0618 expJob[key].terminate()
0619 print(
0620 datetime.now().strftime("%H:%M:%S")
0621 + " Experiment for surface "
0622 + str(key)
0623 + " is over",
0624 flush=True,
0625 )