diff --git a/macro/rich/Import_GDML_Export_ROOT.c b/macro/rich/Import_GDML_Export_ROOT.c index 568df46f26a5ebab0b7e7bc55db065fcfa1cc0ed..af7f7025d3f291d6205a7c84721cea8d768b35bc 100644 --- a/macro/rich/Import_GDML_Export_ROOT.c +++ b/macro/rich/Import_GDML_Export_ROOT.c @@ -1,10 +1,12 @@ void Import_GDML_Export_ROOT() { // we need to use latest root version to work with gdml geometry - system( - string("source /usr/local/Cellar/root/6.16.00/bin/thisroot.sh").c_str()); + // system(string("source /usr/local/Cellar/root/6.20.04_1/bin/thisroot.sh").c_str()); - TString richGeoFilename = "rich_v19d_pipe_v4.gdml"; + // For v17a geometries one need to use root5, one can use installation on lustre + //"source /cvmfs/fairroot.gsi.de/fairsoft/may16_root5/bin/thisroot.sh" + + TString richGeoFilename = "rich_v17a_1e_pcarb_bcarb.gdml"; TGeoManager* gdml = new TGeoManager("gdml", "FAIRGeom"); @@ -20,9 +22,9 @@ void Import_GDML_Export_ROOT() { // Starting from the version v18a position is defined inside the GDML file // Define your position HERE // Z coordinate for v16a = 270, for v17a = 258.75, for v18a = 0. - TGeoRotation* rot = new TGeoRotation("rot", 0., 0., 0.); - TGeoCombiTrans* posrot = - new TGeoCombiTrans(0., 0., 0., rot); // v16a - 270, v17a - 258.75, v18a - 0 + TGeoRotation* rot = new TGeoRotation("rot", 0., 0., 0.); + TGeoCombiTrans* posrot = new TGeoCombiTrans( + 0., 0., 258.75, rot); // v16a - 270, v17a - 258.75, v18a - 0 rootTop->AddNode(gdmlTop, 1, posrot); diff --git a/macro/rich/geotest/run_digi_geotest.C b/macro/rich/geotest/run_digi_geotest.C index 85cebfb4f80db7944a2222bc67d4ed3c4aec1177..ce3316552d0722b5c74bfc0d6672d166bc4a63bb 100644 --- a/macro/rich/geotest/run_digi_geotest.C +++ b/macro/rich/geotest/run_digi_geotest.C @@ -5,7 +5,7 @@ void run_digi_geotest( "/Users/slebedev/Development/cbm/data/sim/rich/geotest/param.00000.root", const string& digiFile = "/Users/slebedev/Development/cbm/data/sim/rich/geotest/digi.00000.root", - int nEvents = 1000) { + int nEvents = 10000) { FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); @@ -31,7 +31,7 @@ void run_digi_geotest( CbmRichDigitizer* richDigitizer = new CbmRichDigitizer(); richDigitizer->SetMaxNofHitsPerPmtCut(65); - run.SetDigitizer(kRich, richDigitizer, "RichPoint", true); + run.SetDigitizer(ECbmModuleId::kRich, richDigitizer, "RichPoint", true); run.Run(nEvents); diff --git a/macro/rich/geotest/run_digi_urqmdtest.C b/macro/rich/geotest/run_digi_urqmdtest.C index 2475836340a571927b3fbd6ba81062241f81e0d0..cbb473a7bd76e342067adec0e01c81238a37784c 100644 --- a/macro/rich/geotest/run_digi_urqmdtest.C +++ b/macro/rich/geotest/run_digi_urqmdtest.C @@ -5,7 +5,7 @@ void run_digi_urqmdtest( "/Users/slebedev/Development/cbm/data/sim/rich/urqmdtest/param.00000.root", const string& digiFile = "/Users/slebedev/Development/cbm/data/sim/rich/urqmdtest/digi.00000.root", - int nEvents = 5) { + int nEvents = 1000) { FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); diff --git a/macro/rich/geotest/run_qa_geotest.C b/macro/rich/geotest/run_qa_geotest.C index 57d8bd06c7914b8d7e47cf32215c9a734f1a33a7..001b083658b67fcf5403df17adc49fd8a5c6f85b 100644 --- a/macro/rich/geotest/run_qa_geotest.C +++ b/macro/rich/geotest/run_qa_geotest.C @@ -8,11 +8,11 @@ void run_qa_geotest( "/Users/slebedev/Development/cbm/data/sim/rich/geotest/digi.00000.root", const string& recoFile = "/Users/slebedev/Development/cbm/data/sim/rich/geotest/reco.00000.root", - const string& qaFile = - "/Users/slebedev/Development/cbm/data/sim/rich/geotest/qa.00000.root", + const string& qaFile = "/Users/slebedev/Development/cbm/data/sim/rich/" + "geotest/qa.g3new_media0.00000.root", const string& geoSetup = "sis100_electron", - const string& resultDir = "results_geotest_geant3_10/", - int nEvents = 1000) { + const string& resultDir = "results_geotest_geant3_new_media0/", + int nEvents = 10000) { FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); TTree::SetMaxTreeSize(90000000000); diff --git a/macro/rich/geotest/run_qa_urqmdtest.C b/macro/rich/geotest/run_qa_urqmdtest.C index 0d62d6a98613741400749420d95fa329ab76a4c2..52390d2f0c7be901efe01e3b8f53a14e84785e7d 100644 --- a/macro/rich/geotest/run_qa_urqmdtest.C +++ b/macro/rich/geotest/run_qa_urqmdtest.C @@ -10,8 +10,8 @@ void run_qa_urqmdtest( const string& qaFile = "/Users/slebedev/Development/cbm/data/sim/rich/urqmdtest/qa.00000.root", const string& geoSetup = "sis100_electron", - const string& resultDir = "results_urqmdtest/", - int nEvents = 5) { + const string& resultDir = "results_urqmdtest_geant4/", + int nEvents = 1000) { FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); diff --git a/macro/rich/geotest/run_reco_geotest.C b/macro/rich/geotest/run_reco_geotest.C index d6281c326b04e727db54e562467753a867b9967e..f84c2dc67ada0153faad8ece9b924bd38d15705a 100644 --- a/macro/rich/geotest/run_reco_geotest.C +++ b/macro/rich/geotest/run_reco_geotest.C @@ -10,7 +10,7 @@ void run_reco_geotest( "/Users/slebedev/Development/cbm/data/sim/rich/geotest/reco.00000.root", const string& geoSetup = "sis100_electron", const string& resultDir = "rich_pipe_v1", // "results_geotest_test/", - int nEvents = 1000) { + int nEvents = 10000) { FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); TTree::SetMaxTreeSize(90000000000); diff --git a/macro/rich/geotest/run_reco_urqmdtest.C b/macro/rich/geotest/run_reco_urqmdtest.C index c078e4b0f5f7def91539aeda5ad2731d93630789..87699756ccf879fec72dca02b3c040e8676b558c 100644 --- a/macro/rich/geotest/run_reco_urqmdtest.C +++ b/macro/rich/geotest/run_reco_urqmdtest.C @@ -9,7 +9,7 @@ void run_reco_urqmdtest( "/Users/slebedev/Development/cbm/data/sim/rich/urqmdtest/reco.00000.root", const string& geoSetup = "sis100_electron", const string& resultDir = "results_urqmdtest/", - int nEvents = 5) { + int nEvents = 1000) { FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); diff --git a/macro/rich/geotest/run_sim_geotest.C b/macro/rich/geotest/run_sim_geotest.C index e2f423f242c1e938a793574c63eb1db03af06e5e..c5e4971de4d0d20c369a34249bc5c4edd6dd002e 100644 --- a/macro/rich/geotest/run_sim_geotest.C +++ b/macro/rich/geotest/run_sim_geotest.C @@ -7,7 +7,7 @@ void run_sim_geotest( const string& geoFile = "/Users/slebedev/Development/cbm/data/sim/rich/geotest/geosim.00000.root", const string& geoSetup = "sis100_electron", //"mirror12_42", - int nEvents = 1000) { + int nEvents = 10) { TTree::SetMaxTreeSize(90000000000); remove(parFile.c_str()); @@ -46,7 +46,7 @@ void run_sim_geotest( run.SetTarget("Gold", 0.025, 2.5); run.SetBeamPosition(0., 0., 0.1, 0.1); //run.SetEngine(kGeant4); - run.StoreTrajectories(true); + //run.StoreTrajectories(true); run.Run(nEvents); timer.Stop(); diff --git a/macro/rich/geotest/run_sim_urqmdtest.C b/macro/rich/geotest/run_sim_urqmdtest.C index 89558a46932cf75ec51b36dc398c48b837f32932..dff5d829ab5cb60a96cde60bda16a420bb4b1b4d 100644 --- a/macro/rich/geotest/run_sim_urqmdtest.C +++ b/macro/rich/geotest/run_sim_urqmdtest.C @@ -8,7 +8,7 @@ void run_sim_urqmdtest( const string& geoFile = "/Users/slebedev/Development/cbm/data/sim/rich/urqmdtest/geosim.00000.root", const string& geoSetup = "sis100_electron", - int nEvents = 5) { + int nEvents = 1000) { TTree::SetMaxTreeSize(90000000000); remove(parFile.c_str()); @@ -26,6 +26,7 @@ void run_sim_urqmdtest( run.LoadSetup(geoSetup.c_str()); run.SetTarget("Gold", 0.025, 2.5); run.SetBeamPosition(0., 0., 0.1, 0.1); + run.SetEngine(kGeant4); run.Run(nEvents); timer.Stop(); diff --git a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis2.C b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis2.C index 9e6a6bc713a4df126a806b015127f255ed2a5859..06dacaf5476f82bd2dc5eba138f689d80c3efa8b 100644 --- a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis2.C +++ b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis2.C @@ -55,7 +55,7 @@ void run_analysis2() { fgCalibrator->SetCalibrationPeriod(100000); fgCalibrator->SetInputFilename( "calibration.root"); // does not actually import data - only defines - // the file that will be used if you specify mode etn_IMPORT + // the file that will be used if you specify mode etn_IMPORT fgCalibrator->SetMode(etn_ONLINE); // Also note the (un)commented line in the end of the macro with export func diff --git a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis3.C b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis3.C index dc9d8bd37c6a512f422b634a1de7c2e3e5bd6055..39b43df21d3a3df9f74c02983c81259a5577ac4d 100644 --- a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis3.C +++ b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis3.C @@ -53,7 +53,7 @@ void run_analysis3() { fgCalibrator->SetCalibrationPeriod(50000000); fgCalibrator->SetInputFilename( "calibration_flib_1040.root"); // does not actually import data - only defines - // the file that will be used if you specidy mode etn_IMPORT + // the file that will be used if you specidy mode etn_IMPORT fgCalibrator->SetMode(etn_ONLINE); // Also note the (un)commented line in the end of the macro with export func diff --git a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu.C b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu.C index be30f0d669dc99d67c6c201a42490832e350827e..faf5dccf3718d7f3da3b5d7415638aa9b9a7fe6a 100644 --- a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu.C +++ b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu.C @@ -196,8 +196,8 @@ void run_analysis_wu(Bool_t generateCalib = kFALSE, Int_t inChannel = 1) { fgCalibrator->SetCalibrationPeriod(50000000); fgCalibrator->SetInputFilename( "calibration_wu.root"); // does not actually import data - only defines - // the file that will be used if you specify mode etn_IMPORT - // Also note the (un)commented line in the end of the macro with export func + // the file that will be used if you specify mode etn_IMPORT + // Also note the (un)commented line in the end of the macro with export func if (generateCalib) { fgCalibrator->SetMode(etn_ONLINE); } else { diff --git a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu_2.C b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu_2.C index ee185c4a408023a5d0b750e87d2bd9b81138713b..d3f4994207d619789f3b8e62d9d8823df70a607c 100644 --- a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu_2.C +++ b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu_2.C @@ -165,8 +165,8 @@ void run_analysis_wu_2(Bool_t generateCalib = kFALSE, fgCalibrator->SetCalibrationPeriod(50000000); fgCalibrator->SetInputFilename( "calibration_wu.root"); // does not actually import data - only defines - // the file that will be used if you specify mode etn_IMPORT - // Also note the (un)commented line in the end of the macro with export func + // the file that will be used if you specify mode etn_IMPORT + // Also note the (un)commented line in the end of the macro with export func if (generateCalib) { fgCalibrator->SetMode(etn_ONLINE); } else { diff --git a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu_2_ip.C b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu_2_ip.C index 35caad619a951b0901f7c0b17b0750b1bb9a2fa2..3b956d5a7690ce7a3315f6ba6261f67bf695fd21 100644 --- a/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu_2_ip.C +++ b/macro/rich/prototype/beamtime/NewUnpacker/run_analysis_wu_2_ip.C @@ -164,8 +164,8 @@ void run_analysis_wu_2_ip(Bool_t generateCalib = kFALSE, fgCalibrator->SetCalibrationPeriod(50000000); fgCalibrator->SetInputFilename( "calibration_wu.root"); // does not actually import data - only defines - // the file that will be used if you specify mode etn_IMPORT - // Also note the (un)commented line in the end of the macro with export func + // the file that will be used if you specify mode etn_IMPORT + // Also note the (un)commented line in the end of the macro with export func if (generateCalib) { fgCalibrator->SetMode(etn_ONLINE); } else { diff --git a/macro/rich/prototype/beamtime/run_analysis.C b/macro/rich/prototype/beamtime/run_analysis.C index ace526634bfeebe40bf4767166e59952ea7d727c..66812c631c0e5201b3888615c150f33f2565791f 100644 --- a/macro/rich/prototype/beamtime/run_analysis.C +++ b/macro/rich/prototype/beamtime/run_analysis.C @@ -63,7 +63,7 @@ void run_analysis() { CbmTrbCalibrator* fgCalibrator = CbmTrbCalibrator::Instance(); fgCalibrator->SetInputFilename( "calibration.root"); // does not actually import data - only defines - // the file that will be used if you specidy mode etn_IMPORT + // the file that will be used if you specidy mode etn_IMPORT fgCalibrator->SetMode(etn_NOCALIB); // Also note the (un)commented line in the end of the macro with export func diff --git a/macro/rich/prototype/beamtime/run_analysis2.C b/macro/rich/prototype/beamtime/run_analysis2.C index f70d1d6f4c0dd178dabdd48e77df712ba56e929e..02d5a11a68e7100e2c2d61418f022d910e840c56 100644 --- a/macro/rich/prototype/beamtime/run_analysis2.C +++ b/macro/rich/prototype/beamtime/run_analysis2.C @@ -49,7 +49,7 @@ void run_analysis2() { //fgCalibrator->SetCalibrationPeriod(20000); fgCalibrator->SetInputFilename( "calibration2.root"); // does not actually import data - only defines - // the file that will be used if you specidy mode etn_IMPORT + // the file that will be used if you specidy mode etn_IMPORT fgCalibrator->SetMode(etn_IMPORT); // Also note the (un)commented line in the end of the macro with export func diff --git a/macro/rich/run/batch_job.py b/macro/rich/run/batch_job.py new file mode 100755 index 0000000000000000000000000000000000000000..0af067e4286c4ba5e37e140aac4286ecc050cc13 --- /dev/null +++ b/macro/rich/run/batch_job.py @@ -0,0 +1,45 @@ +#!/usr/bin/env python3 + +import os +import sys +import shutil + +dataDir = sys.argv[1] +geoSetup = sys.argv[2] +cbmrootConfigPath = "/lustre/nyx/cbm/users/slebedev/cbm/trunk/build/config.sh" +macroDir = "/lustre/nyx/cbm/users/slebedev/cbm/trunk/cbmroot/macro/rich/" +taskId = os.environ.get('SLURM_ARRAY_TASK_ID') +jobId = os.environ.get('SLURM_ARRAY_JOB_ID') + +print("dataDir:" + dataDir) + +os.system(("source {}").format(cbmrootConfigPath)) + +workDir = dataDir + "/workdir/" + jobId + "_" + taskId + "/" +if os.path.exists(workDir): + shutil.rmtree(workDir) +os.makedirs(workDir) +os.chdir(workDir) + + +runRecoqaUrqmd = True + +if runRecoqaUrqmd: + nofEvents = 1000 + dataDirNew = dataDir + "/recoqa_urqmd8gev/" + if not os.path.exists(dataDirNew): os.makedirs(dataDirNew) + urqmdFile = "/lustre/nyx/cbm/prod/gen/urqmd/auau/8gev/centr/urqmd.auau.8gev.centr." + str(taskId).zfill(5) + ".root" + mcFile = dataDirNew + "/mc."+ taskId + ".root" + parFile = dataDirNew + "/param."+ taskId + ".root" + digiFile = dataDirNew + "/digi."+ taskId + ".root" + recoFile = dataDirNew + "/reco."+ taskId + ".root" + qaFile = dataDirNew + "/qa."+ taskId + ".root" + geoSimFile = dataDirNew + "/geosim."+ taskId + ".root" + resultDir = dataDirNew + "/results/results_" + geoSetup + "/" + nofElectrons = 5 + nofPositrons = 5 + plutoFile = "" + os.system( ('root -l -b -q {}/run/run_sim.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},{},\\"{}\\",\\"{}\\",{}\)').format(macroDir, urqmdFile, mcFile, parFile, geoSimFile, nofElectrons, nofPositrons, plutoFile, geoSetup, nofEvents) ) + os.system( ('root -l -b -q {}/run/run_digi.C\(\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, mcFile, parFile, digiFile, nofEvents) ) + os.system( ('root -l -b -q {}/run/run_reco.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, mcFile, parFile, digiFile, recoFile, geoSetup, resultDir, nofEvents) ) + os.system( ('root -l -b -q {}/run/run_qa.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, mcFile, parFile, digiFile, recoFile, qaFile, geoSetup, resultDir, nofEvents) ) diff --git a/macro/rich/run/batch_send.py b/macro/rich/run/batch_send.py new file mode 100755 index 0000000000000000000000000000000000000000..3a49127984c13fcecba1e912fef7fce8652d033f --- /dev/null +++ b/macro/rich/run/batch_send.py @@ -0,0 +1,36 @@ +#!/usr/bin/env python3 + +import os +import shutil + +nofJobs = 100 +timeLimit="08:00:00" +geoSetup="sis100_electron_rich_pcarb_bcarb" +#All data in data dir will be removed +removeData = False +dataDir = "/lustre/nyx/cbm/users/slebedev/cbm/data/rich_carb/" + geoSetup + "/" +logFile = dataDir + "/log/log_slurm-%A_%a.out" +errorFile = dataDir + "/error/error_slurm-%A_%a.out" +jobName = "RICH" + +if removeData: + print("All data in dataDir will be removed. Dir:" + dataDir) + print("Removing...") + if os.path.exists(dataDir): + shutil.rmtree(dataDir) + os.makedirs(dataDir) + +if os.path.exists(os.path.dirname(logFile)): + shutil.rmtree(os.path.dirname(logFile)) +os.makedirs(os.path.dirname(logFile)) + +if os.path.exists(os.path.dirname(errorFile)): + shutil.rmtree(os.path.dirname(errorFile)) +os.makedirs(os.path.dirname(errorFile)) + +#-p debug +commandStr=('sbatch --job-name={} --time={} --output={} --error={} --array=1-{} batch_job.py {} {}').format(jobName, timeLimit, + logFile, errorFile, nofJobs, dataDir, geoSetup) +#print(commandStr) +os.system(commandStr) + diff --git a/macro/rich/run/draw_litqa.C b/macro/rich/run/draw_litqa.C index 77b895418f8c3007a39199068dbfb93682b25574..f4911aaddce522eecb4358448a1895ee83030f78 100644 --- a/macro/rich/run/draw_litqa.C +++ b/macro/rich/run/draw_litqa.C @@ -1,21 +1,19 @@ void draw_litqa() { - gROOT->LoadMacro("$VMCWORKDIR/macro/littrack/loadlibs.C"); - loadlibs(); + //gROOT->LoadMacro("$VMCWORKDIR/macro/littrack/loadlibs.C"); + //loadlibs(); - - //std::string dir = "/Users/slebedev/Development/cbm/data/lmvm/nov13/25gev/trd/1.0field/nomvd/rho0/"; - std::string outputDir = "results_litqa/"; - std::string fileName = "/Users/slebedev/Development/cbm/data/simulations/" - "rich/richreco/25gev.reco.0005.root"; + std::string outputDir = "rich_carb/results_litqa_ac/"; + //std::string fileName = "/Users/slebedev/Development/cbm/data/simulations/rich/richreco/25gev.reco.0005.root"; + std::string fileName = "qa.all.ac.root"; CbmSimulationReport* trackingQaReport = new CbmLitTrackingQaReport(); - trackingQaReport->Create(fileName, ""); //outputDir + "/tracking/"); + trackingQaReport->Create(fileName, outputDir + "/tracking/"); // CbmSimulationReport* fitQaReport = new CbmLitFitQaReport(); // fitQaReport->Create(fileName, outputDir + "/fit/"); - // CbmSimulationReport* clusteringQaReport = new CbmLitClusteringQaReport(); - // clusteringQaReport->Create(fileName, outputDir); + CbmSimulationReport* clusteringQaReport = new CbmLitClusteringQaReport(); + clusteringQaReport->Create(fileName, outputDir + "clustering/"); // CbmLitRadLengthQaReport* radLengthQaReport = new CbmLitRadLengthQaReport(); // radLengthQaReport->Create(fileName, outputDir); diff --git a/macro/rich/run/run_digi.C b/macro/rich/run/run_digi.C index 1d1dd30b65b52766aac7241b270cfc9892afb033..78721071d63b37607c4ba7828f1b0250c9a98d5c 100644 --- a/macro/rich/run/run_digi.C +++ b/macro/rich/run/run_digi.C @@ -5,7 +5,7 @@ void run_digi( "/Users/slebedev/Development/cbm/data/sim/rich/reco/param.00000.root", const string& digiFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/digi.00000.root", - int nEvents = 100) { + int nEvents = 10) { TTree::SetMaxTreeSize(90000000000); FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); @@ -30,7 +30,7 @@ void run_digi( CbmRichDigitizer* richDigitizer = new CbmRichDigitizer(); richDigitizer->SetMaxNofHitsPerPmtCut(65); - run.SetDigitizer(kRich, richDigitizer, "RichPoint", true); + run.SetDigitizer(ECbmModuleId::kRich, richDigitizer, "RichPoint", true); run.Run(nEvents); diff --git a/macro/rich/run/run_qa.C b/macro/rich/run/run_qa.C index f31124145a18bae321487402521f11f0939598f7..ab112466265986eaf953245b1897d11c0a8c8602 100644 --- a/macro/rich/run/run_qa.C +++ b/macro/rich/run/run_qa.C @@ -19,11 +19,7 @@ void run_qa( remove(qaFile.c_str()); - TString setupFile = srcDir + "/geometry/setup/setup_" + geoSetup + ".C"; - TString setupFunct = "setup_" + geoSetup + "()"; - gROOT->LoadMacro(setupFile); - gROOT->ProcessLine(setupFunct); - + CbmSetup::Instance()->LoadSetup(geoSetup.c_str()); std::cout << std::endl << "-I- " << myName << ": Defining parameter files " << std::endl; @@ -31,7 +27,7 @@ void run_qa( TString geoTag; // - TRD digitisation parameters - if (CbmSetup::Instance()->GetGeoTag(kTrd, geoTag)) { + if (CbmSetup::Instance()->GetGeoTag(ECbmModuleId::kTrd, geoTag)) { const Char_t* npar[4] = {"asic", "digi", "gas", "gain"}; TObjString* trdParFile(NULL); for (Int_t i(0); i < 4; i++) { @@ -44,7 +40,7 @@ void run_qa( } // - TOF digitisation parameters - if (CbmSetup::Instance()->GetGeoTag(kTof, geoTag)) { + if (CbmSetup::Instance()->GetGeoTag(ECbmModuleId::kTof, geoTag)) { TObjString* tofFile = new TObjString(srcDir + "/parameters/tof/tof_" + geoTag + ".digi.par"); parFileList->Add(tofFile); @@ -107,7 +103,7 @@ void run_qa( richCat.push_back("ElectronReference"); trackingQa->SetTrackCategories(trackCat); trackingQa->SetRingCategories(richCat); - // run->AddTask(trackingQa); + run->AddTask(trackingQa); CbmLitFitQa* fitQa = new CbmLitFitQa(); fitQa->SetMvdMinNofHits(0); @@ -115,15 +111,15 @@ void run_qa( fitQa->SetMuchMinNofHits(10); fitQa->SetTrdMinNofHits(2); fitQa->SetOutputDir(resultDir); - // run->AddTask(fitQa); + //run->AddTask(fitQa); CbmLitClusteringQa* clusteringQa = new CbmLitClusteringQa(); clusteringQa->SetOutputDir(resultDir); - // run->AddTask(clusteringQa); + run->AddTask(clusteringQa); CbmLitTofQa* tofQa = new CbmLitTofQa(); tofQa->SetOutputDir(std::string(resultDir)); - // run->AddTask(tofQa); + //run->AddTask(tofQa); std::cout << std::endl << std::endl @@ -151,7 +147,7 @@ void run_qa( timer.Stop(); std::cout << std::endl << std::endl; std::cout << "Macro finished succesfully." << std::endl; - std::cout << "Output file is " << recoFile << std::endl; + std::cout << "Output file is " << qaFile << std::endl; std::cout << "Parameter file is " << parFile << std::endl; std::cout << "Real time " << timer.RealTime() << " s, CPU time " << timer.CpuTime() << " s" << std::endl; diff --git a/macro/rich/run/run_reco.C b/macro/rich/run/run_reco.C index abeb7b38b8778da725bb8309b075f49d3e263c6b..7aff89be46f1ae0b599a40768c972cd7885e2f7d 100644 --- a/macro/rich/run/run_reco.C +++ b/macro/rich/run/run_reco.C @@ -17,11 +17,7 @@ void run_reco( remove(recoFile.c_str()); - TString setupFile = srcDir + "/geometry/setup/setup_" + geoSetup + ".C"; - TString setupFunct = "setup_" + geoSetup + "()"; - gROOT->LoadMacro(setupFile); - gROOT->ProcessLine(setupFunct); - + CbmSetup::Instance()->LoadSetup(geoSetup.c_str()); std::cout << std::endl << "-I- " << myName << ": Defining parameter files " << std::endl; @@ -29,7 +25,7 @@ void run_reco( TString geoTag; // - TRD digitisation parameters - if (CbmSetup::Instance()->GetGeoTag(kTrd, geoTag)) { + if (CbmSetup::Instance()->GetGeoTag(ECbmModuleId::kTrd, geoTag)) { const Char_t* npar[4] = {"asic", "digi", "gas", "gain"}; TObjString* trdParFile(NULL); for (Int_t i(0); i < 4; i++) { @@ -42,7 +38,7 @@ void run_reco( } // - TOF digitisation parameters - if (CbmSetup::Instance()->GetGeoTag(kTof, geoTag)) { + if (CbmSetup::Instance()->GetGeoTag(ECbmModuleId::kTof, geoTag)) { TObjString* tofFile = new TObjString(srcDir + "/parameters/tof/tof_" + geoTag + ".digi.par"); parFileList->Add(tofFile); diff --git a/macro/rich/run/run_sim.C b/macro/rich/run/run_sim.C index 1e4041b1f2e356654367c69ffab7f14772da016c..e71783ce444a024e6262921da388b85dbea859d7 100644 --- a/macro/rich/run/run_sim.C +++ b/macro/rich/run/run_sim.C @@ -14,7 +14,7 @@ void run_sim( 5, // number of e- to be generated, if <=0 no positrons are embedded into event const string& plutoFile = "", // if "", no pluto particles are embedded into event - const string& geoSetup = "sis100_electron", + const string& geoSetup = "sis100_electron_rich_pal_bcarb", int nEvents = 10) { TTree::SetMaxTreeSize(90000000000); diff --git a/reco/detectors/rich/utils/CbmRichUtil.h b/reco/detectors/rich/utils/CbmRichUtil.h index f40a9e36841bd8ca9bdd0faa5ea2e1c8a3feb4b5..b67b216af39633948e876e430dc96579de1c8b73 100644 --- a/reco/detectors/rich/utils/CbmRichUtil.h +++ b/reco/detectors/rich/utils/CbmRichUtil.h @@ -33,7 +33,6 @@ public: } static uint16_t GetDirichId(Int_t Address) { - return ((Address >> 16) & 0xFFFF); } diff --git a/sim/detectors/rich/CbmRich.cxx b/sim/detectors/rich/CbmRich.cxx index 7615721162866404096430737a2199e2a2c85d37..6de5fa5b78eb296672ac802a34b333b81a621f09 100644 --- a/sim/detectors/rich/CbmRich.cxx +++ b/sim/detectors/rich/CbmRich.cxx @@ -44,19 +44,12 @@ std::map<TString, TGeoMedium*> CbmRich::fFixedMedia; CbmRich::CbmRich() : FairDetector("RICH", kTRUE, ToIntegralType(ECbmModuleId::kRich)) , fPosIndex(0) - , fRegisterPhotonsOnSensitivePlane(false) + , fRegisterPhotonsOnSensitivePlane(true) , fRichPoints(new TClonesArray("CbmRichPoint")) , fRichRefPlanePoints(new TClonesArray("CbmRichPoint")) , fRichMirrorPoints(new TClonesArray("CbmRichPoint")) , fRotation(NULL) , fPositionRotation(NULL) { - /* - fRichPoints = ; - fRichRefPlanePoints = ; - fRichMirrorPoints = ; - fPosIndex = 0; - */ - fVerboseLevel = 1; } @@ -70,7 +63,7 @@ CbmRich::CbmRich(const char* name, Double_t rz) : FairDetector(name, active, ToIntegralType(ECbmModuleId::kRich)) , fPosIndex(0) - , fRegisterPhotonsOnSensitivePlane(false) + , fRegisterPhotonsOnSensitivePlane(true) , fRichPoints(new TClonesArray("CbmRichPoint")) , fRichRefPlanePoints(new TClonesArray("CbmRichPoint")) , fRichMirrorPoints(new TClonesArray("CbmRichPoint")) @@ -95,16 +88,7 @@ CbmRich::~CbmRich() { } } -void CbmRich::Initialize() { - FairDetector::Initialize(); - FairRun* fRun = FairRun::Instance(); - - fIsGeant4 = std::string(fRun->GetName()) == "TGeant4"; - LOG(info) << "CbmRich::Initialize() Transport engine:" - << std::string(fRun->GetName()); - LOG(info) << "CbmRich::Initialize() fIsGeant4:" - << (fIsGeant4 ? "true" : "false"); -} +void CbmRich::Initialize() { FairDetector::Initialize(); } Bool_t CbmRich::CheckIfSensitive(std::string name) { //return true; @@ -136,8 +120,7 @@ Bool_t CbmRich::ProcessHits(FairVolume* vol) { gMC->TrackPosition(tPos); gMC->TrackMomentum(tMom); - if (pdgCode == 50000050 - && isPhotonReflected(tMom.P())) { // Cherenkovs only + if (pdgCode == 50000050) { // Cherenkovs only AddHit(trackID, iVol, @@ -244,65 +227,6 @@ Bool_t CbmRich::ProcessHits(FairVolume* vol) { return kFALSE; } -Bool_t CbmRich::isPhotonReflected(Double_t photonEnergy) { - // for GEANT3 simulation no need to implement reflectivity - if (!fIsGeant4) return true; - - static std::vector<Double_t> enrg = { - 6.1997, 5.9045, 5.6361, 5.391, 5.1664, 4.9598, 4.769, 4.5924, 4.4284, - 4.2757, 4.1331, 3.9998, 3.8748, 3.7574, 3.6469, 3.5427, 3.4443, 3.3512, - 3.263, 3.1793, 3.0998, 3.0242, 2.9522, 2.8836, 2.818, 2.7554, 2.6955, - 2.6382, 2.5832, 2.5305, 2.4799, 2.4313, 2.3845, 2.3395, 2.2962, 2.2544, - 2.2142, 2.1753, 2.1378, 2.1016, 2.0666, 2.0327, 1.9999, 1.9682, 1.9374, - 1.9076, 1.8787, 1.8507, 1.8234, 1.797, 1.7713, 1.7464, 1.7221, 1.6985, - 1.6756, 1.6533, 1.6315, 1.6103, 1.5897, 1.5695, 1.5499}; - - - static std::vector<Double_t> eff = { - 0.22529, 0.1862, 0.15901, 0.14713, 0.13963, 0.13898, 0.13762, 0.13622, - 0.13868, 0.13951, 0.14031, 0.14073, 0.1409, 0.13977, 0.14205, 0.14072, - 0.1396, 0.13826, 0.13665, 0.13513, 0.13463, 0.13287, 0.13182, 0.13084, - 0.12824, 0.12601, 0.12622, 0.12681, 0.12193, 20.12011, 0.12109, 0.11908, - 0.11526, 0.11364, 0.11385, 0.12015, 0.11935, 0.11712, 0.1208, 0.12021, - 0.11909, 0.11783, 0.12257, 0.1215, 0.12199, 0.12494, 0.13101, 0.12089, - 0.12284, 0.12569, 0.13136, 0.13307, 0.13705, 0.13844, 0.13753, 0.14416, - 0.14761, 0.14953, 0.15218, 0.15315, 0.15719}; - - - // this interpolation causes run time crash - ///static ROOT::Math::Interpolator* interpolator = new ROOT::Math::Interpolator(reflEnergy, reflEff, ROOT::Math::Interpolation::kLINEAR); - // Double_t eff = interpolator->Eval(1.e9 * photonEnergy); - - Double_t phEnrgEv = 1.e9 * photonEnergy; - Double_t effCalc = 0.; - - if (phEnrgEv >= enrg[0]) { - effCalc = eff[0]; - } else if (phEnrgEv <= enrg[enrg.size() - 1]) { - effCalc = eff[eff.size() - 1]; - } else { - Int_t i0 = 0, i1 = 0; - for (size_t i = 0; i < enrg.size(); i++) { - if (phEnrgEv >= enrg[i]) { - i0 = i; - i1 = i0 - 1; - break; - } - } - if (i1 < 0) i1 = 0; - - effCalc = - eff[i0] - + (phEnrgEv - enrg[i0]) * (eff[i1] - eff[i0]) / (enrg[i1] - enrg[i0]); - // LOG(info) << "photonEn:" << phEnrgEv << " i0:" << i0 << ", " << enrg[i0] << ", " << eff[i0] << " i1:" << i1 << ", " << enrg[i1] << ", " << eff[i1] << - // " effCalc:" << effCalc << endl; - } - - Double_t rand = gRandom->Rndm(); - - return (rand > effCalc); -} - void CbmRich::EndOfEvent() { if (fVerboseLevel) Print(""); Reset(); @@ -357,6 +281,7 @@ void CbmRich::CopyClones(TClonesArray* cl1, TClonesArray* cl2, Int_t offset) { void CbmRich::ConstructOpGeometry() { LOG(info) << "CbmRich::ConstructOpGeometry()"; + SetRichGlassPropertiesForGeant4(); } void CbmRich::ConstructGeometry() { @@ -382,6 +307,68 @@ void CbmRich::ConstructGeometry() { } } +void CbmRich::SetRichGlassPropertiesForGeant4() { + FairRun* run = FairRun::Instance(); + Bool_t isGeant4 = std::string(run->GetName()) == "TGeant4"; + LOG(info) << "CbmRich::SetRichGlassPropertiesForGeant4() Transport engine:" + << std::string(run->GetName()); + + // for GEANT3 simulation no need to implement reflectivity + if (!isGeant4) { + LOG(info) << "CbmRich::SetRichGlassPropertiesForGeant4() fIsGeant4 is " + "false. No need to set RICH glass properties. Return."; + return; + } else { + LOG(info) << "CbmRich::SetRichGlassPropertiesForGeant4() fIsGeant4 is " + "true. RICH glass properties will be set."; + } + + std::vector<Double_t> energyAr = { + 6.1997, 5.9045, 5.6361, 5.391, 5.1664, 4.9598, 4.769, 4.5924, 4.4284, + 4.2757, 4.1331, 3.9998, 3.8748, 3.7574, 3.6469, 3.5427, 3.4443, 3.3512, + 3.263, 3.1793, 3.0998, 3.0242, 2.9522, 2.8836, 2.818, 2.7554, 2.6955, + 2.6382, 2.5832, 2.5305, 2.4799, 2.4313, 2.3845, 2.3395, 2.2962, 2.2544, + 2.2142, 2.1753, 2.1378, 2.1016, 2.0666, 2.0327, 1.9999, 1.9682, 1.9374, + 1.9076, 1.8787, 1.8507, 1.8234, 1.797, 1.7713, 1.7464, 1.7221, 1.6985, + 1.6756, 1.6533, 1.6315, 1.6103, 1.5897, 1.5695, 1.5499}; + // transform to GeV + for (size_t i = 0; i < energyAr.size(); i++) { + energyAr[i] = 1.e-9 * energyAr[i]; + } + + std::vector<Double_t> reflectivityAr = { + 0.22529, 0.1862, 0.15901, 0.14713, 0.13963, 0.13898, 0.13762, 0.13622, + 0.13868, 0.13951, 0.14031, 0.14073, 0.1409, 0.13977, 0.14205, 0.14072, + 0.1396, 0.13826, 0.13665, 0.13513, 0.13463, 0.13287, 0.13182, 0.13084, + 0.12824, 0.12601, 0.12622, 0.12681, 0.12193, 0.12011, 0.12109, 0.11908, + 0.11526, 0.11364, 0.11385, 0.12015, 0.11935, 0.11712, 0.1208, 0.12021, + 0.11909, 0.11783, 0.12257, 0.1215, 0.12199, 0.12494, 0.13101, 0.12089, + 0.12284, 0.12569, 0.13136, 0.13307, 0.13705, 0.13844, 0.13753, 0.14416, + 0.14761, 0.14953, 0.15218, 0.15315, 0.15719}; + + for (size_t i = 0; i < reflectivityAr.size(); i++) { + reflectivityAr[i] = 1. - reflectivityAr[i]; + } + + gMC->DefineOpSurface( + "RICHglassSurf", kGlisur, kDielectric_metal, kPolished, 0.0); + gMC->SetMaterialProperty("RICHglassSurf", + "REFLECTIVITY", + energyAr.size(), + energyAr.data(), + reflectivityAr.data()); + + for (int i = 0; i < 10; i++) { + if (gGeoManager->FindVolumeFast( + ("mirror_tile_type" + std::to_string(i)).c_str()) + == nullptr) + continue; + gMC->SetSkinSurface(("RICHglassSkin" + std::to_string(i)).c_str(), + ("mirror_tile_type" + std::to_string(i)).c_str(), + "RICHglassSurf"); + } +} + void CbmRich::ConstructGdmlGeometry(TGeoMatrix* geoMatrix) { TFile* old = gFile; TGDMLParse parser; diff --git a/sim/detectors/rich/CbmRich.h b/sim/detectors/rich/CbmRich.h index 8b610fba721020ffcfa87332d245e0653be98da0..11b95541c3c453d694d2a1cac1b8b17d4e1feb55 100644 --- a/sim/detectors/rich/CbmRich.h +++ b/sim/detectors/rich/CbmRich.h @@ -86,13 +86,6 @@ public: */ virtual Bool_t ProcessHits(FairVolume* vol = 0); - - /** - * \brief Photon reflection efficiency. - * Workaround implementation for GEANT4. - */ - Bool_t isPhotonReflected(Double_t photonEnergy); - /** * \brief If verbosity level is set, print hit collection at the * end of the event and resets it afterwards. @@ -158,6 +151,11 @@ public: */ void ConstructOpGeometry(); + /** + * \brief Set Cherenkov propeties for RICH mirror. + */ + void SetRichGlassPropertiesForGeant4(); + /** Check whether a volume is sensitive. ** The decision is based on the volume name. Only used in case @@ -181,9 +179,6 @@ private: // if false then only charged particles are registered Bool_t fRegisterPhotonsOnSensitivePlane; - // true if GEANT4 simulation is running - Bool_t fIsGeant4; - TClonesArray* fRichPoints; // MC points onto the photodetector plane TClonesArray* fRichRefPlanePoints; // points on the reference plane TClonesArray* fRichMirrorPoints; // mirror points diff --git a/sim/detectors/rich/CbmRichDigitizer.cxx b/sim/detectors/rich/CbmRichDigitizer.cxx index 79f458232d1c07ed37947a7113d439c4f4bb94dc..01414bd8b2c64fdb5647439b1aa8574f05362c44 100644 --- a/sim/detectors/rich/CbmRichDigitizer.cxx +++ b/sim/detectors/rich/CbmRichDigitizer.cxx @@ -47,7 +47,8 @@ CbmRichDigitizer::CbmRichDigitizer() , fTimeResolution(1.) , fDarkRatePerPixel(1000) , fPixelDeadTime(30.) - , fFiredPixelsMap() {} + , fFiredPixelsMap() + , fDoZShift(true) {} CbmRichDigitizer::~CbmRichDigitizer() {} @@ -181,7 +182,8 @@ void CbmRichDigitizer::ProcessPoint(CbmRichPoint* point, // LOG(info) << "ProcessPoint " << pointId; // TGeoNode* node = gGeoManager->FindNode(point->GetX(), point->GetY(), point->GetZ()); // workaround for GEANT4, probably boundary problem - Double_t zNew = point->GetZ() - 0.001; + Double_t zNew = point->GetZ(); + if (fDoZShift) zNew = zNew - 0.001; gGeoManager->FindNode(point->GetX(), point->GetY(), zNew); string path(gGeoManager->GetPath()); Int_t address = diff --git a/sim/detectors/rich/CbmRichDigitizer.h b/sim/detectors/rich/CbmRichDigitizer.h index 5a0bd46bf3632283ac8f28ee7a85cc5228c658d1..9720b1be171e72e999bdf7e96730b7c9a5b5a06c 100644 --- a/sim/detectors/rich/CbmRichDigitizer.h +++ b/sim/detectors/rich/CbmRichDigitizer.h @@ -115,6 +115,11 @@ public: fMaxNofHitsPerPmtCut = nofHits; } + /** + * \brief Set if you want to shift z MC point value (workaround for GEANT4). + */ + void SetDoZShift(Bool_t doZShift) { fDoZShift = doZShift; } + private: Int_t fEventNum; @@ -133,7 +138,7 @@ private: fCrossTalkProbability; // probability of the crosstalk for direct neighbor for one pixel Double_t fNoiseDigiRate; // noise rate per McRichPoint / per pixel / per second : - // hofNoiseDigis = nofRichPoints * nofPixels * dT(50 ns) * (fNoiseDigiRate / 1.e9); + // hofNoiseDigis = nofRichPoints * nofPixels * dT(50 ns) * (fNoiseDigiRate / 1.e9); CbmRichPmtTypeEnum fDetectorType; Int_t fMaxNofHitsPerPmtCut; // maximum number of hits which can be registered per PMT per event. @@ -146,6 +151,8 @@ private: Double_t fPixelDeadTime; // in ns, during this time pixel can not be fired map<Int_t, Double_t> fFiredPixelsMap; // first: pixel address, second: last fired time. + Bool_t + fDoZShift; // Set if you want to shift z MC point value (workaround for GEANT4). Must be set to true if one runs full RICH geoemtry with GEANT4, fot mCBM set to false /* * \brief Add crasstalk digis to the output array for the digi assuming fCrossTalkProbability