From cbb9b74da55c479dec319ccefc1f37203fbe26ed Mon Sep 17 00:00:00 2001 From: Norbert Herrmann <n.herrmann@physi.uni-heidelberg.de> Date: Sun, 6 Dec 2020 10:31:31 +0100 Subject: [PATCH] modifications to include mtof in mcbm common production --- core/detectors/tof/CbmTofCreateDigiPar.cxx | 8 +- macro/beamtime/mcbm2020/ana_trks.C | 13 +- macro/beamtime/mcbm2020/ana_trks_eval.C | 20 +- macro/beamtime/mcbm2020/ana_trksi.C | 17 +- .../mcbm2020/mcbm_build_and_reco_kronos.C | 209 ++++++- macro/beamtime/mcbm2020/mcbm_event_reco.C | 215 +++++++- macro/beamtime/mcbm2020/mcbm_reco.C | 3 +- macro/beamtime/mcbm2020/mtof_reco.C | 27 +- reco/detectors/tof/CbmTofCalibrator.cxx | 78 ++- reco/detectors/tof/CbmTofCalibrator.h | 4 +- reco/detectors/tof/CbmTofEventClusterizer.cxx | 522 +++++++++--------- reco/detectors/tof/CbmTofFindTracks.cxx | 15 +- reco/detectors/tof/CbmTofFindTracks.h | 15 +- 13 files changed, 827 insertions(+), 319 deletions(-) diff --git a/core/detectors/tof/CbmTofCreateDigiPar.cxx b/core/detectors/tof/CbmTofCreateDigiPar.cxx index 57d4186efc..2891b4d9ef 100644 --- a/core/detectors/tof/CbmTofCreateDigiPar.cxx +++ b/core/detectors/tof/CbmTofCreateDigiPar.cxx @@ -140,8 +140,12 @@ InitStatus CbmTofCreateDigiPar::Init() { nodemap.insert(std::pair<Int_t, TGeoNode*>(iAddr, tGeoNode)); LOG(debug) << Form( "Digipar for %d, addr 0x%08x: Node=%p, x %6.2f, y %6.2f, z %6.2f ", - iCell, iAddr, tGeoNode, - fChannelInfo->GetX(), fChannelInfo->GetY(), fChannelInfo->GetZ()); + iCell, + iAddr, + tGeoNode, + fChannelInfo->GetX(), + fChannelInfo->GetY(), + fChannelInfo->GetZ()); } fDigiPar->SetNodeMap(nodemap); diff --git a/macro/beamtime/mcbm2020/ana_trks.C b/macro/beamtime/mcbm2020/ana_trks.C index 091442c90d..1699b4f66f 100644 --- a/macro/beamtime/mcbm2020/ana_trks.C +++ b/macro/beamtime/mcbm2020/ana_trks.C @@ -192,13 +192,14 @@ void ana_trks(Int_t nEvents = 10000, tofFindTracks->SetCalOpt( iCalOpt); // 1 - update offsets, 2 - update walk, 0 - bypass tofFindTracks->SetCorMode(iGenCor); // valid options: 0,1,2,3,4,5,6, 10 - 19 + //tofFindTracks->SetTtTarg(0.065); // target value for Mar2020 triple stack -> betapeak ~ 0.95 + //tofFindTracks->SetTtTarg(0.0605); // target value for Mar2020 triple stack -> betapeak ~ 0.95 + //tofFindTracks->SetTtTarg(0.058); // target value for Mar2020 double stack + //tofFindTracks->SetTtTarg(0.055); // target value Nov2019 (triple stack run 831) + //tofFindTracks->SetTtTarg(0.048); // target value Nov2019 (double stack run 714) tofFindTracks->SetTtTarg( - 0.065); // target value for Mar2020 triple stack -> betapeak ~ 0.95 - //tofFindTracks->SetTtTarg(0.0605); // target value for Mar2020 triple stack -> betapeak ~ 0.95 - //tofFindTracks->SetTtTarg(0.058); // target value for Mar2020 double stack - //tofFindTracks->SetTtTarg(0.055); // target value Nov2019 (triple stack run 831) - //tofFindTracks->SetTtTarg(0.048); // target value Nov2019 (double stack run 714) - //tofFindTracks->SetTtTarg(0.035); // target value for inverse velocity, > 0.033 ns/cm! + 0.044); // target value Mar2020, after T0 fix (double stack run 714) + //tofFindTracks->SetTtTarg(0.035); // target value for inverse velocity, > 0.033 ns/cm! tofFindTracks->SetCalParFileName( cTrkFile); // Tracker parameter value file name tofFindTracks->SetBeamCounter(5, 0, 0); // default beam counter diff --git a/macro/beamtime/mcbm2020/ana_trks_eval.C b/macro/beamtime/mcbm2020/ana_trks_eval.C index 051803ef97..49b4384d15 100644 --- a/macro/beamtime/mcbm2020/ana_trks_eval.C +++ b/macro/beamtime/mcbm2020/ana_trks_eval.C @@ -689,11 +689,11 @@ void ana_trks_eval(Int_t nEvents = 10000, switch (iRSelin) { case 500: if (iMc == 0) { - tofAnaTestbeam->SetTShift(4.8); // Shift DTD4 to 0 - tofAnaTestbeam->SetTOffD4(11.); // Shift DTD4 to physical value - } else { // MC - tofAnaTestbeam->SetTShift(-12.); // Shift DTD4 to 0 - tofAnaTestbeam->SetTOffD4(15.); // Shift DTD4 to physical value + tofAnaTestbeam->SetTShift(4.8); // Shift DTD4 to 0 + tofAnaTestbeam->SetTOffD4(11.); // Shift DTD4 to physical value + } else { // MC + tofAnaTestbeam->SetTShift(-2.); // Shift DTD4 to 0 + tofAnaTestbeam->SetTOffD4(15.); // Shift DTD4 to physical value } switch (iSel2in) { case 30: @@ -702,15 +702,19 @@ void ana_trks_eval(Int_t nEvents = 10000, case 31: if (iMc == 0) { switch (iRun) { - case 717: + case 727: + case 726: + case 723: + case 721: tofAnaTestbeam->SetTShift(6.5); // Shift DTD4 to 0 tofAnaTestbeam->SetSel2TOff( 0.6); // Shift Sel2 time peak to 0 break; + case 717: default: // 714 - // tofAnaTestbeam->SetSel2TOff(-1.3); // Shift Sel2 time peak to 0 + //tofAnaTestbeam->SetSel2TOff(-1.3); // Shift Sel2 time peak to 0 tofAnaTestbeam->SetSel2TOff( - 0.3); // Shift Sel2 time peak to 0 + 0.); // Shift Sel2 time peak to 0 } } else { // MC tofAnaTestbeam->SetSel2TOff( diff --git a/macro/beamtime/mcbm2020/ana_trksi.C b/macro/beamtime/mcbm2020/ana_trksi.C index 07036e1cff..eca5d8b347 100644 --- a/macro/beamtime/mcbm2020/ana_trksi.C +++ b/macro/beamtime/mcbm2020/ana_trksi.C @@ -200,7 +200,7 @@ void ana_trksi(Int_t nEvents = 10000, tofFindTracks->SetCalParFileName( cTrkFile); // Tracker parameter value file name tofFindTracks->SetBeamCounter(5, 0, 0); // default beam counter - tofFindTracks->SetR0Lim(100.); + tofFindTracks->SetR0Lim(2.); tofFindTracks->SetStationMaxHMul( 30); // Max Hit Multiplicity in any used station @@ -214,6 +214,7 @@ void ana_trksi(Int_t nEvents = 10000, tofTrackFinder->SetSIGLIM(dChi2Lim2 * 2.); // matching window in multiples of chi2 tofTrackFinder->SetChiMaxAccept(dChi2Lim2); // max tracklet chi2 + tofTrackFinder->SetSIGLIMMOD(5.); // max deviation for last hit Int_t iMinNofHits = -1; Int_t iNStations = 0; @@ -468,7 +469,7 @@ void ana_trksi(Int_t nEvents = 10000, } //CbmTofAnaTestbeam defaults tofAnaTestbeam->SetR0LimFit( - 20.); // limit distance of fitted track to nominal vertex + 2.); // limit distance of fitted track to nominal vertex tofAnaTestbeam->SetDXMean(0.); tofAnaTestbeam->SetDYMean(0.); tofAnaTestbeam->SetDTMean(0.); // in ns @@ -593,22 +594,26 @@ void ana_trksi(Int_t nEvents = 10000, case 901041: switch (iRSelin) { case 500: - tofAnaTestbeam->SetTShift(-12.); // Shift DTD4 to 0 - tofAnaTestbeam->SetTOffD4(15.); // Shift DTD4 to physical value + tofAnaTestbeam->SetTShift(-2.); // Shift DTD4 to 0 + tofAnaTestbeam->SetTOffD4(15.); // Shift DTD4 to physical value switch (iSel2in) { case 30: tofAnaTestbeam->SetSel2TOff(-0.3); // Shift Sel2 time peak to 0 break; case 31: switch (iRun) { - case 717: + case 727: + case 726: + case 723: + case 721: tofAnaTestbeam->SetTShift(6.5); // Shift DTD4 to 0 tofAnaTestbeam->SetSel2TOff( 0.6); // Shift Sel2 time peak to 0 break; + case 717: default: //714 tofAnaTestbeam->SetSel2TOff( - 0.3); // Shift Sel2 time peak to 0 + 0.); // Shift Sel2 time peak to 0 } break; case 600: diff --git a/macro/beamtime/mcbm2020/mcbm_build_and_reco_kronos.C b/macro/beamtime/mcbm2020/mcbm_build_and_reco_kronos.C index f4d5108754..fd8becba10 100644 --- a/macro/beamtime/mcbm2020/mcbm_build_and_reco_kronos.C +++ b/macro/beamtime/mcbm2020/mcbm_build_and_reco_kronos.C @@ -39,7 +39,7 @@ void mcbm_build_and_reco_kronos(UInt_t uRunIdx = 28, inFile = Form("/lustre/cbm/users/ploizeau/mcbm2020/" "unp_evt_data_7f229b3f_20201103/unp_mcbm_%u.root", uRunId); - parFileIn = Form("/lustre/cbm/users/ploizeau/mcbm2020/ + parFileIn = Form("/lustre/cbm/users/ploizeau/mcbm2020/" "unp_evt_data_7f229b3f_20201103/unp_mcbm_params_%u.root", uRunId); } // if( 99999 != uRunIdx ) @@ -277,10 +277,214 @@ void mcbm_build_and_reco_kronos(UInt_t uRunIdx = 28, std::cout << "-I- : Added task " << trdHit->GetName() << std::endl; // ------------------------------------------------------------------------ - // ----- Local reconstruction in TOF ---------------------------------- // ------------------------------------------------------------------------ + // TOF defaults + Int_t calMode = 93; + Int_t calSel = 1; + Int_t calSm = 0; + Int_t RefSel = 0; + Double_t dDeadtime = 50.; + Int_t iSel2 = 500; + TString TofGeoTag = "v20f_mcbm"; + TString cCalId = "831.50.3.0"; + Int_t iCalSet = 12022500; // calibration settings + + TObjString* tofBdfFile = + new TObjString(parDir + "/tof/tof_" + TofGeoTag + ".digibdf.par"); + parFileList->Add(tofBdfFile); + std::cout << "-I- Using parameter file " << tofBdfFile->GetString() + << std::endl; + CbmTofEventClusterizer* tofCluster = + new CbmTofEventClusterizer("TOF Event Clusterizer", 0, 1); + TString cFname = parDir + "/tof/" + + Form("%s_set%09d_%02d_%01dtofClust.hst.root", + cCalId.Data(), + iCalSet, + calMode, + calSel); + tofCluster->SetCalParFileName(cFname); + tofCluster->SetCalMode(calMode); + tofCluster->SetCalSel(calSel); + tofCluster->PosYMaxScal(0.75); //in % of 2*length + tofCluster->SetChannelDeadtime(dDeadtime); // artificial deadtime in ns + + run->AddTask(tofCluster); + std::cout << "-I- Added task " << tofCluster->GetName() << std::endl; + + // ----- Track reconstruction ------------------------------------------ + Int_t iTrackMode = 2; + switch (iTrackMode) { + case 2: { + Int_t iGenCor = 1; + Double_t dScalFac = 1.; + Double_t dChi2Lim2 = 3.5; + TString cTrkFile = + parDir + "/tof/" + Form("%s_tofFindTracks.hst.root", cCalId.Data()); + Int_t iTrackingSetup = 1; + Int_t iCalOpt = 0; + + CbmTofTrackFinder* tofTrackFinder = new CbmTofTrackFinderNN(); + tofTrackFinder->SetMaxTofTimeDifference(0.2); // in ns/cm + tofTrackFinder->SetTxLIM(0.3); // max slope dx/dz + tofTrackFinder->SetTyLIM(0.3); // max dev from mean slope dy/dz + tofTrackFinder->SetTyMean(0.); // mean slope dy/dz + CbmTofTrackFitter* tofTrackFitter = new CbmTofTrackFitterKF(0, 211); + TFitter* MyFit = new TFitter(1); // initialize Minuit + tofTrackFinder->SetFitter(tofTrackFitter); + + CbmTofFindTracks* tofFindTracks = + new CbmTofFindTracks("TOF Track Finder"); + tofFindTracks->UseFinder(tofTrackFinder); + tofFindTracks->UseFitter(tofTrackFitter); + tofFindTracks->SetCalOpt(iCalOpt); + // 1 - update offsets, 2 - update walk, 0 - bypass + tofFindTracks->SetCorMode( + iGenCor); // valid options: 0,1,2,3,4,5,6, 10 - 19 + tofFindTracks->SetTtTarg( + 0.065); // target value for Mar2020 triple stack -> betapeak ~ 0.95 + //tofFindTracks->SetTtTarg(0.041); // target value for inverse velocity, > 0.033 ns/cm! + //tofFindTracks->SetTtTarg(0.035); // target value for inverse velocity, > 0.033 ns/cm! + tofFindTracks->SetCalParFileName( + cTrkFile); // Tracker parameter value file name + tofFindTracks->SetBeamCounter(5, 0, 0); // default beam counter + tofFindTracks->SetStationMaxHMul( + 30); // Max Hit Multiplicity in any used station + + tofFindTracks->SetT0MAX(dScalFac); // in ns + tofFindTracks->SetSIGT(0.08); // default in ns + tofFindTracks->SetSIGX(0.3); // default in cm + tofFindTracks->SetSIGY(0.45); // default in cm + tofFindTracks->SetSIGZ(0.05); // default in cm + tofFindTracks->SetUseSigCalib( + kFALSE); // ignore resolutions in CalPar file + tofTrackFinder->SetSIGLIM(dChi2Lim2 + * 2.); // matching window in multiples of chi2 + tofTrackFinder->SetChiMaxAccept(dChi2Lim2); // max tracklet chi2 + + Int_t iMinNofHits = -1; + Int_t iNStations = 0; + Int_t iNReqStations = 3; + switch (iTrackingSetup) { + case 0: // bypass mode + iMinNofHits = -1; + iNStations = 1; + tofFindTracks->SetStation(0, 5, 0, 0); // Diamond + break; + + case 1: // for calibration mode of full setup + iMinNofHits = 3; + iNStations = 28; + iNReqStations = 4; + tofFindTracks->SetStation(0, 5, 0, 0); + tofFindTracks->SetStation(1, 0, 2, 2); + tofFindTracks->SetStation(2, 0, 1, 2); + tofFindTracks->SetStation(3, 0, 0, 2); + tofFindTracks->SetStation(4, 0, 2, 1); + tofFindTracks->SetStation(5, 0, 1, 1); + tofFindTracks->SetStation(6, 0, 0, 1); + tofFindTracks->SetStation(7, 0, 2, 3); + tofFindTracks->SetStation(8, 0, 1, 3); + tofFindTracks->SetStation(9, 0, 0, 3); + tofFindTracks->SetStation(10, 0, 2, 0); + tofFindTracks->SetStation(11, 0, 1, 0); + tofFindTracks->SetStation(12, 0, 0, 0); + tofFindTracks->SetStation(13, 0, 2, 4); + tofFindTracks->SetStation(14, 0, 1, 4); + tofFindTracks->SetStation(15, 0, 0, 4); + tofFindTracks->SetStation(16, 0, 4, 0); + tofFindTracks->SetStation(17, 0, 3, 0); + tofFindTracks->SetStation(18, 0, 4, 1); + tofFindTracks->SetStation(19, 0, 3, 1); + tofFindTracks->SetStation(20, 0, 4, 2); + tofFindTracks->SetStation(21, 0, 3, 2); + tofFindTracks->SetStation(22, 0, 4, 3); + tofFindTracks->SetStation(23, 0, 3, 3); + tofFindTracks->SetStation(24, 0, 4, 4); + tofFindTracks->SetStation(25, 0, 3, 4); + tofFindTracks->SetStation(26, 9, 0, 0); + tofFindTracks->SetStation(27, 9, 0, 1); + break; + + case 2: // for geometry check mode of full setup + iMinNofHits = 3; + iNStations = 27; + iNReqStations = 4; + tofFindTracks->SetStation(0, 0, 2, 2); + tofFindTracks->SetStation(1, 0, 1, 2); + tofFindTracks->SetStation(2, 0, 0, 2); + tofFindTracks->SetStation(3, 0, 2, 1); + tofFindTracks->SetStation(4, 0, 1, 1); + tofFindTracks->SetStation(5, 0, 0, 1); + tofFindTracks->SetStation(6, 0, 2, 3); + tofFindTracks->SetStation(7, 0, 1, 3); + tofFindTracks->SetStation(8, 0, 0, 3); + tofFindTracks->SetStation(9, 0, 2, 0); + tofFindTracks->SetStation(10, 0, 1, 0); + tofFindTracks->SetStation(11, 0, 0, 0); + tofFindTracks->SetStation(12, 0, 2, 4); + tofFindTracks->SetStation(13, 0, 1, 4); + tofFindTracks->SetStation(14, 0, 0, 4); + tofFindTracks->SetStation(15, 0, 4, 0); + tofFindTracks->SetStation(16, 0, 3, 0); + tofFindTracks->SetStation(17, 0, 4, 1); + tofFindTracks->SetStation(18, 0, 3, 1); + tofFindTracks->SetStation(19, 0, 4, 2); + tofFindTracks->SetStation(20, 0, 3, 2); + tofFindTracks->SetStation(21, 0, 4, 3); + tofFindTracks->SetStation(22, 0, 3, 3); + tofFindTracks->SetStation(23, 0, 4, 4); + tofFindTracks->SetStation(24, 0, 3, 4); + tofFindTracks->SetStation(25, 9, 0, 0); + tofFindTracks->SetStation(26, 9, 0, 1); + break; + + case 3: // for reduced bias tracking of full setup + iMinNofHits = 3; + iNStations = 28; + iNReqStations = 4; + tofFindTracks->SetStation(0, 0, 2, 2); + tofFindTracks->SetStation(1, 0, 1, 2); + tofFindTracks->SetStation(2, 0, 0, 2); + tofFindTracks->SetStation(3, 0, 2, 1); + tofFindTracks->SetStation(4, 0, 1, 1); + tofFindTracks->SetStation(5, 0, 0, 1); + tofFindTracks->SetStation(6, 0, 2, 3); + tofFindTracks->SetStation(7, 0, 1, 3); + tofFindTracks->SetStation(8, 0, 0, 3); + tofFindTracks->SetStation(9, 0, 2, 0); + tofFindTracks->SetStation(10, 0, 1, 0); + tofFindTracks->SetStation(11, 0, 0, 0); + tofFindTracks->SetStation(12, 0, 2, 4); + tofFindTracks->SetStation(13, 0, 1, 4); + tofFindTracks->SetStation(14, 0, 0, 4); + tofFindTracks->SetStation(15, 0, 4, 0); + tofFindTracks->SetStation(16, 0, 3, 0); + tofFindTracks->SetStation(17, 0, 4, 1); + tofFindTracks->SetStation(18, 0, 3, 1); + tofFindTracks->SetStation(19, 0, 4, 2); + tofFindTracks->SetStation(20, 0, 3, 2); + tofFindTracks->SetStation(21, 0, 4, 3); + tofFindTracks->SetStation(22, 0, 3, 3); + tofFindTracks->SetStation(23, 0, 4, 4); + tofFindTracks->SetStation(24, 0, 3, 4); + tofFindTracks->SetStation(25, 9, 0, 0); + tofFindTracks->SetStation(26, 9, 0, 1); + tofFindTracks->SetStation(27, 5, 0, 0); + break; + } + tofFindTracks->SetMinNofHits(iMinNofHits); + tofFindTracks->SetNStations(iNStations); + tofFindTracks->SetNReqStations(iNReqStations); + //tofFindTracks->PrintSetup(); + run->AddTask(tofFindTracks); + } break; + case 1: { + } + case 0: + default:; + } // ----- Local reconstruction of RICH Hits ------------------------------ CbmRichMCbmHitProducer* hitProdRich = new CbmRichMCbmHitProducer(); @@ -317,6 +521,7 @@ void mcbm_build_and_reco_kronos(UInt_t uRunIdx = 28, parIo2->open(parFileList, "in"); rtdb->setSecondInput(parIo2); parIo3->open(parFileOut.Data(), "RECREATE"); + // ------------------------------------------------------------------------ diff --git a/macro/beamtime/mcbm2020/mcbm_event_reco.C b/macro/beamtime/mcbm2020/mcbm_event_reco.C index ff8aa4d97f..96019068ca 100644 --- a/macro/beamtime/mcbm2020/mcbm_event_reco.C +++ b/macro/beamtime/mcbm2020/mcbm_event_reco.C @@ -6,7 +6,7 @@ // // -------------------------------------------------------------------------- -void mcbm_event_reco(Int_t runId = 831, Int_t nTimeslices = 300) { +void mcbm_event_reco(Int_t runId = 831, Int_t nTimeslices = -1) { // --- Logger settings ---------------------------------------------------- TString logLevel = "INFO"; @@ -15,9 +15,10 @@ void mcbm_event_reco(Int_t runId = 831, Int_t nTimeslices = 300) { // ----- Environment -------------------------------------------------- - TString myName = "mcbm_reco"; // this macro's name for screen output + TString myName = "mcbm_event_reco"; // this macro's name for screen output TString srcDir = gSystem->Getenv("VMCWORKDIR"); // top source directory TString paramDir = srcDir + "/macro/beamtime/mcbm2020/"; + TString parDir = srcDir + "/parameters"; // ------------------------------------------------------------------------ @@ -26,7 +27,6 @@ void mcbm_event_reco(Int_t runId = 831, Int_t nTimeslices = 300) { TString parFile = Form("./data/unp_mcbm_params_%i.root", runId); TString outFile = Form("./data/reco_mcbm_%i.root", runId); // ------------------------------------------------------------------------ - // --- Load the geometry setup ---- // This is currently only required by the TRD (parameters) std::string geoSetupTag = "mcbm_beam_2020_03"; @@ -37,6 +37,7 @@ void mcbm_event_reco(Int_t runId = 831, Int_t nTimeslices = 300) { TList* parFileList = new TList(); // ------------------------------------------------------------------------ + // ----- Timer -------------------------------------------------------- TStopwatch timer; timer.Start(); @@ -84,8 +85,6 @@ void mcbm_event_reco(Int_t runId = 831, Int_t nTimeslices = 300) { // ----- Local reconstruction in MUCH --------------------------------- Int_t flag = 1; - TString parDir = - TString(gSystem->Getenv("VMCWORKDIR")) + TString("/parameters"); TString muchDigiFile( parDir + "/much/much_v19c_mcbm_digi_sector.root"); // MUCH digi file CbmMuchFindHitsGem* muchFindHits = @@ -177,7 +176,212 @@ void mcbm_event_reco(Int_t runId = 831, Int_t nTimeslices = 300) { // ----- Local reconstruction in TOF ---------------------------------- // ------------------------------------------------------------------------ + // TOF defaults + Int_t calMode = 93; + Int_t calSel = 1; + Int_t calSm = 0; + Int_t RefSel = 0; + Double_t dDeadtime = 50.; + Int_t iSel2 = 500; + TString TofGeoTag = "v20f_mcbm"; + TString cCalId = "831.50.3.0"; + Int_t iCalSet = 12022500; // calibration settings + + TObjString* tofBdfFile = + new TObjString(parDir + "/tof/tof_" + TofGeoTag + ".digibdf.par"); + parFileList->Add(tofBdfFile); + std::cout << "-I- Using parameter file " << tofBdfFile->GetString() + << std::endl; + CbmTofEventClusterizer* tofCluster = + new CbmTofEventClusterizer("TOF Event Clusterizer", 0, 1); + TString cFname = parDir + "/tof/" + + Form("%s_set%09d_%02d_%01dtofClust.hst.root", + cCalId.Data(), + iCalSet, + calMode, + calSel); + tofCluster->SetCalParFileName(cFname); + tofCluster->SetCalMode(calMode); + tofCluster->SetCalSel(calSel); + tofCluster->PosYMaxScal(0.75); //in % of 2*length + tofCluster->SetChannelDeadtime(dDeadtime); // artificial deadtime in ns + + run->AddTask(tofCluster); + std::cout << "-I- Added task " << tofCluster->GetName() << std::endl; + + // ----- Track reconstruction ------------------------------------------ + Int_t iTrackMode = 2; + switch (iTrackMode) { + case 2: { + Int_t iGenCor = 1; + Double_t dScalFac = 1.; + Double_t dChi2Lim2 = 3.5; + TString cTrkFile = + parDir + "/tof/" + Form("%s_tofFindTracks.hst.root", cCalId.Data()); + Int_t iTrackingSetup = 1; + Int_t iCalOpt = 0; + + CbmTofTrackFinder* tofTrackFinder = new CbmTofTrackFinderNN(); + tofTrackFinder->SetMaxTofTimeDifference(0.2); // in ns/cm + tofTrackFinder->SetTxLIM(0.3); // max slope dx/dz + tofTrackFinder->SetTyLIM(0.3); // max dev from mean slope dy/dz + tofTrackFinder->SetTyMean(0.); // mean slope dy/dz + CbmTofTrackFitter* tofTrackFitter = new CbmTofTrackFitterKF(0, 211); + TFitter* MyFit = new TFitter(1); // initialize Minuit + tofTrackFinder->SetFitter(tofTrackFitter); + + CbmTofFindTracks* tofFindTracks = + new CbmTofFindTracks("TOF Track Finder"); + tofFindTracks->UseFinder(tofTrackFinder); + tofFindTracks->UseFitter(tofTrackFitter); + tofFindTracks->SetCalOpt(iCalOpt); + // 1 - update offsets, 2 - update walk, 0 - bypass + tofFindTracks->SetCorMode( + iGenCor); // valid options: 0,1,2,3,4,5,6, 10 - 19 + tofFindTracks->SetTtTarg( + 0.065); // target value for Mar2020 triple stack -> betapeak ~ 0.95 + //tofFindTracks->SetTtTarg(0.041); // target value for inverse velocity, > 0.033 ns/cm! + //tofFindTracks->SetTtTarg(0.035); // target value for inverse velocity, > 0.033 ns/cm! + tofFindTracks->SetCalParFileName( + cTrkFile); // Tracker parameter value file name + tofFindTracks->SetBeamCounter(5, 0, 0); // default beam counter + tofFindTracks->SetStationMaxHMul( + 30); // Max Hit Multiplicity in any used station + + tofFindTracks->SetT0MAX(dScalFac); // in ns + tofFindTracks->SetSIGT(0.08); // default in ns + tofFindTracks->SetSIGX(0.3); // default in cm + tofFindTracks->SetSIGY(0.45); // default in cm + tofFindTracks->SetSIGZ(0.05); // default in cm + tofFindTracks->SetUseSigCalib( + kFALSE); // ignore resolutions in CalPar file + tofTrackFinder->SetSIGLIM(dChi2Lim2 + * 2.); // matching window in multiples of chi2 + tofTrackFinder->SetChiMaxAccept(dChi2Lim2); // max tracklet chi2 + + Int_t iMinNofHits = -1; + Int_t iNStations = 0; + Int_t iNReqStations = 3; + switch (iTrackingSetup) { + case 0: // bypass mode + iMinNofHits = -1; + iNStations = 1; + tofFindTracks->SetStation(0, 5, 0, 0); // Diamond + break; + + case 1: // for calibration mode of full setup + iMinNofHits = 3; + iNStations = 28; + iNReqStations = 4; + tofFindTracks->SetStation(0, 5, 0, 0); + tofFindTracks->SetStation(1, 0, 2, 2); + tofFindTracks->SetStation(2, 0, 1, 2); + tofFindTracks->SetStation(3, 0, 0, 2); + tofFindTracks->SetStation(4, 0, 2, 1); + tofFindTracks->SetStation(5, 0, 1, 1); + tofFindTracks->SetStation(6, 0, 0, 1); + tofFindTracks->SetStation(7, 0, 2, 3); + tofFindTracks->SetStation(8, 0, 1, 3); + tofFindTracks->SetStation(9, 0, 0, 3); + tofFindTracks->SetStation(10, 0, 2, 0); + tofFindTracks->SetStation(11, 0, 1, 0); + tofFindTracks->SetStation(12, 0, 0, 0); + tofFindTracks->SetStation(13, 0, 2, 4); + tofFindTracks->SetStation(14, 0, 1, 4); + tofFindTracks->SetStation(15, 0, 0, 4); + tofFindTracks->SetStation(16, 0, 4, 0); + tofFindTracks->SetStation(17, 0, 3, 0); + tofFindTracks->SetStation(18, 0, 4, 1); + tofFindTracks->SetStation(19, 0, 3, 1); + tofFindTracks->SetStation(20, 0, 4, 2); + tofFindTracks->SetStation(21, 0, 3, 2); + tofFindTracks->SetStation(22, 0, 4, 3); + tofFindTracks->SetStation(23, 0, 3, 3); + tofFindTracks->SetStation(24, 0, 4, 4); + tofFindTracks->SetStation(25, 0, 3, 4); + tofFindTracks->SetStation(26, 9, 0, 0); + tofFindTracks->SetStation(27, 9, 0, 1); + break; + + case 2: // for geometry check mode of full setup + iMinNofHits = 3; + iNStations = 27; + iNReqStations = 4; + tofFindTracks->SetStation(0, 0, 2, 2); + tofFindTracks->SetStation(1, 0, 1, 2); + tofFindTracks->SetStation(2, 0, 0, 2); + tofFindTracks->SetStation(3, 0, 2, 1); + tofFindTracks->SetStation(4, 0, 1, 1); + tofFindTracks->SetStation(5, 0, 0, 1); + tofFindTracks->SetStation(6, 0, 2, 3); + tofFindTracks->SetStation(7, 0, 1, 3); + tofFindTracks->SetStation(8, 0, 0, 3); + tofFindTracks->SetStation(9, 0, 2, 0); + tofFindTracks->SetStation(10, 0, 1, 0); + tofFindTracks->SetStation(11, 0, 0, 0); + tofFindTracks->SetStation(12, 0, 2, 4); + tofFindTracks->SetStation(13, 0, 1, 4); + tofFindTracks->SetStation(14, 0, 0, 4); + tofFindTracks->SetStation(15, 0, 4, 0); + tofFindTracks->SetStation(16, 0, 3, 0); + tofFindTracks->SetStation(17, 0, 4, 1); + tofFindTracks->SetStation(18, 0, 3, 1); + tofFindTracks->SetStation(19, 0, 4, 2); + tofFindTracks->SetStation(20, 0, 3, 2); + tofFindTracks->SetStation(21, 0, 4, 3); + tofFindTracks->SetStation(22, 0, 3, 3); + tofFindTracks->SetStation(23, 0, 4, 4); + tofFindTracks->SetStation(24, 0, 3, 4); + tofFindTracks->SetStation(25, 9, 0, 0); + tofFindTracks->SetStation(26, 9, 0, 1); + break; + + case 3: // for reduced bias tracking of full setup + iMinNofHits = 3; + iNStations = 28; + iNReqStations = 4; + tofFindTracks->SetStation(0, 0, 2, 2); + tofFindTracks->SetStation(1, 0, 1, 2); + tofFindTracks->SetStation(2, 0, 0, 2); + tofFindTracks->SetStation(3, 0, 2, 1); + tofFindTracks->SetStation(4, 0, 1, 1); + tofFindTracks->SetStation(5, 0, 0, 1); + tofFindTracks->SetStation(6, 0, 2, 3); + tofFindTracks->SetStation(7, 0, 1, 3); + tofFindTracks->SetStation(8, 0, 0, 3); + tofFindTracks->SetStation(9, 0, 2, 0); + tofFindTracks->SetStation(10, 0, 1, 0); + tofFindTracks->SetStation(11, 0, 0, 0); + tofFindTracks->SetStation(12, 0, 2, 4); + tofFindTracks->SetStation(13, 0, 1, 4); + tofFindTracks->SetStation(14, 0, 0, 4); + tofFindTracks->SetStation(15, 0, 4, 0); + tofFindTracks->SetStation(16, 0, 3, 0); + tofFindTracks->SetStation(17, 0, 4, 1); + tofFindTracks->SetStation(18, 0, 3, 1); + tofFindTracks->SetStation(19, 0, 4, 2); + tofFindTracks->SetStation(20, 0, 3, 2); + tofFindTracks->SetStation(21, 0, 4, 3); + tofFindTracks->SetStation(22, 0, 3, 3); + tofFindTracks->SetStation(23, 0, 4, 4); + tofFindTracks->SetStation(24, 0, 3, 4); + tofFindTracks->SetStation(25, 9, 0, 0); + tofFindTracks->SetStation(26, 9, 0, 1); + tofFindTracks->SetStation(27, 5, 0, 0); + break; + } + tofFindTracks->SetMinNofHits(iMinNofHits); + tofFindTracks->SetNStations(iNStations); + tofFindTracks->SetNReqStations(iNReqStations); + //tofFindTracks->PrintSetup(); + run->AddTask(tofFindTracks); + } break; + case 1: { + } + case 0: + default:; + } // ----- Local reconstruction of RICH Hits ------------------------------ CbmRichMCbmHitProducer* hitProdRich = new CbmRichMCbmHitProducer(); @@ -209,6 +413,7 @@ void mcbm_event_reco(Int_t runId = 831, Int_t nTimeslices = 300) { FairParRootFileIo* parIo1 = new FairParRootFileIo(); FairParAsciiFileIo* parIo2 = new FairParAsciiFileIo(); parIo1->open(parFile.Data(), "UPDATE"); + parIo2->open(parFileList, "in"); rtdb->setFirstInput(parIo1); parIo2->open(parFileList, "in"); rtdb->setSecondInput(parIo2); diff --git a/macro/beamtime/mcbm2020/mcbm_reco.C b/macro/beamtime/mcbm2020/mcbm_reco.C index e17246e511..4087e633be 100644 --- a/macro/beamtime/mcbm2020/mcbm_reco.C +++ b/macro/beamtime/mcbm2020/mcbm_reco.C @@ -192,6 +192,7 @@ void mcbm_reco(Int_t runId = 831, Int_t nTimeslices = 0) { FairParRootFileIo* parIo1 = new FairParRootFileIo(); FairParAsciiFileIo* parIo2 = new FairParAsciiFileIo(); parIo1->open(parFile.Data(), "UPDATE"); + parIo2->open(parFileList, "in"); rtdb->setFirstInput(parIo1); parIo2->open(parFileList, "in"); rtdb->setSecondInput(parIo2); @@ -223,7 +224,7 @@ void mcbm_reco(Int_t runId = 831, Int_t nTimeslices = 0) { std::cout << std::endl << std::endl; std::cout << "Macro finished successfully." << std::endl; std::cout << "Output file is " << outFile << std::endl; - std::cout << "Parameter file is " << parFile << std::endl; + std::cout << "Parameter file is " << parFileOut << std::endl; std::cout << "Real time " << rtime << " s, CPU time " << ctime << " s" << std::endl; std::cout << std::endl; diff --git a/macro/beamtime/mcbm2020/mtof_reco.C b/macro/beamtime/mcbm2020/mtof_reco.C index 7a9fa316dd..a130313784 100644 --- a/macro/beamtime/mcbm2020/mtof_reco.C +++ b/macro/beamtime/mcbm2020/mtof_reco.C @@ -3,7 +3,7 @@ // // -------------------------------------------------------------------------- -void mtof_reco(Int_t nEvents = 100, // number of Timeslices +void mtof_reco(Int_t nEvents = -1, // number of Timeslices TString dataset = "data/unp_mcbm_831", TString setup = "mcbm_beam_2020_03", TString cCalId = "831.50.3.0", @@ -25,6 +25,7 @@ void mtof_reco(Int_t nEvents = 100, // number of Timeslices // ----- Environment -------------------------------------------------- TString myName = "mtof_reco"; // this macro's name for screen output TString srcDir = gSystem->Getenv("VMCWORKDIR"); // top source directory + TString parDir = srcDir + "/parameters"; // ------------------------------------------------------------------------ @@ -152,11 +153,12 @@ void mtof_reco(Int_t nEvents = 100, // number of Timeslices case 1: { CbmTofEventClusterizer* tofCluster = new CbmTofEventClusterizer("TOF Event Clusterizer", 0, 1); - TString cFname = Form("%s_set%09d_%02d_%01dtofClust.hst.root", - cCalId.Data(), - iCalSet, - calMode, - calSel); + TString cFname = parDir + "/tof/" + + Form("%s_set%09d_%02d_%01dtofClust.hst.root", + cCalId.Data(), + iCalSet, + calMode, + calSel); tofCluster->SetCalParFileName(cFname); tofCluster->SetCalMode(calMode); tofCluster->SetCalSel(calSel); @@ -223,7 +225,7 @@ void mtof_reco(Int_t nEvents = 100, // number of Timeslices Int_t iDutSm = iDut % 10; iDut = (iDut - iDutSm) / 10; - tofCluster->SetDutId(iDut); + //tofCluster->SetDutId(iDut); tofCluster->SetDutSm(iDutSm); tofCluster->SetDutRpc(iDutRpc); @@ -253,12 +255,13 @@ void mtof_reco(Int_t nEvents = 100, // number of Timeslices Double_t beamWidthY = 0.1; switch (iTrackMode) { case 2: { - Int_t iGenCor = 1; - Double_t dScalFac = 1.; - Double_t dChi2Lim2 = 3.5; - TString cTrkFile = Form("%s_tofFindTracks.hst.root", cCalId.Data()); + Int_t iGenCor = 1; + Double_t dScalFac = 1.; + Double_t dChi2Lim2 = 3.5; + TString cTrkFile = + parDir + "/tof/" + Form("%s_tofFindTracks.hst.root", cCalId.Data()); Int_t iTrackingSetup = 1; - Int_t iCalOpt = 0; + Int_t iCalOpt = 1; CbmTofTrackFinder* tofTrackFinder = new CbmTofTrackFinderNN(); tofTrackFinder->SetMaxTofTimeDifference(0.2); // in ns/cm diff --git a/reco/detectors/tof/CbmTofCalibrator.cxx b/reco/detectors/tof/CbmTofCalibrator.cxx index 7032e00e5e..9ff57c3cca 100644 --- a/reco/detectors/tof/CbmTofCalibrator.cxx +++ b/reco/detectors/tof/CbmTofCalibrator.cxx @@ -6,6 +6,7 @@ // CBMroot classes and includes #include "CbmTofCalibrator.h" +#include "CbmEvent.h" #include "CbmMatch.h" #include "CbmTofAddress.h" // in cbmdata/tof #include "CbmTofCell.h" // in tof/TofData @@ -62,6 +63,12 @@ CbmTofCalibrator::~CbmTofCalibrator() {} InitStatus CbmTofCalibrator::Init() { FairRootManager* fManager = FairRootManager::Instance(); + // Get access to TofCalDigis + fTofCalDigiVec = + fManager->InitObjectAs<std::vector<CbmTofDigi> const*>("TofCalDigi"); + //dynamic_cast<std::vector<CbmTofDigi> const*>(fManager->GetObject("TofCalDigi")); + if (NULL == fTofCalDigiVec) LOG(fatal) << "No access to TofCalDigis!"; + // check for availability of digis fDigiMan = CbmDigiManager::Instance(); if (NULL == fDigiMan) { @@ -88,9 +95,9 @@ InitStatus CbmTofCalibrator::Init() { } // Get Access to MatchCollection - fTofDigiMatchColl = (TClonesArray*) fManager->GetObject("TofDigiMatch"); + fTofDigiMatchColl = (TClonesArray*) fManager->GetObject("TofCalDigiMatch"); if (NULL == fTofDigiMatchColl) - fTofDigiMatchColl = (TClonesArray*) fManager->GetObject("TofCalDigiMatch"); + fTofDigiMatchColl = (TClonesArray*) fManager->GetObject("TofDigiMatch"); if (NULL == fTofDigiMatchColl) { LOG(error) << "CbmTofCalibrator: no access to DigiMatch "; @@ -243,7 +250,9 @@ Bool_t CbmTofCalibrator::CreateCalHist() { return kTRUE; } -void CbmTofCalibrator::FillCalHist(CbmTofTracklet* pTrk, Int_t iOpt) { +void CbmTofCalibrator::FillCalHist(CbmTofTracklet* pTrk, + Int_t iOpt, + CbmEvent* tEvent) { // fill deviation histograms on walk level if (pTrk->GetTt() < 0) return; // take tracks with positive velocity only if (fbBeam @@ -304,12 +313,14 @@ void CbmTofCalibrator::FillCalHist(CbmTofTracklet* pTrk, Int_t iOpt) { - pTrk->GetFitT(pHit->GetZ())); // residuals transformed into LRF //fhCalTOff[iDetIndx]->Fill((Double_t)iCh,fTrackletTools->GetTdif(pTrk, iDetId, pHit)); // prediction by other hits - Int_t iMA = pTrk->GetTofHitIndex(iHit); - if (iMA > fTofDigiMatchColl->GetEntries()) { - LOG(error) << " Inconsistent DigiMatches for Hitind " << iMA + Int_t iEA = pTrk->GetTofHitIndex(iHit); + Int_t iTSA = fTofFindTracks->GetTofHitIndex(iEA); + + if (iTSA > fTofDigiMatchColl->GetEntries()) { + LOG(error) << " Inconsistent DigiMatches for Hitind " << iTSA << ", TClonesArraySize: " << fTofDigiMatchColl->GetEntries(); } - CbmMatch* digiMatch = (CbmMatch*) fTofDigiMatchColl->At(iMA); + CbmMatch* digiMatch = (CbmMatch*) fTofDigiMatchColl->At(iTSA); Double_t hlocal_d[3]; for (Int_t iLink = 0; iLink < digiMatch->GetNofLinks(); @@ -318,11 +329,43 @@ void CbmTofCalibrator::FillCalHist(CbmTofTracklet* pTrk, Int_t iOpt) { Int_t iDigInd0 = L0.GetIndex(); Int_t iDigInd1 = (digiMatch->GetLink(iLink + 1)).GetIndex(); - const CbmTofDigi* tDigi0 = fDigiMan->Get<CbmTofDigi>(iDigInd0); - Int_t iCh0 = tDigi0->GetChannel(); - Int_t iSide0 = tDigi0->GetSide(); - LOG(debug) << "Fill Walk for " << iDetIndx << ", TSRCS " << iSmType << iSm - << iRpc << iCh0 << iSide0 << ", " << tDigi0 << ", " << pTrk; + const CbmTofDigi* tDigi0 = NULL; + const CbmTofDigi* tDigi1 = NULL; + if (tEvent != NULL) { //disable + LOG(debug) << "Locate MatchDigiInd " << iDigInd0 << " and " << iDigInd1 + << " in CalDigiVec of size " << fTofCalDigiVec->size(); + // <<" in current event not implemented"; + //continue; + tDigi0 = &(fTofCalDigiVec->at(iDigInd0)); + tDigi1 = &(fTofCalDigiVec->at(iDigInd1)); + } else { // event wise entries + tDigi0 = fDigiMan->Get<CbmTofDigi>(iDigInd0); + tDigi1 = fDigiMan->Get<CbmTofDigi>(iDigInd1); + } + + Int_t iCh0 = tDigi0->GetChannel(); + Int_t iSide0 = tDigi0->GetSide(); + + LOG(debug) << "Fill Walk for Hit Ind " << iEA << ", " << iTSA + << Form(", TSRC %d%d%d%2d, DigiInd %2d, %2d", + iSmType, + iSm, + iRpc, + iCh, + iDigInd0, + iDigInd1) + << Form(", TSRCS %d%d%d%2d%d %d%d%d%2d%d", + (Int_t) tDigi0->GetType(), + (Int_t) tDigi0->GetSm(), + (Int_t) tDigi0->GetRpc(), + (Int_t) tDigi0->GetChannel(), + (Int_t) tDigi0->GetSide(), + (Int_t) tDigi1->GetType(), + (Int_t) tDigi1->GetSm(), + (Int_t) tDigi1->GetRpc(), + (Int_t) tDigi1->GetChannel(), + (Int_t) tDigi1->GetSide()); + if (iDetIndx > (Int_t) fhCalWalk.size()) { LOG(error) << "Invalid DetIndx " << iDetIndx; continue; @@ -336,11 +379,8 @@ void CbmTofCalibrator::FillCalHist(CbmTofTracklet* pTrk, Int_t iOpt) { continue; } - const CbmTofDigi* tDigi1 = fDigiMan->Get<CbmTofDigi>(iDigInd1); - Int_t iCh1 = tDigi1->GetChannel(); - Int_t iSide1 = tDigi1->GetSide(); - LOG(debug) << "Fill Walk for " << iDetIndx << ", TSRCS " << iSmType << iSm - << iRpc << iCh1 << iSide1 << ", " << tDigi1 << ", " << pTrk; + Int_t iCh1 = tDigi1->GetChannel(); + Int_t iSide1 = tDigi1->GetSide(); if (iCh1 > (Int_t) fhCalWalk[iDetIndx].size()) { LOG(error) << "Invalid Ch1 " << iCh1; continue; @@ -352,8 +392,8 @@ void CbmTofCalibrator::FillCalHist(CbmTofTracklet* pTrk, Int_t iOpt) { if (iCh0 != iCh1 || iSide0 == iSide1) { LOG(error) << "Invalid digi pair for TSR " << iSmType << iSm << iRpc - << " Ch " << iCh0 << " " << iCh1 << ", Side " << iSide0 - << " " << iSide1; + << Form( + " Ch %2d %2d side %d %d", iCh0, iCh1, iSide0, iSide1); continue; } diff --git a/reco/detectors/tof/CbmTofCalibrator.h b/reco/detectors/tof/CbmTofCalibrator.h index ce596f2982..efc15c2b88 100644 --- a/reco/detectors/tof/CbmTofCalibrator.h +++ b/reco/detectors/tof/CbmTofCalibrator.h @@ -53,7 +53,7 @@ public: InitStatus Init(); Bool_t InitParameters(); Bool_t CreateCalHist(); - void FillCalHist(CbmTofTracklet* pTrk, Int_t iOpt); + void FillCalHist(CbmTofTracklet* pTrk, Int_t iOpt, CbmEvent* pEvent = NULL); Bool_t UpdateCalHist(Int_t iOpt); void ReadHist(TFile* fhFile); void WriteHist(TFile* fhFile); @@ -67,6 +67,8 @@ private: CbmTofFindTracks* fTofFindTracks; CbmTofTrackletTools* fTrackletTools; + const std::vector<CbmTofDigi>* fTofCalDigiVec = nullptr; + CbmTofDigiPar* fDigiPar; CbmTofDigiBdfPar* fDigiBdfPar; TClonesArray* fTofDigiMatchColl; // TOF Digi Links diff --git a/reco/detectors/tof/CbmTofEventClusterizer.cxx b/reco/detectors/tof/CbmTofEventClusterizer.cxx index abd11b6220..cc32cc309d 100644 --- a/reco/detectors/tof/CbmTofEventClusterizer.cxx +++ b/reco/detectors/tof/CbmTofEventClusterizer.cxx @@ -104,19 +104,13 @@ CbmTofEventClusterizer::CbmTofEventClusterizer(const char* name, , fTrbHeader(NULL) , fTofPointsColl(NULL) , fMcTracksColl(NULL) - , - //fTofDigisColl(NULL), - fDigiMan(nullptr) + , fDigiMan(nullptr) , fEventsColl(nullptr) , fbWriteHitsInOut(writeDataInOut) , fbWriteDigisInOut(writeDataInOut) - , - //fTofCalDigisColl(NULL), - fTofHitsColl(NULL) + , fTofHitsColl(NULL) , fTofDigiMatchColl(NULL) - , - //fTofCalDigisCollOut(NULL), - fTofHitsCollOut(NULL) + , fTofHitsCollOut(NULL) , fTofDigiMatchCollOut(NULL) , fiNbHits(0) , fVerbose(verbose) @@ -226,7 +220,7 @@ CbmTofEventClusterizer::CbmTofEventClusterizer(const char* name, , fTRefMode(0) , fTRefHits(0) , fIdMode(0) - , fDutId(0) + , fDutId(-1) , fDutSm(0) , fDutRpc(0) , fDutAddr(0) @@ -390,7 +384,7 @@ void CbmTofEventClusterizer::Exec(Option_t* option) { // --- In event-by-event mode: copy caldigis, hits and matches to output array and register them to event - // Int_t iDigi0=iNbCalDigis; //starting index of current event (VF) not used + Int_t iDigi0 = iNbCalDigis; //starting index of current event for (UInt_t index = 0; index < fTofCalDigiVec->size(); index++) { // for (Int_t index = 0; index < fTofCalDigisColl->GetEntriesFast(); index++){ CbmTofDigi* tDigi = &(fTofCalDigiVec->at(index)); @@ -407,14 +401,25 @@ void CbmTofEventClusterizer::Exec(Option_t* option) { tEvent->AddData(ECbmDataType::kTofHit, iNbHits); CbmMatch* pDigiMatch = (CbmMatch*) fTofDigiMatchColl->At(index); - // update content of match object, not necessary if event definition is kept ! - /* - for (Int_t iLink=0; iLink<pDigiMatch->GetNofLinks(); iLink++) { // loop over digis - CbmLink Link = pDigiMatch->GetLink(iLink); - Link.SetIndex(Link.GetIndex()+iDigi0); - } - */ - new ((*fTofDigiMatchCollOut)[iNbHits]) CbmMatch(*pDigiMatch); + + TString cstr = Form("Modify for hit %d, (%d", index, iNbHits); + cstr += Form(") %d links by iDigi0 %d: ", + (Int_t) pDigiMatch->GetNofLinks(), + iDigi0); + + // update content of match object to TS array + CbmMatch* mDigiMatch = new CbmMatch(); + for (Int_t iLink = 0; iLink < pDigiMatch->GetNofLinks(); iLink++) { + CbmLink Link = pDigiMatch->GetLink(iLink); + Link.SetIndex(Link.GetIndex() + iDigi0); + cstr += Form(" %d", (Int_t) Link.GetIndex()); + mDigiMatch->AddLink(Link); + } + //delete pDigiMatch; + + LOG(debug) << cstr; + + new ((*fTofDigiMatchCollOut)[iNbHits]) CbmMatch(*mDigiMatch); iNbHits++; } @@ -1171,8 +1176,18 @@ Bool_t CbmTofEventClusterizer::CreateHistos() { 4000, 0.0, 4.0); + // RPC related distributions + Int_t iNbDet = fDigiBdfPar->GetNbDet(); + fviDetId.resize(iNbDet); - /* TH2* hTemp = 0;*/ + fDetIdIndexMap.clear(); + for (Int_t iDetIndx = 0; iDetIndx < iNbDet; iDetIndx++) { + Int_t iUniqueId = fDigiBdfPar->GetDetUId(iDetIndx); + fDetIdIndexMap[iUniqueId] = iDetIndx; + fviDetId[iDetIndx] = iUniqueId; + } + + if (fDutId < 0) return kTRUE; // Sm related distributions fhSmCluPosition.resize(fDigiBdfPar->GetNbSmTypes()); @@ -1312,8 +1327,6 @@ Bool_t CbmTofEventClusterizer::CreateHistos() { } } - // RPC related distributions - Int_t iNbDet = fDigiBdfPar->GetNbDet(); LOG(info) << " Define Clusterizer histos for " << iNbDet << " detectors "; fhRpcDigiCor.resize(iNbDet); @@ -1342,17 +1355,12 @@ Bool_t CbmTofEventClusterizer::CreateHistos() { fhRpcDTLastHits.resize(iNbDet); fhRpcDTLastHits_Tot.resize(iNbDet); fhRpcDTLastHits_CluSize.resize(iNbDet); - fviDetId.resize(iNbDet); - - fDetIdIndexMap.clear(); if (fTotMax != 0.) fdTOTMax = fTotMax; if (fTotMin != 0.) fdTOTMin = fTotMin; for (Int_t iDetIndx = 0; iDetIndx < iNbDet; iDetIndx++) { - Int_t iUniqueId = fDigiBdfPar->GetDetUId(iDetIndx); - fDetIdIndexMap[iUniqueId] = iDetIndx; - fviDetId[iDetIndx] = iUniqueId; + Int_t iUniqueId = fDigiBdfPar->GetDetUId(iDetIndx); Int_t iSmType = CbmTofAddress::GetSmType(iUniqueId); Int_t iSmId = CbmTofAddress::GetSmId(iUniqueId); @@ -2388,6 +2396,9 @@ Bool_t CbmTofEventClusterizer::CreateHistos() { Bool_t CbmTofEventClusterizer::FillHistos() { fhClustBuildTime->Fill(fStop.GetSec() - fStart.GetSec() + (fStop.GetNanoSec() - fStart.GetNanoSec()) / 1e9); + + if (fDutId < 0) return kTRUE; + Int_t iNbTofHits = fTofHitsColl->GetEntries(); CbmTofHit* pHit; //gGeoManager->SetTopVolume( gGeoManager->FindVolumeFast("tof_v14a") ); @@ -3586,6 +3597,7 @@ Bool_t CbmTofEventClusterizer::FillHistos() { } Bool_t CbmTofEventClusterizer::WriteHistos() { + if (fDutId < 0) return kTRUE; TDirectory* oldir = gDirectory; TFile* fHist; fHist = new TFile(fOutHstFileName, "RECREATE"); @@ -4600,14 +4612,14 @@ Bool_t CbmTofEventClusterizer::WriteHistos() { TFitResultPtr fRes = hTy->Fit("gaus", "SQM0", "", dFMean - dFLim, dFMean + dFLim); LOG(debug) << "CalibF " - << Form("TSRC %d%d%d%d gaus %8.2f %8.2f %8.2f ", - iSmType, - iSm, - iRpc, - iCh, - fRes->Parameter(0), - fRes->Parameter(1), - fRes->Parameter(2)); + << Form("TSRC %d%d%d%d gaus %8.2f %8.2f %8.2f ", + iSmType, + iSm, + iRpc, + iCh, + fRes->Parameter(0), + fRes->Parameter(1), + fRes->Parameter(2)); TMean = fRes->Parameter(1); //overwrite mean } } @@ -5490,7 +5502,6 @@ Bool_t CbmTofEventClusterizer::BuildClusters() { } } iNbTofDigi = fTofDigiVec.size(); - //iNbTofDigi = fTofDigisColl->GetEntries(); // Update } if (kTRUE) { @@ -5505,7 +5516,7 @@ Bool_t CbmTofEventClusterizer::BuildClusters() { CbmTofDigi* pDigi = &(fTofDigiVec.at(iDigInd)); Int_t iDetIndx = fDigiBdfPar->GetDetInd(pDigi->GetAddress()); - LOG(debug) << iDigInd << " " << pDigi + LOG(debug) << "RawDigi" << iDigInd << " " << pDigi << Form(" Address : 0x%08x ", pDigi->GetAddress()) << " SmT " << pDigi->GetType() << " Sm " << pDigi->GetSm() << " Rpc " << pDigi->GetRpc() << " Ch " << pDigi->GetChannel() << " S " @@ -5521,143 +5532,150 @@ Bool_t CbmTofEventClusterizer::BuildClusters() { break; } - if (NULL == fhRpcDigiCor[iDetIndx]) { - if (100 < iMess++) - LOG(warning) << Form( - " DigiCor Histo for DetIndx %d derived from 0x%08x not found", - iDetIndx, - pDigi->GetAddress()); - continue; - } - - size_t iDigiCh = pDigi->GetChannel() * 2 + pDigi->GetSide(); - if (iDigiCh < fvTimeLastDigi[iDetIndx].size()) { + CbmTofDigi* pDigi2Min = NULL; + Double_t dTDifMin = dDoubleMax; - if (fvTimeLastDigi[iDetIndx][iDigiCh] > 0) { - if (fdStartAna10s > 0.) { - Double_t dTimeAna10s = (pDigi->GetTime() - fdStartAna10s) / 1.E9; - if (dTimeAna10s < fdSpillDuration) - fhRpcDigiDTLD[iDetIndx]->Fill( - iDigiCh, - (pDigi->GetTime() - fvTimeLastDigi[iDetIndx][iDigiCh]) / 1.E9); + if (fDutId > -1) { + if (fhRpcDigiCor.size() > 0) { + if (NULL == fhRpcDigiCor[iDetIndx]) { + if (100 < iMess++) + LOG(warning) << Form( + " DigiCor Histo for DetIndx %d derived from 0x%08x not found", + iDetIndx, + pDigi->GetAddress()); + continue; } - } - fvTimeLastDigi[iDetIndx][iDigiCh] = pDigi->GetTime(); + } else + break; - if (fvTimeFirstDigi[iDetIndx][iDigiCh] != 0.) { - fhRpcDigiDTFD[iDetIndx]->Fill( - iDigiCh, (pDigi->GetTime() - fvTimeFirstDigi[iDetIndx][iDigiCh])); - fvMulDigi[iDetIndx][iDigiCh]++; - } else { - fvTimeFirstDigi[iDetIndx][iDigiCh] = pDigi->GetTime(); - fvMulDigi[iDetIndx][iDigiCh]++; - } - } + size_t iDigiCh = pDigi->GetChannel() * 2 + pDigi->GetSide(); + if (iDigiCh < fvTimeLastDigi[iDetIndx].size()) { + if (fvTimeLastDigi[iDetIndx][iDigiCh] > 0) { + if (fdStartAna10s > 0.) { + Double_t dTimeAna10s = (pDigi->GetTime() - fdStartAna10s) / 1.E9; + if (dTimeAna10s < fdSpillDuration) + fhRpcDigiDTLD[iDetIndx]->Fill( + iDigiCh, + (pDigi->GetTime() - fvTimeLastDigi[iDetIndx][iDigiCh]) + / 1.E9); + } + } + fvTimeLastDigi[iDetIndx][iDigiCh] = pDigi->GetTime(); - Double_t dTDifMin = dDoubleMax; - CbmTofDigi* pDigi2Min = NULL; - // for (Int_t iDigI2 =iDigInd+1; iDigI2<iNbTofDigi;iDigI2++){ - for (Int_t iDigI2 = 0; iDigI2 < iNbTofDigi; iDigI2++) { - CbmTofDigi* pDigi2 = &(fTofDigiVec.at(iDigI2)); - // CbmTofDigi *pDigi2 = (CbmTofDigi*) fTofDigisColl->At( iDigI2 ); - // Fill digi correlation histogram per counter - if (iDetIndx == fDigiBdfPar->GetDetInd(pDigi2->GetAddress())) { - if (0. == pDigi->GetSide() && 1. == pDigi2->GetSide()) { - fhRpcDigiCor[iDetIndx]->Fill(pDigi->GetChannel(), - pDigi2->GetChannel()); + if (fvTimeFirstDigi[iDetIndx][iDigiCh] != 0.) { + fhRpcDigiDTFD[iDetIndx]->Fill( + iDigiCh, (pDigi->GetTime() - fvTimeFirstDigi[iDetIndx][iDigiCh])); + fvMulDigi[iDetIndx][iDigiCh]++; } else { - if (1. == pDigi->GetSide() && 0. == pDigi2->GetSide()) { - fhRpcDigiCor[iDetIndx]->Fill(pDigi2->GetChannel(), - pDigi->GetChannel()); - } + fvTimeFirstDigi[iDetIndx][iDigiCh] = pDigi->GetTime(); + fvMulDigi[iDetIndx][iDigiCh]++; } - if (pDigi->GetSide() != pDigi2->GetSide()) { - if (pDigi->GetChannel() == pDigi2->GetChannel()) { - Double_t dTDif = TMath::Abs(pDigi->GetTime() - pDigi2->GetTime()); - if (dTDif < dTDifMin) { - dTDifMin = dTDif; - pDigi2Min = pDigi2; - } - } else if ( - TMath::Abs(pDigi->GetChannel() - pDigi2->GetChannel()) - == 1) { // opposite side missing, neighbouring channel has hit on opposite side // FIXME - // check that same side digi of neighbouring channel is absent - Int_t iDigI3 = 0; - for (; iDigI3 < iNbTofDigi; iDigI3++) { - CbmTofDigi* pDigi3 = &(fTofDigiVec.at(iDigI3)); - // CbmTofDigi *pDigi3 = (CbmTofDigi*) fTofDigisColl->At( iDigI3 ); - if (pDigi3->GetSide() == pDigi->GetSide() - && pDigi2->GetChannel() == pDigi3->GetChannel()) - break; + } + + // for (Int_t iDigI2 =iDigInd+1; iDigI2<iNbTofDigi;iDigI2++){ + for (Int_t iDigI2 = 0; iDigI2 < iNbTofDigi; iDigI2++) { + CbmTofDigi* pDigi2 = &(fTofDigiVec.at(iDigI2)); + // CbmTofDigi *pDigi2 = (CbmTofDigi*) fTofDigisColl->At( iDigI2 ); + // Fill digi correlation histogram per counter + if (iDetIndx == fDigiBdfPar->GetDetInd(pDigi2->GetAddress())) { + if (0. == pDigi->GetSide() && 1. == pDigi2->GetSide()) { + fhRpcDigiCor[iDetIndx]->Fill(pDigi->GetChannel(), + pDigi2->GetChannel()); + } else { + if (1. == pDigi->GetSide() && 0. == pDigi2->GetSide()) { + fhRpcDigiCor[iDetIndx]->Fill(pDigi2->GetChannel(), + pDigi->GetChannel()); } - if (iDigI3 == iNbTofDigi) // same side neighbour did not fire - { - Int_t iCorMode = 0; // Missing hit correction mode - switch (iCorMode) { - case 0: // no action + } + if (pDigi->GetSide() != pDigi2->GetSide()) { + if (pDigi->GetChannel() == pDigi2->GetChannel()) { + Double_t dTDif = + TMath::Abs(pDigi->GetTime() - pDigi2->GetTime()); + if (dTDif < dTDifMin) { + dTDifMin = dTDif; + pDigi2Min = pDigi2; + } + } else if ( + TMath::Abs(pDigi->GetChannel() - pDigi2->GetChannel()) + == 1) { // opposite side missing, neighbouring channel has hit on opposite side // FIXME + // check that same side digi of neighbouring channel is absent + Int_t iDigI3 = 0; + for (; iDigI3 < iNbTofDigi; iDigI3++) { + CbmTofDigi* pDigi3 = &(fTofDigiVec.at(iDigI3)); + // CbmTofDigi *pDigi3 = (CbmTofDigi*) fTofDigisColl->At( iDigI3 ); + if (pDigi3->GetSide() == pDigi->GetSide() + && pDigi2->GetChannel() == pDigi3->GetChannel()) break; - case 1: // shift found hit - LOG(debug2) << Form("shift channel %d%d%d%d%d and ", - (Int_t) pDigi->GetType(), - (Int_t) pDigi->GetSm(), - (Int_t) pDigi->GetRpc(), - (Int_t) pDigi->GetChannel(), - (Int_t) pDigi->GetSide()) - << Form(" %d%d%d%d%d ", - (Int_t) pDigi2->GetType(), - (Int_t) pDigi2->GetSm(), - (Int_t) pDigi2->GetRpc(), - (Int_t) pDigi2->GetChannel(), - (Int_t) pDigi2->GetSide()); - //if(pDigi->GetTime() < pDigi2->GetTime()) - if (pDigi->GetSide() == 0) - pDigi2->SetAddress(pDigi->GetSm(), + } + if (iDigI3 == iNbTofDigi) // same side neighbour did not fire + { + Int_t iCorMode = 0; // Missing hit correction mode + switch (iCorMode) { + case 0: // no action + break; + case 1: // shift found hit + LOG(debug2) << Form("shift channel %d%d%d%d%d and ", + (Int_t) pDigi->GetType(), + (Int_t) pDigi->GetSm(), + (Int_t) pDigi->GetRpc(), + (Int_t) pDigi->GetChannel(), + (Int_t) pDigi->GetSide()) + << Form(" %d%d%d%d%d ", + (Int_t) pDigi2->GetType(), + (Int_t) pDigi2->GetSm(), + (Int_t) pDigi2->GetRpc(), + (Int_t) pDigi2->GetChannel(), + (Int_t) pDigi2->GetSide()); + //if(pDigi->GetTime() < pDigi2->GetTime()) + if (pDigi->GetSide() == 0) + pDigi2->SetAddress(pDigi->GetSm(), + pDigi->GetRpc(), + pDigi->GetChannel(), + 1 - pDigi->GetSide(), + pDigi->GetType()); + else + pDigi->SetAddress(pDigi2->GetSm(), + pDigi2->GetRpc(), + pDigi2->GetChannel(), + 1 - pDigi2->GetSide(), + pDigi2->GetType()); + + LOG(debug2) << Form("resultchannel %d%d%d%d%d and ", + (Int_t) pDigi->GetType(), + (Int_t) pDigi->GetSm(), + (Int_t) pDigi->GetRpc(), + (Int_t) pDigi->GetChannel(), + (Int_t) pDigi->GetSide()) + << Form(" %d%d%d%d%d ", + (Int_t) pDigi2->GetType(), + (Int_t) pDigi2->GetSm(), + (Int_t) pDigi2->GetRpc(), + (Int_t) pDigi2->GetChannel(), + (Int_t) pDigi2->GetSide()); + break; + case 2: // insert missing hits + fTofDigiVec.push_back(CbmTofDigi(*pDigi)); + CbmTofDigi* pDigiN = &(fTofDigiVec.back()); + // CbmTofDigi *pDigiN = new((*fTofDigisColl)[iNbTofDigi++]) CbmTofDigi( *pDigi ); + pDigiN->SetAddress(pDigi->GetSm(), pDigi->GetRpc(), - pDigi->GetChannel(), - 1 - pDigi->GetSide(), + pDigi2->GetChannel(), + pDigi->GetSide(), pDigi->GetType()); - else - pDigi->SetAddress(pDigi2->GetSm(), - pDigi2->GetRpc(), - pDigi2->GetChannel(), - 1 - pDigi2->GetSide(), - pDigi2->GetType()); - - LOG(debug2) << Form("resultchannel %d%d%d%d%d and ", - (Int_t) pDigi->GetType(), - (Int_t) pDigi->GetSm(), - (Int_t) pDigi->GetRpc(), - (Int_t) pDigi->GetChannel(), - (Int_t) pDigi->GetSide()) - << Form(" %d%d%d%d%d ", - (Int_t) pDigi2->GetType(), - (Int_t) pDigi2->GetSm(), - (Int_t) pDigi2->GetRpc(), - (Int_t) pDigi2->GetChannel(), - (Int_t) pDigi2->GetSide()); - break; - case 2: // insert missing hits - fTofDigiVec.push_back(CbmTofDigi(*pDigi)); - CbmTofDigi* pDigiN = &(fTofDigiVec.back()); - // CbmTofDigi *pDigiN = new((*fTofDigisColl)[iNbTofDigi++]) CbmTofDigi( *pDigi ); - pDigiN->SetAddress(pDigi->GetSm(), - pDigi->GetRpc(), - pDigi2->GetChannel(), - pDigi->GetSide(), - pDigi->GetType()); - pDigiN->SetTot(pDigi2->GetTot()); - - fTofDigiVec.push_back(CbmTofDigi(*pDigi2)); - CbmTofDigi* pDigiN2 = &(fTofDigiVec.back()); - // CbmTofDigi *pDigiN2 = new((*fTofDigisColl)[iNbTofDigi++]) CbmTofDigi( *pDigi2 ); - pDigiN2->SetAddress(pDigi2->GetSm(), - pDigi2->GetRpc(), - pDigi->GetChannel(), - pDigi2->GetSide(), - pDigi2->GetType()); - pDigiN2->SetTot(pDigi->GetTot()); + pDigiN->SetTot(pDigi2->GetTot()); + + fTofDigiVec.push_back(CbmTofDigi(*pDigi2)); + CbmTofDigi* pDigiN2 = &(fTofDigiVec.back()); + // CbmTofDigi *pDigiN2 = new((*fTofDigisColl)[iNbTofDigi++]) CbmTofDigi( *pDigi2 ); + pDigiN2->SetAddress(pDigi2->GetSm(), + pDigi2->GetRpc(), + pDigi->GetChannel(), + pDigi2->GetSide(), + pDigi2->GetType()); + pDigiN2->SetTot(pDigi->GetTot()); - break; + break; + } } } } @@ -5719,7 +5737,7 @@ Bool_t CbmTofEventClusterizer::BuildClusters() { } // kTRUE end // Calibrate RawDigis - if (kTRUE == fDigiBdfPar->UseExpandedDigi()) { + if (kTRUE) { CbmTofDigi* pDigi; // CbmTofDigi *pCalDigi=NULL; (VF) not used CalibRawDigis(); @@ -5766,67 +5784,68 @@ Bool_t CbmTofEventClusterizer::BuildClusters() { } // for( Int_t iDigInd = 0; iDigInd < nTofDigi; iDigInd++ ) // inspect digi array - Int_t iNbDet = fDigiBdfPar->GetNbDet(); - for (Int_t iDetIndx = 0; iDetIndx < iNbDet; iDetIndx++) { - Int_t iDetId = fviDetId[iDetIndx]; - Int_t iSmType = CbmTofAddress::GetSmType(iDetId); - Int_t iSm = CbmTofAddress::GetSmId(iDetId); - Int_t iRpc = CbmTofAddress::GetRpcId(iDetId); - Int_t iNbStrips = fDigiBdfPar->GetNbChan(iSmType, iRpc); - for (Int_t iStrip = 0; iStrip < iNbStrips; iStrip++) { - Int_t iDigiMul = fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) - + iRpc][iStrip] - .size(); - //LOG(info)<<"Inspect TSRC "<<iSmType<<iSm<<iRpc<<iStrip<<" with "<<iNbStrips<<" strips: Mul "<<iDigiMul; - if (iDigiMul > 0) { - fhRpcDigiMul[iDetIndx]->Fill(iStrip, iDigiMul); - if (iDigiMul == 1) { - fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 0); - if (iStrip > 0) - if (fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) - + iRpc][iStrip - 1] - .size() - > 1) { - fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 1); - if (TMath::Abs( - fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) - + iRpc][iStrip][0] - ->GetTime() - - fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) - + iRpc][iStrip - 1][0] - ->GetTime()) - < fMaxTimeDist) - fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 3); - } - if (iStrip < iNbStrips - 2) { - if (fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) - + iRpc][iStrip + 1] - .size() - > 1) { - fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 2); - if (TMath::Abs( - fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) - + iRpc][iStrip][0] - ->GetTime() - - fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) - + iRpc][iStrip + 1][0] - ->GetTime()) - < fMaxTimeDist) - fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 4); + if (fDutId > -1) { + Int_t iNbDet = fDigiBdfPar->GetNbDet(); + for (Int_t iDetIndx = 0; iDetIndx < iNbDet; iDetIndx++) { + Int_t iDetId = fviDetId[iDetIndx]; + Int_t iSmType = CbmTofAddress::GetSmType(iDetId); + Int_t iSm = CbmTofAddress::GetSmId(iDetId); + Int_t iRpc = CbmTofAddress::GetRpcId(iDetId); + Int_t iNbStrips = fDigiBdfPar->GetNbChan(iSmType, iRpc); + for (Int_t iStrip = 0; iStrip < iNbStrips; iStrip++) { + Int_t iDigiMul = + fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) + iRpc] + [iStrip] + .size(); + //LOG(info)<<"Inspect TSRC "<<iSmType<<iSm<<iRpc<<iStrip<<" with "<<iNbStrips<<" strips: Mul "<<iDigiMul; + if (iDigiMul > 0) { + fhRpcDigiMul[iDetIndx]->Fill(iStrip, iDigiMul); + if (iDigiMul == 1) { + fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 0); + if (iStrip > 0) + if (fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) + + iRpc][iStrip - 1] + .size() + > 1) { + fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 1); + if (TMath::Abs( + fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) + + iRpc][iStrip][0] + ->GetTime() + - fStorDigi[iSmType] + [iSm * fDigiBdfPar->GetNbRpc(iSmType) + iRpc] + [iStrip - 1][0] + ->GetTime()) + < fMaxTimeDist) + fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 3); + } + if (iStrip < iNbStrips - 2) { + if (fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) + + iRpc][iStrip + 1] + .size() + > 1) { + fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 2); + if (TMath::Abs( + fStorDigi[iSmType][iSm * fDigiBdfPar->GetNbRpc(iSmType) + + iRpc][iStrip][0] + ->GetTime() + - fStorDigi[iSmType] + [iSm * fDigiBdfPar->GetNbRpc(iSmType) + iRpc] + [iStrip + 1][0] + ->GetTime()) + < fMaxTimeDist) + fhRpcDigiStatus[iDetIndx]->Fill(iStrip, 4); + } } } } } } } - BuildHits(); - } // if( kTRUE == fDigiBdfPar->UseExpandedDigi() ) - else { - LOG(error) << " Compressed Digis not implemented ... "; - return kFALSE; // not implemented properly yet - } + } // if( kTRUE ) obsolete ) + return kTRUE; } @@ -5841,7 +5860,7 @@ Bool_t CbmTofEventClusterizer::MergeClusters() { CbmTofHit* pHit = (CbmTofHit*) fTofHitsColl->At(iHitInd); if (NULL == pHit) continue; - Int_t iDetId = (pHit->GetAddress() & SelMask); + Int_t iDetId = (pHit->GetAddress() & DetMask); Int_t iSmType = CbmTofAddress::GetSmType(iDetId); Int_t iNbRpc = fDigiBdfPar->GetNbRpc(iSmType); if (iSmType != 5 && iSmType != 8) continue; // only merge diamonds and Pad @@ -5865,7 +5884,7 @@ Bool_t CbmTofEventClusterizer::MergeClusters() { iHitInd2++) { CbmTofHit* pHit2 = (CbmTofHit*) fTofHitsColl->At(iHitInd2); if (NULL == pHit2) continue; - Int_t iDetId2 = (pHit2->GetAddress() & SelMask); + Int_t iDetId2 = (pHit2->GetAddress() & DetMask); Int_t iSmType2 = CbmTofAddress::GetSmType(iDetId2); if (iSmType2 == iSmType) { Int_t iSm2 = CbmTofAddress::GetSmId(iDetId2); @@ -6208,7 +6227,9 @@ Bool_t CbmTofEventClusterizer::AddNextChan(Int_t iSmType, if (iCh == iNbCh) return kFALSE; if (0 == fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()) return kFALSE; if (0 < fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()) - fhNbDigiPerChan->Fill(fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()); + if (fDutId > -1) + fhNbDigiPerChan->Fill( + fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()); if (1 < fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()) { Bool_t AddedHit = kFALSE; for (size_t i1 = 0; @@ -6493,7 +6514,7 @@ Bool_t CbmTofEventClusterizer::BuildHits() { Double_t dTimeDif = 0.0; Double_t dTotS = 0.0; fiNbSameSide = 0; - if (kTRUE == fDigiBdfPar->UseExpandedDigi()) { + if (kTRUE) { for (Int_t iSmType = 0; iSmType < iNbSmTypes; iSmType++) { Int_t iNbSm = fDigiBdfPar->GetNbSm(iSmType); Int_t iNbRpc = fDigiBdfPar->GetNbRpc(iSmType); @@ -6549,8 +6570,9 @@ Bool_t CbmTofEventClusterizer::BuildHits() { if (fvDeadStrips[iDetIndx] & (1 << iCh)) continue; // skip over dead channels if (0 < fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()) - fhNbDigiPerChan->Fill( - fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()); + if (fDutId > -1) + fhNbDigiPerChan->Fill( + fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()); while (1 < fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].size()) { @@ -6813,6 +6835,9 @@ Bool_t CbmTofEventClusterizer::BuildHits() { if (TMath::Abs(dPosY) > fChannelInfo->GetSizey() * fPosYMaxScal) { // remove both digis + LOG(debug1) << "Remove digis on TSRC " << iSmType << iSm + << iRpc << iCh << " with dPosY " << dPosY + << " > " << fChannelInfo->GetSizey(); fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].erase( fStorDigi[iSmType][iSm * iNbRpc + iRpc][iCh].begin()); fStorDigiInd[iSmType][iSm * iNbRpc + iRpc][iCh].erase( @@ -6852,10 +6877,12 @@ Bool_t CbmTofEventClusterizer::BuildHits() { // Now check if a hit/cluster is already started if (0 < iNbChanInHit) { if (iLastChan == iCh - 1) { - fhDigTimeDifClust->Fill(dTime - dLastTime); - fhDigSpacDifClust->Fill(dPosY - dLastPosY); - fhDigDistClust->Fill(dPosY - dLastPosY, - dTime - dLastTime); + if (fDutId > -1) { + fhDigTimeDifClust->Fill(dTime - dLastTime); + fhDigSpacDifClust->Fill(dPosY - dLastPosY); + fhDigDistClust->Fill(dPosY - dLastPosY, + dTime - dLastTime); + } } // if( iLastChan == iCh - 1 ) // a cluster is already started => check distance in space/time @@ -7077,10 +7104,7 @@ Bool_t CbmTofEventClusterizer::BuildHits() { } else { pHit->Delete(); } - /* - new((*fTofDigiMatchColl)[fiNbHits]) CbmMatch(); - CbmMatch* digiMatch = (CbmMatch *)fTofDigiMatchColl->At(fiNbHits); - */ + CbmMatch* digiMatch = new ((*fTofDigiMatchColl)[fiNbHits]) CbmMatch(); for (size_t i = 0; i < vDigiIndRef.size(); i++) { @@ -7311,25 +7335,30 @@ Bool_t CbmTofEventClusterizer::BuildHits() { if (iChm > iNbCh - 1) iChm = iNbCh - 1; iDetId = CbmTofAddress::GetUniqueAddress(iSm, iRpc, iChm, 0, iSmType); - Int_t iRefId = 0; // Index of the correspondng TofPoint + //Int_t iRefId = 0; // Index of the correspondng TofPoint //if(NULL != fTofPointsColl) iRefId = fTofPointsColl->IndexOf( vPtsRef[0] ); TString cstr = "Save V-Hit "; - cstr += Form(" %3d %3d 0x%08x %3d 0x%08x %8.2f %6.2f", // %3d %3d - fiNbHits, - iNbChanInHit, - iDetId, - iLastChan, - iRefId, //vPtsRef.size(),vPtsRef[0]) - dWeightedTime, - dWeightedPosY); - - cstr += Form(", DigiSize: %lu ", vDigiIndRef.size()); + cstr += Form( + " %3d %3d 0x%08x TSR %d%d%d Ch %2d %8.2f %6.2f", // %3d %3d + fiNbHits, + iNbChanInHit, + iDetId, + iSmType, + iSm, + iRpc, + iChm, + dWeightedTime, + dWeightedPosY); + + cstr += Form(", DigiSize: %lu (%3lu)", + vDigiIndRef.size(), + fTofCalDigiVec->size()); cstr += ", DigiInds: "; fviClusterMul[iSmType][iSm][iRpc]++; for (UInt_t i = 0; i < vDigiIndRef.size(); i++) { - cstr += Form(" %d (M,%d)", + cstr += Form(" %3d (M,%2d)", vDigiIndRef.at(i), fviClusterMul[iSmType][iSm][iRpc]); } @@ -7367,10 +7396,7 @@ Bool_t CbmTofEventClusterizer::BuildHits() { } else { pHit->Delete(); } - /* - new((*fTofDigiMatchColl)[fiNbHits]) CbmMatch(); - CbmMatch* digiMatch = (CbmMatch *)fTofDigiMatchColl->At(fiNbHits); - */ + CbmMatch* digiMatch = new ((*fTofDigiMatchColl)[fiNbHits]) CbmMatch(); diff --git a/reco/detectors/tof/CbmTofFindTracks.cxx b/reco/detectors/tof/CbmTofFindTracks.cxx index ec1493a06f..f3cbf06794 100644 --- a/reco/detectors/tof/CbmTofFindTracks.cxx +++ b/reco/detectors/tof/CbmTofFindTracks.cxx @@ -82,6 +82,7 @@ CbmTofFindTracks::CbmTofFindTracks(const char* name, , fTofHitArrayIn(NULL) , fTofMatchArrayIn(NULL) , fTofHitArray(NULL) + , fTofHitIndexArray() , fTrackArray(NULL) , fTrackArrayOut(nullptr) , fTofUHitArray(NULL) @@ -246,7 +247,7 @@ InitStatus CbmTofFindTracks::Init() { return kERROR; } - // Get TOF hit Array + // Get TOF DigiMatch Array fTofMatchArrayIn = (TClonesArray*) ioman->GetObject("TofCalDigiMatch"); if (!fTofMatchArrayIn) { LOG(fatal) << "CbmTofFindTracks::Init: No TofDigiMatch array!"; @@ -1058,16 +1059,18 @@ void CbmTofFindTracks::Exec(Option_t* opt) { if (fTofHitArray) fTofHitArray->Clear("C"); Int_t iNbHits = 0; + fTofHitIndexArray.resize(tEvent->GetNofData(ECbmDataType::kTofHit)); for (Int_t iHit = 0; iHit < tEvent->GetNofData(ECbmDataType::kTofHit); iHit++) { Int_t iHitIndex = static_cast<Int_t>(tEvent->GetIndex(ECbmDataType::kTofHit, iHit)); CbmTofHit* tHit = dynamic_cast<CbmTofHit*>(fTofHitArrayIn->At(iHitIndex)); + fTofHitIndexArray[iNbHits] = iHitIndex; new ((*fTofHitArray)[iNbHits++]) CbmTofHit(*tHit); } - ExecFind(opt); + ExecFind(opt, tEvent); // --- In event-by-event mode: copy tracks to output array and register them to event for (Int_t iTrk = 0; iTrk < fTrackArray->GetEntries(); iTrk++) { @@ -1081,7 +1084,7 @@ void CbmTofFindTracks::Exec(Option_t* opt) { } } -void CbmTofFindTracks::ExecFind(Option_t* /*opt*/) { +void CbmTofFindTracks::ExecFind(Option_t* /*opt*/, CbmEvent* tEvent) { fiEvent++; ResetStationsFired(); if (NULL != fTofUHitArray) fTofUHitArray->Clear("C"); @@ -1220,7 +1223,7 @@ void CbmTofFindTracks::ExecFind(Option_t* /*opt*/) { FindVertex(); - FillHistograms(); + FillHistograms(tEvent); } FillUHits(); // put unused hits into TClonesArray @@ -1780,7 +1783,7 @@ void CbmTofFindTracks::FindVertex() { static Int_t iWarnNotDefined = 0; -void CbmTofFindTracks::FillHistograms() { +void CbmTofFindTracks::FillHistograms(CbmEvent* tEvent) { // Locate reference ("beam counter") hit CbmTofHit* pRefHit = NULL; Double_t RefMinTime = 1.E300; @@ -1854,7 +1857,7 @@ void CbmTofFindTracks::FillHistograms() { iTMul++; fhTrklChi2->Fill(pTrk->GetNofHits(), pTrk->GetChiSq()); - if (fiCalOpt > 0) fTofCalibrator->FillCalHist(pTrk, fiCalOpt); + if (fiCalOpt > 0) fTofCalibrator->FillCalHist(pTrk, fiCalOpt, tEvent); CbmTofTrackletParam* tPar = pTrk->GetTrackParameter(); Double_t dTt = pTrk->GetTt(); diff --git a/reco/detectors/tof/CbmTofFindTracks.h b/reco/detectors/tof/CbmTofFindTracks.h index b912459754..6a1bfe1f3e 100644 --- a/reco/detectors/tof/CbmTofFindTracks.h +++ b/reco/detectors/tof/CbmTofFindTracks.h @@ -39,6 +39,7 @@ class CbmTofDigiBdfPar; class CbmTofAddress; class CbmTofHit; class CbmMatch; +class CbmEvent; class CbmTofFindTracks : public FairTask { friend class CbmTofTrackFinderNN; @@ -72,7 +73,7 @@ public: /** Task execution **/ virtual void Exec(Option_t* opt); - virtual void ExecFind(Option_t* opt); + virtual void ExecFind(Option_t* opt, CbmEvent* tEvent = NULL); /** Finish at the end of each event **/ virtual void Finish(); @@ -88,7 +89,7 @@ public: virtual void FindVertex(); - virtual void FillHistograms(); + virtual void FillHistograms(CbmEvent* tEvent = NULL); /** Accessors **/ CbmTofTrackFinder* GetFinder() { return fFinder; }; @@ -173,6 +174,13 @@ public: inline Double_t GetVertexY() const { return fVTX_Y; } inline Double_t GetVertexZ() const { return fVTX_Z; } + inline Int_t GetTofHitIndex(Int_t iHit) { + if (fTofHitIndexArray.size() < 1) + return iHit; + else + return fTofHitIndexArray[iHit]; + } + private: static CbmTofFindTracks* fInstance; CbmTofTrackFinder* fFinder; // Pointer to TrackFinder concrete class @@ -183,7 +191,8 @@ private: TClonesArray* fTofHitArrayIn; // Input array of TOF hits TClonesArray* fTofMatchArrayIn; // Input array of TOF hit matches TClonesArray* fTofHitArray; // Output array of recalibrated TOF hits - TClonesArray* fTrackArray; // Output array of CbmTofTracks + std::vector<Int_t> fTofHitIndexArray; // Index of hit in TS + TClonesArray* fTrackArray; // Output array of CbmTofTracks TClonesArray* fTrackArrayOut; // Output array of CbmTofTracks in CbmEvent mode TClonesArray* fTofUHitArray; // Output array of unused TOF hits -- GitLab