diff --git a/.vscode/settings.json b/.vscode/settings.json index 1c692d7..5e66083 100644 --- a/.vscode/settings.json +++ b/.vscode/settings.json @@ -126,5 +126,6 @@ "QQQStage2.C": "cpp", "inspect.C": "cpp" }, - "github-enterprise.uri": "https://fsunuc.physics.fsu.edu" + "github-enterprise.uri": "https://fsunuc.physics.fsu.edu", + "C_Cpp.default.compilerPath": "/usr/bin/gcc" } \ No newline at end of file diff --git a/MakeVertex.C b/MakeVertex.C index cb0cb5e..c88f65e 100755 --- a/MakeVertex.C +++ b/MakeVertex.C @@ -26,28 +26,28 @@ PW pw_contr; PW pwinstance; TVector3 hitPos; double qqqenergy, qqqtimestamp; -class Event { +class Event +{ public: - Event(TVector3 p, double e1, double e2, double t1, double t2) : pos(p), Energy1(e1), Energy2(e2), Time1(t1), Time2(t2) {} - Event(TVector3 p, double e1, double e2, double t1, double t2, int c1, int c2) : pos(p), Energy1(e1), Energy2(e2), Time1(t1), Time2(t2), ch1(c1), ch2(c2) {} + Event(TVector3 p, double e1, double e2, double t1, double t2) : pos(p), Energy1(e1), Energy2(e2), Time1(t1), Time2(t2) {} + Event(TVector3 p, double e1, double e2, double t1, double t2, int c1, int c2) : pos(p), Energy1(e1), Energy2(e2), Time1(t1), Time2(t2), ch1(c1), ch2(c2) {} - TVector3 pos; - int ch1=-1; //int(ch1/16) gives qqq id, ch1%16 gives ring# - int ch2=-1; //int(ch2/16) gives qqq id, ch2%16 gives wedge# - double Energy1=-1; //Front for QQQ, Anode for PC - double Energy2=-1; //Back for QQQ, Cathode for PC - double Time1=-1; - double Time2=-1; - + TVector3 pos; + int ch1 = -1; // int(ch1/16) gives qqq id, ch1%16 gives ring# + int ch2 = -1; // int(ch2/16) gives qqq id, ch2%16 gives wedge# + double Energy1 = -1; // Front for QQQ, Anode for PC + double Energy2 = -1; // Back for QQQ, Cathode for PC + double Time1 = -1; + double Time2 = -1; }; // Calibration globals const int MAX_QQQ = 4; const int MAX_RING = 16; const int MAX_WEDGE = 16; -const double qqqpos=100.0; -const double vertexpos=14.2; -const double pcrad=37.0; +const double qqqpos = 100.0; +const double vertexpos = 14.2; +const double pcrad = 37.0; double qqqGain[MAX_QQQ][MAX_RING][MAX_WEDGE] = {{{0}}}; bool qqqGainValid[MAX_QQQ][MAX_RING][MAX_WEDGE] = {{{false}}}; double qqqCalib[MAX_QQQ][MAX_RING][MAX_WEDGE] = {{{0}}}; @@ -109,7 +109,7 @@ void MakeVertex::Begin(TTree * /*tree*/) Crossover[i][j][0].y = pwinstance.An[i].first.Y() + alpha * a.Y(); Crossover[i][j][0].z = pwinstance.An[i].first.Z() + alpha * a.Z(); - if (Crossover[i][j][0].z < -190 || Crossover[i][j][0].z > 190 || (i+j)%24 == 12) + if (Crossover[i][j][0].z < -190 || Crossover[i][j][0].z > 190 || (i + j) % 24 == 12) { Crossover[i][j][0].z = 9999999; } @@ -195,30 +195,30 @@ void MakeVertex::Begin(TTree * /*tree*/) infile.close(); } } - + { - std::ifstream infile("sx3cal/backgains.dat"); - std::string temp; - int backpos, frontpos, clkpos; - std::cout << "foo" << std::endl; - if (infile.is_open()) - while(infile>>clkpos>>temp>>frontpos>>temp>>backpos>>sx3BackGain[clkpos][frontpos][backpos]) - std::cout << sx3BackGain[clkpos][frontpos][backpos] << std::endl; - infile.close(); + std::ifstream infile("sx3cal/backgains.dat"); + std::string temp; + int backpos, frontpos, clkpos; + std::cout << "foo" << std::endl; + if (infile.is_open()) + while (infile >> clkpos >> temp >> frontpos >> temp >> backpos >> sx3BackGain[clkpos][frontpos][backpos]) + std::cout << sx3BackGain[clkpos][frontpos][backpos] << std::endl; + infile.close(); - infile.open("sx3cal/frontgains.dat"); - if (infile.is_open()) - while(infile>>clkpos>>temp>>temp>>frontpos>>sx3FrontOffset[clkpos][frontpos]>>sx3FrontGain[clkpos][frontpos]) - std::cout << sx3FrontOffset[clkpos][frontpos] << " " << sx3FrontGain[clkpos][frontpos] << std::endl; - infile.close(); - - infile.open("sx3cal/rightgains.dat"); - if (infile.is_open()) - while(infile>>clkpos>>frontpos>>temp>>sx3RightGain[clkpos][frontpos]) { - sx3RightGain[clkpos][frontpos]=TMath::Abs(sx3RightGain[clkpos][frontpos]); - } - infile.close(); + infile.open("sx3cal/frontgains.dat"); + if (infile.is_open()) + while (infile >> clkpos >> temp >> temp >> frontpos >> sx3FrontOffset[clkpos][frontpos] >> sx3FrontGain[clkpos][frontpos]) + std::cout << sx3FrontOffset[clkpos][frontpos] << " " << sx3FrontGain[clkpos][frontpos] << std::endl; + infile.close(); + infile.open("sx3cal/rightgains.dat"); + if (infile.is_open()) + while (infile >> clkpos >> frontpos >> temp >> sx3RightGain[clkpos][frontpos]) + { + sx3RightGain[clkpos][frontpos] = TMath::Abs(sx3RightGain[clkpos][frontpos]); + } + infile.close(); } std::cout << "aaa" << std::endl; } @@ -227,7 +227,7 @@ Bool_t MakeVertex::Process(Long64_t entry) { hitPos.Clear(); qqqenergy = -1; - qqqtimestamp=-1; + qqqtimestamp = -1; HitNonZero = false; bool qqq1000cut = false; b_sx3Multi->GetEntry(entry); @@ -252,7 +252,7 @@ Bool_t MakeVertex::Process(Long64_t entry) std::vector sx3Events; // if(sx3.multi>1) { - // std::array Fsx3; + // std::array Fsx3; // //std::cout << "-----" << std::endl; // for(int i=0; i=12) continue; @@ -298,7 +298,7 @@ Bool_t MakeVertex::Process(Long64_t entry) // } // } // } - //return kTRUE; + // return kTRUE; // QQQ Processing int qqqCount = 0; @@ -322,45 +322,59 @@ Bool_t MakeVertex::Process(Long64_t entry) std::vector QQQ_Events, PC_Events; std::vector QQQ_Events_Raw, PC_Events_Raw; - std::vector QQQ_Events2; //clustering done + std::vector QQQ_Events2; // clustering done - std::unordered_map> qvecr[4], qvecw[4]; - if(qqq.multi>1) { - //if(qqq.multi>=3) std::cout << "-----" << std::endl; - for(int i=0; i=3) std::cout << std::setprecision(16) << "qqq"<< qqq.id[i] << " " << std::string(qqq.ch[i]/16?"ring":"wedge") << qqq.ch[i]%16 << " " << qqq.e[i] << " " << qqq.t[i] - qqq.t[0] << std::endl; - if(qqq.ch[i]/16) { - if(qvecr[qqq.id[i]].find(qqq.ch[i])!=qvecr[qqq.id[i]].end()) std::cout << "mayday!" << std::endl; - qvecr[qqq.id[i]][qqq.ch[i]] = std::tuple(qqq.id[i],qqq.ch[i],qqq.e[i],qqq.t[i]); - } else { - if(qvecw[qqq.id[i]].find(qqq.ch[i])!=qvecw[qqq.id[i]].end()) std::cout << "mayday!" << std::endl; - qvecw[qqq.id[i]][qqq.ch[i]] = std::tuple(qqq.id[i],qqq.ch[i],qqq.e[i],qqq.t[i]); - } - } + std::unordered_map> qvecr[4], qvecw[4]; + if (qqq.multi > 1) + { + // if(qqq.multi>=3) std::cout << "-----" << std::endl; + for (int i = 0; i < qqq.multi; i++) + { + // if(qqq.multi>=3) std::cout << std::setprecision(16) << "qqq"<< qqq.id[i] << " " << std::string(qqq.ch[i]/16?"ring":"wedge") << qqq.ch[i]%16 << " " << qqq.e[i] << " " << qqq.t[i] - qqq.t[0] << std::endl; + if (qqq.ch[i] / 16) + { + if (qvecr[qqq.id[i]].find(qqq.ch[i]) != qvecr[qqq.id[i]].end()) + std::cout << "mayday!" << std::endl; + qvecr[qqq.id[i]][qqq.ch[i]] = std::tuple(qqq.id[i], qqq.ch[i], qqq.e[i], qqq.t[i]); + } + else + { + if (qvecw[qqq.id[i]].find(qqq.ch[i]) != qvecw[qqq.id[i]].end()) + std::cout << "mayday!" << std::endl; + qvecw[qqq.id[i]][qqq.ch[i]] = std::tuple(qqq.id[i], qqq.ch[i], qqq.e[i], qqq.t[i]); + } + } } - + bool PCQQQTimeCut = false; - for (int i = 0; i < qqq.multi; i++) { + for (int i = 0; i < qqq.multi; i++) + { plotter->Fill2D("QQQ_Index_Vs_Energy", 16 * 8, 0, 16 * 8, 2000, 0, 16000, qqq.index[i], qqq.e[i], "hRawQQQ"); - for (int j = 0; j < qqq.multi; j++) { + for (int j = 0; j < qqq.multi; j++) + { if (j == i) continue; plotter->Fill2D("QQQ_Coincidence_Matrix", 16 * 8, 0, 16 * 8, 16 * 8, 0, 16 * 8, qqq.index[i], qqq.index[j], "hRawQQQ"); } - for (int k = 0; k < pc.multi; k++) { - if (pc.index[k] < 24 && pc.e[k] > 50) { + for (int k = 0; k < pc.multi; k++) + { + if (pc.index[k] < 24 && pc.e[k] > 50) + { plotter->Fill2D("QQQ_Vs_Anode_Energy", 400, 0, 4000, 1000, 0, 16000, qqq.e[i], pc.e[k], "hRawQQQ"); plotter->Fill2D("QQQ_Vs_PC_Index", 16 * 8, 0, 16 * 8, 24, 0, 24, qqq.index[i], pc.index[k], "hRawQQQ"); } - else if (pc.index[k] >= 24 && pc.e[k] > 50) { + else if (pc.index[k] >= 24 && pc.e[k] > 50) + { plotter->Fill2D("QQQ_Vs_Cathode_Energy", 400, 0, 4000, 1000, 0, 16000, qqq.e[i], pc.e[k], "hRawQQQ"); } } - for (int j = i + 1; j < qqq.multi; j++) { - if (qqq.id[i] == qqq.id[j]) { + for (int j = i + 1; j < qqq.multi; j++) + { + if (qqq.id[i] == qqq.id[j]) + { qqqCount++; int chWedge = -1; @@ -372,7 +386,8 @@ Bool_t MakeVertex::Process(Long64_t entry) double tRing = 0.0; double tWedge = 0.0; - if (qqq.ch[i] < 16 && qqq.ch[j] >= 16 && qqqGainValid[qqq.id[i]][qqq.ch[i]][qqq.ch[j] - 16]) { + if (qqq.ch[i] < 16 && qqq.ch[j] >= 16 && qqqGainValid[qqq.id[i]][qqq.ch[i]][qqq.ch[j] - 16]) + { chWedge = qqq.ch[i]; eWedge = qqq.e[i] * qqqGain[qqq.id[i]][qqq.ch[i]][qqq.ch[j] - 16]; chRing = qqq.ch[j] - 16; @@ -380,7 +395,8 @@ Bool_t MakeVertex::Process(Long64_t entry) tRing = static_cast(qqq.t[j]); tWedge = static_cast(qqq.t[i]); } - else if (qqq.ch[j] < 16 && qqq.ch[i] >= 16 && qqqGainValid[qqq.id[j]][qqq.ch[j]][qqq.ch[i] - 16]) { + else if (qqq.ch[j] < 16 && qqq.ch[i] >= 16 && qqqGainValid[qqq.id[j]][qqq.ch[j]][qqq.ch[i] - 16]) + { chWedge = qqq.ch[j]; eWedge = qqq.e[j] * qqqGain[qqq.id[j]][qqq.ch[j]][qqq.ch[i] - 16]; chRing = qqq.ch[i] - 16; @@ -395,24 +411,27 @@ Bool_t MakeVertex::Process(Long64_t entry) plotter->Fill2D("RingE_vs_Index", 16 * 4, 0, 16 * 4, 1000, 0, 16000, chRing + qqq.id[i] * 16, eRing, "hRawQQQ"); plotter->Fill2D("WedgeE_vs_Index", 16 * 4, 0, 16 * 4, 1000, 0, 16000, chWedge + qqq.id[i] * 16, eWedge, "hRawQQQ"); - if (qqqCalibValid[qqq.id[i]][chWedge][chRing]) { + if (qqqCalibValid[qqq.id[i]][chWedge][chRing]) + { eWedgeMeV = eWedge * qqqCalib[qqq.id[i]][chWedge][chRing] / 1000; eRingMeV = eRing * qqqCalib[qqq.id[i]][chWedge][chRing] / 1000; - if(eRingMeV/eWedgeMeV > 3.0 || eRingMeV/eWedgeMeV<1.0/3.0) continue; - //if(eRingMeV<4.0 || eWedgeMeV<4.0) continue; - + if (eRingMeV / eWedgeMeV > 3.0 || eRingMeV / eWedgeMeV < 1.0 / 3.0) + continue; + // if(eRingMeV<4.0 || eWedgeMeV<4.0) continue; + double theta = -TMath::Pi() / 2 + 2 * TMath::Pi() / 16 / 4. * (qqq.id[i] * 16 + chWedge + 0.5); double rho = 50. + (50. / 16.) * (chRing + 0.5); //"?" - //z used to be 75+30+23=128 - //we found a 12mm shift towards the vertex later --> 116 - Event qqqevent(TVector3(rho*TMath::Cos(theta),rho*TMath::Sin(theta),qqqpos), eRingMeV, eWedgeMeV, tRing, tWedge,chRing+qqq.id[i]*16, chWedge+qqq.id[i]*16); - Event qqqeventr(TVector3(rho*TMath::Cos(theta),rho*TMath::Sin(theta),qqqpos), eRing, eWedge, tRing, tWedge,chRing+qqq.id[i]*16, chWedge+qqq.id[i]*16); + // z used to be 75+30+23=128 + // we found a 12mm shift towards the vertex later --> 116 + Event qqqevent(TVector3(rho * TMath::Cos(theta), rho * TMath::Sin(theta), qqqpos), eRingMeV, eWedgeMeV, tRing, tWedge, chRing + qqq.id[i] * 16, chWedge + qqq.id[i] * 16); + Event qqqeventr(TVector3(rho * TMath::Cos(theta), rho * TMath::Sin(theta), qqqpos), eRing, eWedge, tRing, tWedge, chRing + qqq.id[i] * 16, chWedge + qqq.id[i] * 16); QQQ_Events.push_back(qqqevent); QQQ_Events_Raw.push_back(qqqeventr); plotter->Fill2D("QQQCartesianPlot", 200, -100, 100, 200, -100, 100, rho * TMath::Cos(theta), rho * TMath::Sin(theta), "hCalQQQ"); plotter->Fill2D("QQQCartesianPlot" + std::to_string(qqq.id[i]), 200, -100, 100, 200, -100, 100, rho * TMath::Cos(theta), rho * TMath::Sin(theta), "hCalQQQ"); - if (PCQQQTimeCut) { + if (PCQQQTimeCut) + { plotter->Fill2D("PC_XY_Projection_QQQ_TimeCut" + std::to_string(qqq.id[i]), 400, -100, 100, 400, -100, 100, rho * TMath::Cos(theta), rho * TMath::Sin(theta), "hPCQQQ"); } plotter->Fill2D("PC_XY_Projection_QQQ" + std::to_string(qqq.id[i]), 400, -100, 100, 400, -100, 100, rho * TMath::Cos(theta), rho * TMath::Sin(theta), "hPCQQQ"); @@ -438,19 +457,21 @@ Bool_t MakeVertex::Process(Long64_t entry) plotter->Fill2D("CalibratedQQQEvsPCE_R", 1000, 0, 10, 2000, 0, 30000, eRingMeV, pc.e[k], "hPCQQQ"); plotter->Fill2D("CalibratedQQQEvsPCE_W", 1000, 0, 10, 2000, 0, 30000, eWedgeMeV, pc.e[k], "hPCQQQ"); if (tRing - static_cast(pc.t[k]) < -150) // proton tests, 27Al - //if (tRing - static_cast(pc.t[k]) < -150 && tRing - static_cast(pc.t[k]) > -450) // 27Al - //if (tRing - static_cast(pc.t[k]) < -70 && tRing - static_cast(pc.t[k]) > -150) // 17F + // if (tRing - static_cast(pc.t[k]) < -150 && tRing - static_cast(pc.t[k]) > -450) // 27Al + // if (tRing - static_cast(pc.t[k]) < -70 && tRing - static_cast(pc.t[k]) > -150) // 17F { PCQQQTimeCut = true; } } - if (pc.index[k] >= 24 && pc.e[k] > 50) { + if (pc.index[k] >= 24 && pc.e[k] > 50) + { plotter->Fill2D("Timing_Difference_QQQ_PC_Cathode", 500, -2000, 2000, 16, 0, 16, tRing - static_cast(pc.t[k]), chRing, "hTiming"); } - } //end of pc k loop + } // end of pc k loop - if (!HitNonZero) { + if (!HitNonZero) + { double theta = -TMath::Pi() / 2 + 2 * TMath::Pi() / 16 / 4. * (qqq.id[i] * 16 + chWedge + 0.5); double rho = 50. + (50. / 16.) * (chRing + 0.5); //"?" double x = rho * TMath::Cos(theta); @@ -460,22 +481,21 @@ Bool_t MakeVertex::Process(Long64_t entry) qqqtimestamp = tRing; HitNonZero = true; } - } // if j==i - } //j loop end - } //i loop end + } // if j==i + } // j loop end + } // i loop end plotter->Fill1D("QQQ_Multiplicity", 10, 0, 10, qqqCount, "hRawQQQ"); /*if(QQQ_Events.size()>=1) { std::cout<< " ---->" << std::endl; - for(auto qe: QQQ_Events) { - std::cout << qe.ch1/16 << " " <> WireEvent; //this stores nearest neighbour wire events, or a 'cluster' - WireEvent aWireEvents, cWireEvents; //naming for book keeping + typedef std::unordered_map> WireEvent; // this stores nearest neighbour wire events, or a 'cluster' + WireEvent aWireEvents, cWireEvents; // naming for book keeping aWireEvents.clear(); aWireEvents.reserve(24); @@ -489,8 +509,10 @@ Bool_t MakeVertex::Process(Long64_t entry) if (pc.e[i] > 50) { plotter->Fill2D("PC_Index_Vs_Energy", 48, 0, 48, 2000, 0, 30000, pc.index[i], static_cast(pc.e[i]), "hRawPC"); - } else { - continue; + } + else + { + continue; } if (pc.index[i] < 48) @@ -503,13 +525,13 @@ Bool_t MakeVertex::Process(Long64_t entry) { anodeT = static_cast(pc.t[i]); anodeIndex = pc.index[i]; - aWireEvents[pc.index[i]] = std::tuple(pc.index[i],pc.e[i],static_cast(pc.t[i])); + aWireEvents[pc.index[i]] = std::tuple(pc.index[i], pc.e[i], static_cast(pc.t[i])); } else { cathodeT = static_cast(pc.t[i]); cathodeIndex = pc.index[i] - 24; - cWireEvents[pc.index[i]-24] = std::tuple(pc.index[i]-24,pc.e[i],static_cast(pc.t[i])); + cWireEvents[pc.index[i] - 24] = std::tuple(pc.index[i] - 24, pc.e[i], static_cast(pc.t[i])); } if (anodeT != -99999 && cathodeT != 99999) @@ -552,83 +574,121 @@ Bool_t MakeVertex::Process(Long64_t entry) double aEMax = 0; int aIDMax = 0; - - - for (int i = 0; i < pc.multi; i++) { + for (int i = 0; i < pc.multi; i++) + { // if (pc.e[i] > 100) { - if (pc.index[i] < 24) { + if (pc.index[i] < 24) + { anodeHits.push_back(std::pair(pc.index[i], pc.e[i])); } - else if (pc.index[i] >= 24) { + else if (pc.index[i] >= 24) + { cathodeHits.push_back(std::pair(pc.index[i] - 24, pc.e[i])); } } } - std::sort(anodeHits.begin(),anodeHits.end(),[](std::pair a, std::pair b){ return a.first < b.first;}); - std::sort(cathodeHits.begin(),cathodeHits.end(),[](std::pair a, std::pair b){ return a.first < b.first;}); + std::sort(anodeHits.begin(), anodeHits.end(), [](std::pair a, std::pair b) + { return a.first < b.first; }); + std::sort(cathodeHits.begin(), cathodeHits.end(), [](std::pair a, std::pair b) + { return a.first < b.first; }); - //clusters = collection of (collection of wires) where each wire is (index, energy, timestamp) - std::vector>> aClusters = pwinstance.Make_Clusters(aWireEvents); - std::vector>> cClusters = pwinstance.Make_Clusters(cWireEvents); - - std::vector> sumE_AC; - for(auto aCluster: aClusters) { - for(auto cCluster: cClusters) { - if(aCluster.size()<=1 && cCluster.size()<=1) continue; - auto [crossover,alpha,apSumE,cpSumE,apMaxE,cpMaxE,apTSMaxE,cpTSMaxE] = pwinstance.FindCrossoverProperties(aCluster, cCluster); - if(alpha!=9999999 && apSumE!=-1) { - //Event PCEvent(crossover,apMaxE,cpMaxE,apTSMaxE,cpTSMaxE); - //Event PCEvent(crossover,apSumE,cpSumE,apTSMaxE,cpTSMaxE); - Event PCEvent(crossover,apSumE,cpMaxE,apTSMaxE,cpTSMaxE); //run12 shows cathode-max and anode-sum provide best dE signals. - //std::cout << apSumE << " " << crossover.Perp() << " " << apMaxE << " " << apTSMaxE << std::endl; - PC_Events.push_back(PCEvent); - sumE_AC.push_back(std::pair(apSumE,cpSumE)); - } - } - } - if(QQQ_Events.size() && PC_Events.size()) - plotter->Fill2D("PCEv_vs_QQQEv",20,0,20,20,0,20,QQQ_Events.size(),PC_Events.size()); + // clusters = collection of (collection of wires) where each wire is (index, energy, timestamp) + std::vector>> aClusters = pwinstance.Make_Clusters(aWireEvents); + std::vector>> cClusters = pwinstance.Make_Clusters(cWireEvents); - for(auto pcevent:PC_Events) { - for(auto sx3event:sx3Events) { - plotter->Fill1D("dt_pcA_sx3B"+std::to_string(sx3event.ch2),640,-2000,2000,sx3event.Time1 - pcevent.Time1); - plotter->Fill1D("dt_pcC_sx3B"+std::to_string(sx3event.ch2),640,-2000,2000,sx3event.Time1 - pcevent.Time2); - plotter->Fill2D("dE_E_Anodesx3B",400,0,10,800,0,40000,sx3event.Energy1*0.001,pcevent.Energy1); + std::vector> sumE_AC; + for (auto aCluster : aClusters) + { + for (auto cCluster : cClusters) + { + // if (aCluster.size() <= 1 && cCluster.size() <= 1) + // continue; + auto [crossover, alpha, apSumE, cpSumE, apMaxE, cpMaxE, apTSMaxE, cpTSMaxE] = pwinstance.FindCrossoverProperties(aCluster, cCluster); + if (alpha != 9999999 && apSumE != -1) + { + // Event PCEvent(crossover,apMaxE,cpMaxE,apTSMaxE,cpTSMaxE); + // Event PCEvent(crossover,apSumE,cpSumE,apTSMaxE,cpTSMaxE); + Event PCEvent(crossover, apSumE, cpMaxE, apTSMaxE, cpTSMaxE); // run12 shows cathode-max and anode-sum provide best dE signals. + // std::cout << apSumE << " " << crossover.Perp() << " " << apMaxE << " " << apTSMaxE << std::endl; + PC_Events.push_back(PCEvent); + sumE_AC.push_back(std::pair(apSumE, cpSumE)); + } + } + } + if (QQQ_Events.size() && PC_Events.size()) + plotter->Fill2D("PCEv_vs_QQQEv", 20, 0, 20, 20, 0, 20, QQQ_Events.size(), PC_Events.size()); - plotter->Fill2D("dE_E_Cathodesx3B",400,0,10,800,0,10000,sx3event.Energy1*0.001,pcevent.Energy2); - double sx3z = sx3event.pos.Z()+(75.0/2.0)-3.0; //w.r.t target origin at 90 for run12 - double sx3rho = 88.0;//approximate barrel radius - double sx3theta = TMath::ATan2(sx3rho,sx3z-vertexpos); - double pczguess = pcrad/TMath::Tan(sx3theta) + vertexpos; - plotter->Fill2D("pcz_vs_sx3pczguess",300,0,200,150,0,200,pczguess,pcevent.pos.Z()); - plotter->Fill2D("pcz_vs_sx3pczguess"+std::to_string(sx3event.ch2),300,0,200,150,0,200,pczguess,pcevent.pos.Z()); - plotter->Fill2D("pcz_vs_sx3z",300,0,200,150,0,200,sx3z,pcevent.pos.Z()); - } - } + for (auto pcevent : PC_Events) + { + for (auto sx3event : sx3Events) + { + plotter->Fill1D("dt_pcA_sx3B" + std::to_string(sx3event.ch2), 640, -2000, 2000, sx3event.Time1 - pcevent.Time1); + plotter->Fill1D("dt_pcC_sx3B" + std::to_string(sx3event.ch2), 640, -2000, 2000, sx3event.Time1 - pcevent.Time2); + plotter->Fill2D("dE_E_Anodesx3B", 400, 0, 10, 800, 0, 40000, sx3event.Energy1 * 0.001, pcevent.Energy1); - for(auto pcevent: PC_Events) { - for(auto qqqevent: QQQ_Events) { - plotter->Fill1D("dt_pcA_qqqR",640,-2000,2000,qqqevent.Time1 - pcevent.Time1); - plotter->Fill1D("dt_pcC_qqqW",640,-2000,2000,qqqevent.Time2 - pcevent.Time2); - plotter->Fill2D("dE_E_AnodeQQQR",400,0,10,800,0,40000,qqqevent.Energy1,pcevent.Energy1); - plotter->Fill2D("dE_E_CathodeQQQR",400,0,10,800,0,10000,qqqevent.Energy2,pcevent.Energy2); - double sinTheta = TMath::Sin((qqqevent.pos - TVector3(0,0,vertexpos)).Theta())/TMath::Sin((TVector3(51.5,0,qqqpos) - TVector3(0,0,vertexpos)).Theta()); - plotter->Fill2D("dE2_E_AnodeQQQR",400,0,10,800,0,40000,qqqevent.Energy1,pcevent.Energy1*sinTheta); - plotter->Fill2D("dE2_E_CathodeQQQR",400,0,10,800,0,10000,qqqevent.Energy2,pcevent.Energy2*sinTheta); + plotter->Fill2D("dE_E_Cathodesx3B", 400, 0, 10, 800, 0, 10000, sx3event.Energy1 * 0.001, pcevent.Energy2); + double sx3z = sx3event.pos.Z() + (75.0 / 2.0) - 3.0; // w.r.t target origin at 90 for run12 + double sx3rho = 88.0; // approximate barrel radius + double sx3theta = TMath::ATan2(sx3rho, sx3z - vertexpos); + double pczguess = pcrad / TMath::Tan(sx3theta) + vertexpos; + plotter->Fill2D("pcz_vs_sx3pczguess", 300, 0, 200, 150, 0, 200, pczguess, pcevent.pos.Z()); + plotter->Fill2D("pcz_vs_sx3pczguess" + std::to_string(sx3event.ch2), 300, 0, 200, 150, 0, 200, pczguess, pcevent.pos.Z()); + plotter->Fill2D("pcz_vs_sx3z", 300, 0, 200, 150, 0, 200, sx3z, pcevent.pos.Z()); + } + } - if(qqqevent.pos.Phi() <= pcevent.pos.Phi()+TMath::Pi()/4. && qqqevent.pos.Phi() >= pcevent.pos.Phi()-TMath::Pi()/4.) { - plotter->Fill1D("PCZ",800,-200,200,pcevent.pos.Z(),"phicut"); - double pcz_guess = pcrad/TMath::Tan((qqqevent.pos-TVector3(0,0,vertexpos)).Theta()) + vertexpos; //this is ideally kept to be all QQQ+userinput for calibration of pcz - plotter->Fill2D("pczguess_vs_pc",300,0,200,150,0,200,pcz_guess,pcevent.pos.Z(),"phicut"); - plotter->Fill2D("pczguess_vs_pc_phi="+std::to_string(qqqevent.pos.Phi()*180./M_PI),300,0,200,150,0,200,pcz_guess,pcevent.pos.Z(),"phicut"); - //plotter->Fill1D("PCZ",800,-200,200,pcevent.pos.Z(),"phicut"); - } - } - } - //HALFTIME! Can stop here in future versions - //return kTRUE; + for (auto pcevent : PC_Events) + { + for (auto qqqevent : QQQ_Events) + { + plotter->Fill1D("dt_pcA_qqqR", 640, -2000, 2000, qqqevent.Time1 - pcevent.Time1); + plotter->Fill1D("dt_pcC_qqqW", 640, -2000, 2000, qqqevent.Time2 - pcevent.Time2); + plotter->Fill2D("dE_E_AnodeQQQR", 400, 0, 10, 800, 0, 40000, qqqevent.Energy1, pcevent.Energy1); + plotter->Fill2D("dE_E_CathodeQQQR", 400, 0, 10, 800, 0, 10000, qqqevent.Energy2, pcevent.Energy2); + double sinTheta = TMath::Sin((qqqevent.pos - TVector3(0, 0, vertexpos)).Theta()) / TMath::Sin((TVector3(51.5, 0, qqqpos) - TVector3(0, 0, vertexpos)).Theta()); + plotter->Fill2D("dE2_E_AnodeQQQR", 400, 0, 10, 800, 0, 40000, qqqevent.Energy1, pcevent.Energy1 * sinTheta); + plotter->Fill2D("dE2_E_CathodeQQQR", 400, 0, 10, 800, 0, 10000, qqqevent.Energy2, pcevent.Energy2 * sinTheta); + + if (qqqevent.pos.Phi() <= pcevent.pos.Phi() + TMath::Pi() / 4. && qqqevent.pos.Phi() >= pcevent.pos.Phi() - TMath::Pi() / 4.) + { + plotter->Fill1D("PCZ", 800, -200, 200, pcevent.pos.Z(), "phicut"); + double pcz_guess = pcrad / TMath::Tan((qqqevent.pos - TVector3(0, 0, vertexpos)).Theta()) + vertexpos; // this is ideally kept to be all QQQ+userinput for calibration of pcz + plotter->Fill2D("pczguess_vs_pc", 300, 0, 200, 150, 0, 200, pcz_guess, pcevent.pos.Z(), "phicut"); + plotter->Fill2D("pczguess_vs_pc_phi=" + std::to_string(qqqevent.pos.Phi() * 180. / M_PI), 300, 0, 200, 150, 0, 200, pcz_guess, pcevent.pos.Z(), "phicut"); + // plotter->Fill1D("PCZ",800,-200,200,pcevent.pos.Z(),"phicut"); + } + + if (cClusters.size() == 1) + { + plotter->Fill1D("pcz_1cCluster", 600, -300, 300, pcevent.pos.Z(), "hPCzQQQ"); + } + else if (cClusters.size() == 2) + { + plotter->Fill1D("pcz_2cCluster", 600, -300, 300, pcevent.pos.Z(), "hPCzQQQ"); + } + else if (cClusters.size() >= 3) + { + plotter->Fill1D("pcz_ncCluster", 600, -300, 300, pcevent.pos.Z(), "hPCzQQQ"); + } + + if (aClusters.size() == 1) + { + plotter->Fill1D("pcz_1aCluster", 600, -300, 300, pcevent.pos.Z(), "hPCzQQQ"); + } + else if (aClusters.size() == 2) + { + plotter->Fill1D("pcz_2aCluster", 600, -300, 300, pcevent.pos.Z(), "hPCzQQQ"); + } + else if (aClusters.size() >= 3) + { + plotter->Fill1D("pcz_naCluster", 600, -300, 300, pcevent.pos.Z(), "hPCzQQQ"); + } + } + } + // HALFTIME! Can stop here in future versions + // return kTRUE; if (anodeHits.size() >= 1 && cathodeHits.size() >= 1) { @@ -670,7 +730,7 @@ Bool_t MakeVertex::Process(Long64_t entry) } } - TVector3 anodeIntersection,vector_closest_to_z; + TVector3 anodeIntersection, vector_closest_to_z; anodeIntersection.Clear(); vector_closest_to_z.Clear(); if (corrcatMax.size() > 0) @@ -696,26 +756,28 @@ Bool_t MakeVertex::Process(Long64_t entry) } bool PCQQQPhiCut = false; // flip the algorithm for cathode 1 multi anode events - if ((hitPos.Phi() > (anodeIntersection.Phi() - TMath::PiOver4())) && (hitPos.Phi() < (anodeIntersection.Phi() + TMath::PiOver4()))) { + if ((hitPos.Phi() > (anodeIntersection.Phi() - TMath::PiOver4())) && (hitPos.Phi() < (anodeIntersection.Phi() + TMath::PiOver4()))) + { PCQQQPhiCut = true; } for (double Tz = 60; Tz <= 100; Tz += 1.0) { - TVector3 TargetPos(0, 0, Tz); - if(PCQQQPhiCut && anodeIntersection.Perp()>0 && anodeIntersection.Z()!=0 && cathodeHits.size()>=2) { - plotter->Fill2D("Inttheta_vs_QQQtheta_TC" + std::to_string(PCQQQTimeCut) + "_TZ" + std::to_string(Tz), 400, 0, 180, 90, 0, 90, (anodeIntersection - TargetPos).Theta() * 180. / TMath::Pi(), (hitPos - TargetPos).Theta() * 180. / TMath::Pi(), "TPosVariation"); - //plotter->Fill2D("R_ratio_to_Z_ratio" + std::to_string(PCQQQTimeCut) + "_TZ" + std::to_string(Tz), 100, -2, 2, 100, -2, 2, (anodeIntersection - TargetPos).Z()/(hitPos-TargetPos).Z(), ((anodeIntersection - TargetPos).Perp()+2.5)/(hitPos-TargetPos).Perp(), "TPosVariation"); - } + TVector3 TargetPos(0, 0, Tz); + if (PCQQQPhiCut && anodeIntersection.Perp() > 0 && anodeIntersection.Z() != 0 && cathodeHits.size() >= 2) + { + plotter->Fill2D("Inttheta_vs_QQQtheta_TC" + std::to_string(PCQQQTimeCut) + "_TZ" + std::to_string(Tz), 400, 0, 180, 90, 0, 90, (anodeIntersection - TargetPos).Theta() * 180. / TMath::Pi(), (hitPos - TargetPos).Theta() * 180. / TMath::Pi(), "TPosVariation"); + // plotter->Fill2D("R_ratio_to_Z_ratio" + std::to_string(PCQQQTimeCut) + "_TZ" + std::to_string(Tz), 100, -2, 2, 100, -2, 2, (anodeIntersection - TargetPos).Z()/(hitPos-TargetPos).Z(), ((anodeIntersection - TargetPos).Perp()+2.5)/(hitPos-TargetPos).Perp(), "TPosVariation"); + } } - if (anodeIntersection.Z() != 0 && anodeIntersection.Perp()>0 && HitNonZero) + if (anodeIntersection.Z() != 0 && anodeIntersection.Perp() > 0 && HitNonZero) { plotter->Fill1D("PC_Z_Projection", 600, -300, 300, anodeIntersection.Z(), "hPCzQQQ"); plotter->Fill2D("Z_Proj_VsDelTime", 600, -300, 300, 200, -2000, 2000, anodeIntersection.Z(), anodeT - cathodeT, "hPCzQQQ"); plotter->Fill2D("IntPhi_vs_QQQphi", 100, -200, 200, 80, -200, 200, anodeIntersection.Phi() * 180. / TMath::Pi(), hitPos.Phi() * 180. / TMath::Pi(), "hPCQQQ"); - //plotter->Fill2D("Inttheta_vs_QQQtheta", 90, 0, 180, 20, 0, 45, anodeIntersection.Theta() * 180. / TMath::Pi(), hitPos.Theta() * 180. / TMath::Pi(), "hPCQQQ"); - //plotter->Fill2D("Inttheta_vs_QQQtheta_TC" + std::to_string(PCQQQTimeCut)+ "_PC"+std::to_string(PCQQQPhiCut), 90, 0, 180, 20, 0, 45, anodeIntersection.Theta() * 180. / TMath::Pi(), hitPos.Theta() * 180. / TMath::Pi(), "hPCQQQ"); + // plotter->Fill2D("Inttheta_vs_QQQtheta", 90, 0, 180, 20, 0, 45, anodeIntersection.Theta() * 180. / TMath::Pi(), hitPos.Theta() * 180. / TMath::Pi(), "hPCQQQ"); + // plotter->Fill2D("Inttheta_vs_QQQtheta_TC" + std::to_string(PCQQQTimeCut)+ "_PC"+std::to_string(PCQQQPhiCut), 90, 0, 180, 20, 0, 45, anodeIntersection.Theta() * 180. / TMath::Pi(), hitPos.Theta() * 180. / TMath::Pi(), "hPCQQQ"); plotter->Fill2D("IntPhi_vs_QQQphi_TC" + std::to_string(PCQQQTimeCut) + "PhiC" + std::to_string(PCQQQPhiCut), 100, -200, 200, 80, -200, 200, anodeIntersection.Phi() * 180. / TMath::Pi(), hitPos.Phi() * 180. / TMath::Pi(), "hPCQQQ"); } if (anodeIntersection.Z() != 0 && cathodeHits.size() >= 2) @@ -796,42 +858,46 @@ Bool_t MakeVertex::Process(Long64_t entry) { pw_contr.CalTrack2(hitPos, anodeIntersection); plotter->Fill1D("VertexRecon", 600, -1300, 1300, pw_contr.GetZ0()); - plotter->Fill1D("VertexRecon_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 600, -1300, 1300, pw_contr.GetZ0()); + plotter->Fill1D("VertexRecon_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 600, -1300, 1300, pw_contr.GetZ0()); if (cathodeHits.size() == 2) - plotter->Fill1D("VertexRecon_2c_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 600, -1300, 1300, pw_contr.GetZ0()); + plotter->Fill1D("VertexRecon_2c_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 600, -1300, 1300, pw_contr.GetZ0()); TVector3 x2(anodeIntersection), x1(hitPos); - TVector3 v = x2-x1; - double t_minimum = -1.0*(x1.X()*v.X()+x1.Y()*v.Y())/(v.X()*v.X()+v.Y()*v.Y()); - vector_closest_to_z = x1 + t_minimum*v; - - plotter->Fill1D("VertexRecon_Z_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 600, -1300, 1300, vector_closest_to_z.Z() ,"customVertex"); - if(vector_closest_to_z.Perp() < 20) { - plotter->Fill1D("VertexRecon_RadialCut_Z_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 600, -1300, 1300, vector_closest_to_z.Z() ,"customVertex"); - } - - plotter->Fill2D("VertexRecon_XY_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 100, -100, 100, 100,-100,100, vector_closest_to_z.X(), vector_closest_to_z.Y() ,"customVertex"); - if(cathodeHits.size()==2) { - plotter->Fill1D("VertexRecon2C_Z_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 600, -1300, 1300, vector_closest_to_z.Z() ,"customVertex"); - if(vector_closest_to_z.Perp() < 20) { - plotter->Fill1D("VertexRecon2C_RadialCut_Z_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 600, -1300, 1300, vector_closest_to_z.Z() ,"customVertex"); - } - plotter->Fill2D("VertexRecon2C_XY_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 100, -100, 100, 100,-100,100, vector_closest_to_z.X(), vector_closest_to_z.Y() ,"customVertex"); - plotter->Fill2D("VertexRecon2C_RhoZ_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 100, -100, 100, 600,-1300,1300, vector_closest_to_z.Perp(), vector_closest_to_z.Z() ,"customVertex"); - plotter->Fill2D("VertexRecon2C_Z_vs_QQQE_TC"+std::to_string(PCQQQTimeCut)+"_PhiC"+std::to_string(PCQQQPhiCut), 600, -1300, 1300, 800,0,20000, vector_closest_to_z.Z(), qqqenergy ,"customVertex"); - } - + TVector3 v = x2 - x1; + double t_minimum = -1.0 * (x1.X() * v.X() + x1.Y() * v.Y()) / (v.X() * v.X() + v.Y() * v.Y()); + vector_closest_to_z = x1 + t_minimum * v; + + plotter->Fill1D("VertexRecon_Z_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 600, -1300, 1300, vector_closest_to_z.Z(), "customVertex"); + if (vector_closest_to_z.Perp() < 20) + { + plotter->Fill1D("VertexRecon_RadialCut_Z_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 600, -1300, 1300, vector_closest_to_z.Z(), "customVertex"); + } + + plotter->Fill2D("VertexRecon_XY_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 100, -100, 100, 100, -100, 100, vector_closest_to_z.X(), vector_closest_to_z.Y(), "customVertex"); + if (cathodeHits.size() == 2) + { + plotter->Fill1D("VertexRecon2C_Z_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 600, -1300, 1300, vector_closest_to_z.Z(), "customVertex"); + if (vector_closest_to_z.Perp() < 20) + { + plotter->Fill1D("VertexRecon2C_RadialCut_Z_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 600, -1300, 1300, vector_closest_to_z.Z(), "customVertex"); + } + plotter->Fill2D("VertexRecon2C_XY_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 100, -100, 100, 100, -100, 100, vector_closest_to_z.X(), vector_closest_to_z.Y(), "customVertex"); + plotter->Fill2D("VertexRecon2C_RhoZ_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 100, -100, 100, 600, -1300, 1300, vector_closest_to_z.Perp(), vector_closest_to_z.Z(), "customVertex"); + plotter->Fill2D("VertexRecon2C_Z_vs_QQQE_TC" + std::to_string(PCQQQTimeCut) + "_PhiC" + std::to_string(PCQQQPhiCut), 600, -1300, 1300, 800, 0, 20000, vector_closest_to_z.Z(), qqqenergy, "customVertex"); + } } for (int i = 0; i < qqq.multi; i++) { - if(anodeIntersection.Perp() > 0) { //suppress x,y=0,0 events - if (PCQQQTimeCut) { - plotter->Fill2D("PC_XY_Projection_QQQ_TimeCut" + std::to_string(qqq.id[i]), 400, -100, 100, 400, -100, 100, anodeIntersection.X(), anodeIntersection.Y(), "hPCQQQ"); - } - plotter->Fill2D("PC_XY_Projection_QQQ" + std::to_string(qqq.id[i]), 400, -100, 100, 400, -100, 100, anodeIntersection.X(), anodeIntersection.Y(), "hPCQQQ"); + if (anodeIntersection.Perp() > 0) + { // suppress x,y=0,0 events + if (PCQQQTimeCut) + { + plotter->Fill2D("PC_XY_Projection_QQQ_TimeCut" + std::to_string(qqq.id[i]), 400, -100, 100, 400, -100, 100, anodeIntersection.X(), anodeIntersection.Y(), "hPCQQQ"); + } + plotter->Fill2D("PC_XY_Projection_QQQ" + std::to_string(qqq.id[i]), 400, -100, 100, 400, -100, 100, anodeIntersection.X(), anodeIntersection.Y(), "hPCQQQ"); } for (int j = i + 1; j < qqq.multi; j++) { @@ -891,48 +957,49 @@ Bool_t MakeVertex::Process(Long64_t entry) // plotter->Fill2D("EdE_PC_vs_QQQ_timegate_ls1000"+std::to_string()) plotter->Fill2D("PC_Z_vs_QQQRing_Det" + std::to_string(qqqID), 600, -300, 300, 16, 0, 16, anodeIntersection.Z(), chRing, "hPCQQQ"); - //double theta = -TMath::Pi() / 2 + 2 * TMath::Pi() / 16 / 4. * (qqq.id[i] * 16 + chWedge + 0.5); - //double rho = 50. + 40. / 16. * (chRing + 0.5); + // double theta = -TMath::Pi() / 2 + 2 * TMath::Pi() / 16 / 4. * (qqq.id[i] * 16 + chWedge + 0.5); + // double rho = 50. + 40. / 16. * (chRing + 0.5); for (int k = 0; k < pc.multi; k++) { - if(pc.index[k] >= 24) + if (pc.index[k] >= 24) continue; -// double sinTheta = TMath::Sin((hitPos-vector_closest_to_z).Theta()); - double sinTheta = TMath::Sin((anodeIntersection-TVector3(0,0,90.0)).Theta()); -// double sinTheta = TMath::Sin((anodeIntersection-vector_closest_to_z).Theta()); -// double sinTheta = TMath::Sin((hitPos-TVector3(0,0,30.0)).Theta()); -// double sinTheta = TMath::Sin(hitPos.Theta()); + // double sinTheta = TMath::Sin((hitPos-vector_closest_to_z).Theta()); + double sinTheta = TMath::Sin((anodeIntersection - TVector3(0, 0, 90.0)).Theta()); + // double sinTheta = TMath::Sin((anodeIntersection-vector_closest_to_z).Theta()); + // double sinTheta = TMath::Sin((hitPos-TVector3(0,0,30.0)).Theta()); + // double sinTheta = TMath::Sin(hitPos.Theta()); - if(cathodeHits.size()==2 && PCQQQPhiCut) { - plotter->Fill2D("CalibratedQQQE_RvsCPCE_TC" + std::to_string(PCQQQTimeCut) , 400, 0, 10, 400, 0, 30000, eRingMeV, pc.e[k]*sinTheta, "hPCQQQ"); - plotter->Fill2D("CalibratedQQQE_WvsCPCE_TC" + std::to_string(PCQQQTimeCut) , 400, 0, 10, 400, 0, 30000, eWedgeMeV, pc.e[k]*sinTheta, "hPCQQQ"); - plotter->Fill2D("CalibratedQQQE_RvsPCE_TC" + std::to_string(PCQQQTimeCut) , 400, 0, 10, 400, 0, 30000, eRingMeV, pc.e[k], "hPCQQQ"); - plotter->Fill2D("CalibratedQQQE_WvsPCE_TC" + std::to_string(PCQQQTimeCut) , 400, 0, 10, 400, 0, 30000, eWedgeMeV, pc.e[k], "hPCQQQ"); - plotter->Fill2D("PCQQQ_dTimevsdPhi", 200, -2000, 2000, 80, -200, 200, tRing - static_cast(pc.t[k]), (hitPos.Phi()-anodeIntersection.Phi()) * 180. / TMath::Pi(), "hTiming"); + if (cathodeHits.size() == 2 && PCQQQPhiCut) + { + plotter->Fill2D("CalibratedQQQE_RvsCPCE_TC" + std::to_string(PCQQQTimeCut), 400, 0, 10, 400, 0, 30000, eRingMeV, pc.e[k] * sinTheta, "hPCQQQ"); + plotter->Fill2D("CalibratedQQQE_WvsCPCE_TC" + std::to_string(PCQQQTimeCut), 400, 0, 10, 400, 0, 30000, eWedgeMeV, pc.e[k] * sinTheta, "hPCQQQ"); + plotter->Fill2D("CalibratedQQQE_RvsPCE_TC" + std::to_string(PCQQQTimeCut), 400, 0, 10, 400, 0, 30000, eRingMeV, pc.e[k], "hPCQQQ"); + plotter->Fill2D("CalibratedQQQE_WvsPCE_TC" + std::to_string(PCQQQTimeCut), 400, 0, 10, 400, 0, 30000, eWedgeMeV, pc.e[k], "hPCQQQ"); + plotter->Fill2D("PCQQQ_dTimevsdPhi", 200, -2000, 2000, 80, -200, 200, tRing - static_cast(pc.t[k]), (hitPos.Phi() - anodeIntersection.Phi()) * 180. / TMath::Pi(), "hTiming"); } - } - }///qqq i==j case end - } //j loop end + } /// qqq i==j case end + } // j loop end } // qqq i loop end - - TVector3 guessVertex(0,0,90.); //for run12, subtract anodeIntersection.Z() by ~74.0 seems to work - //rho=40.0 mm is halfway between the cathodes(rho=42) and anodes(rho=37) - double pcz_guess = 42.0/TMath::Tan((hitPos-guessVertex).Theta()) + guessVertex.Z(); //this is ideally kept to be all QQQ+userinput for calibration of pcz - if(PCQQQTimeCut && PCQQQPhiCut && hitPos.Perp()>0 && anodeIntersection.Perp()>0 && cathodeHits.size()>=2) { - plotter->Fill2D("pczguess_vs_qqqE",100,0,200,800,0,20,pcz_guess,qqqenergy,"pczguess"); - double pczoffset=30.0; - //plotter->Fill2D("pczguess_vs_pcz_rad="+std::to_string(hitPos.Perp()),100,0,200,150,0,200,pcz_guess,anodeIntersection.Z(),"pczguess"); //entirely qqq-derived position vs entirely PC derived position - plotter->Fill2D("pczguess_vs_pcz_phi="+std::to_string(hitPos.Phi()*180./M_PI),100,0,200,150,0,200,pcz_guess,anodeIntersection.Z()+pczoffset,"pczguess"); //entirely qqq-derived position vs entirely PC derived position - plotter->Fill2D("pczguess_vs_pcz",100,0,200,150,0,200,pcz_guess,anodeIntersection.Z()+pczoffset); - plotter->Fill2D("pcz_vs_pcPhi_rad="+std::to_string(hitPos.Perp()),360,0,360,150,0,200,anodeIntersection.Phi()*180./M_PI,anodeIntersection.Z()+pczoffset,"pczguess"); + + TVector3 guessVertex(0, 0, 90.); // for run12, subtract anodeIntersection.Z() by ~74.0 seems to work + // rho=40.0 mm is halfway between the cathodes(rho=42) and anodes(rho=37) + double pcz_guess = 42.0 / TMath::Tan((hitPos - guessVertex).Theta()) + guessVertex.Z(); // this is ideally kept to be all QQQ+userinput for calibration of pcz + if (PCQQQTimeCut && PCQQQPhiCut && hitPos.Perp() > 0 && anodeIntersection.Perp() > 0 && cathodeHits.size() >= 2) + { + plotter->Fill2D("pczguess_vs_qqqE", 100, 0, 200, 800, 0, 20, pcz_guess, qqqenergy, "pczguess"); + double pczoffset = 30.0; + // plotter->Fill2D("pczguess_vs_pcz_rad="+std::to_string(hitPos.Perp()),100,0,200,150,0,200,pcz_guess,anodeIntersection.Z(),"pczguess"); //entirely qqq-derived position vs entirely PC derived position + plotter->Fill2D("pczguess_vs_pcz_phi=" + std::to_string(hitPos.Phi() * 180. / M_PI), 100, 0, 200, 150, 0, 200, pcz_guess, anodeIntersection.Z() + pczoffset, "pczguess"); // entirely qqq-derived position vs entirely PC derived position + plotter->Fill2D("pczguess_vs_pcz", 100, 0, 200, 150, 0, 200, pcz_guess, anodeIntersection.Z() + pczoffset); + plotter->Fill2D("pcz_vs_pcPhi_rad=" + std::to_string(hitPos.Perp()), 360, 0, 360, 150, 0, 200, anodeIntersection.Phi() * 180. / M_PI, anodeIntersection.Z() + pczoffset, "pczguess"); } for (int i = 0; i < sx3.multi; i++) { // plotting sx3 strip hits vs anode phi - if (sx3.ch[i] < 8 && anodeIntersection.Perp()>0) + if (sx3.ch[i] < 8 && anodeIntersection.Perp() > 0) plotter->Fill2D("PCPhi_vs_SX3Strip", 100, -200, 200, 8 * 24, 0, 8 * 24, anodeIntersection.Phi() * 180. / TMath::Pi(), sx3.id[i] * 8 + sx3.ch[i]); }