diff --git a/analysis/fillChain.cxx b/analysis/fillChain.cxx index 2f5094eae729c25b79cb9e116bf4292f143ac377..06a8ece4ac5482ddbc12c6f1646fe4cf8718ba1b 100644 --- a/analysis/fillChain.cxx +++ b/analysis/fillChain.cxx @@ -8,19 +8,52 @@ using std::cout; using std::endl; -void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { +Int_t GetClustersMWPC(unsigned short n, unsigned short *x) +{ + Int_t clusters = 1; + + if(n == 0) clusters = 0.; + else if(n > 1 && n<=32) + { + for(Int_t i = 1; i < n; ++i) + { +// if (x[i]>50 || x[i-1]>50) continue; + if( (abs(x[i] - x[i-1]) > 1) && (abs(x[i] - x[i-1]) < 32) ) clusters += 1.; + } + } + return clusters; +} //-------------------------------------------------------------------- + +Double_t GetClusterPositionMWPC(unsigned short n, unsigned short *x, Float_t planeOffset=0.) +{ + Double_t position = -100.; + + + + +// if(n == 0) clusters = 0.; + if(n > 0 && n<=32) + { + position = (x[0]+n/2.+0.5)*1.25-20. + planeOffset; +// cout << n << endl; + } + return position; +} //-------------------------------------------------------------------- + + +void fillTree(const TString beam = "he", Int_t nofile = 0, const Int_t noevents = 0) { TString inFile; - inFile.Form("~/data/exp1804/h5_14_%02d.root", nofile); -// inFile.Form("~/data/exp1804/be10_03_%d0.root", nofile); //files 00,10,...,90 + if (beam.Contains("he")) inFile.Form("~/data/exp1804/h5_14_%02d.root", nofile); + if (beam.Contains("be")) inFile.Form("~/data/exp1804/be10_03_%d0.root", nofile); //files 00,10,...,90 //where 70 is run 01 // 80 is run 02 // 90 is run 05 // inFile.Form("~/data/exp1804/calib/si_20_03.root"); TString outFile; - outFile.Form("~/data/exp1804/h5_14_%02d_calib.root", nofile); -// outFile.Form("~/data/exp1804/be10_03_%d0_calib.root", nofile); //files 00,10,...,60 + if (beam.Contains("he")) outFile.Form("~/data/exp1804/h5_14_%02d_calib.root", nofile); + if (beam.Contains("be")) outFile.Form("~/data/exp1804/be10_03_%d0_calib.root", nofile); //files 00,10,...,60 //where 70 is run 01 // 80 is run 02 // 90 is run 05 @@ -39,9 +72,12 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { TFile *fw = new TFile(outFile, "RECREATE"); TTree *tw = new TTree("cal", "Calibrated information"); + Int_t trigger; + Float_t SQ20E[16]; Float_t SQ20Ecorr[16]; Float_t SQ20EcorrHit[16]; + Float_t SQ20EcorrHitC[16]; Float_t SQ20Esum; Float_t SQ20time[16]; @@ -68,13 +104,26 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { Int_t SQLYmult; Int_t SQ20EcorrMult; Int_t SQ20EcorrHitMult; + Int_t SQ20EcorrHitCMult; + + Bool_t CsI_L_veto; + Int_t expectedThinStrip; + Float_t SQRXE[32]; Float_t SQRXEsum; Int_t SQRXmult; Float_t TOF, dEbeam; + + Float_t tF5; + + //MWPC, wire multiplicity Float_t x1p, y1p, x2p, y2p, xt, yt; + //MWPC, cluster multiplicity + Int_t x1MultC, x2MultC, y1MultC, y2MultC; + Float_t x1c, y1c, x2c, y2c, xtc, ytc; + //left 1 mm position Float_t x1mm, y1mm; @@ -91,6 +140,8 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { TVector3 vNorm(0.,0.,1.); Double_t angleLeft; + tw->Branch("trigger",&trigger,"trigger/I"); + tw->Branch("angleLeft",&angleLeft,"angleLeft/D"); tw->Branch("x1mm",&x1mm,"x1mm/F"); tw->Branch("y1mm",&y1mm,"y1mm/F"); @@ -103,6 +154,7 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { tw->Branch("SQ20E",SQ20E,"SQ20E[16]/F"); tw->Branch("SQ20Ecorr",SQ20Ecorr,"SQ20Ecorr[16]/F"); tw->Branch("SQ20EcorrHit",SQ20EcorrHit,"SQ20EcorrHit[16]/F"); + tw->Branch("SQ20EcorrHitC",SQ20EcorrHitC,"SQ20EcorrHitC[16]/F"); tw->Branch("SQ20Esum",&SQ20Esum,"SQ20Esum/F"); tw->Branch("SQ20time",SQ20time,"SQ20time[16]/F"); @@ -130,6 +182,11 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { tw->Branch("SQLXmult",&SQLXmult,"SQLXmult/I"); tw->Branch("SQLYmult",&SQLYmult,"SQLYmult/I"); tw->Branch("SQ20EcorrMult",&SQ20EcorrMult,"SQ20EcorrMult/I"); + tw->Branch("SQ20EcorrHitCMult",&SQ20EcorrHitCMult,"SQ20EcorrHitCMult/I"); +// tw->Branch("SQ20EcorrMult",&SQ20EcorrMult,"SQ20EcorrMult/I"); + + tw->Branch("CsI_L_veto", &CsI_L_veto, "CsI_L_veto/O"); + tw->Branch("expectedThinStrip",&expectedThinStrip,"expectedThinStrip/I"); tw->Branch("SQRXE",SQRXE,"SQRXE[32]/F"); tw->Branch("SQRXEsum",&SQRXEsum,"SQRXE/F"); @@ -137,6 +194,7 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { tw->Branch("TOF", &TOF, "TOF/F"); tw->Branch("dEbeam", &dEbeam, "dEbeam/F"); + tw->Branch("tF5", &tF5, "tF5/F"); tw->Branch("x1p", &x1p, "x1p/F"); tw->Branch("y1p", &y1p, "y1p/F"); @@ -146,6 +204,19 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { tw->Branch("xt", &xt, "xt/F"); tw->Branch("yt", &yt, "yt/F"); + tw->Branch("x1MultC",&x1MultC,"x1MultC/I"); + tw->Branch("x2MultC",&x2MultC,"x2MultC/I"); + tw->Branch("y1MultC",&y1MultC,"y1MultC/I"); + tw->Branch("y2MultC",&y2MultC,"y2MultC/I"); + + tw->Branch("x1c", &x1c, "x1c/F"); + tw->Branch("y1c", &y1c, "y1c/F"); + tw->Branch("x2c", &x2c, "x2c/F"); + tw->Branch("y2c", &y2c, "y2c/F"); + + tw->Branch("xtc", &xtc, "xtc/F"); + tw->Branch("ytc", &ytc, "ytc/F"); + Long64_t nevents; if (noevents == 0) nevents = tr->GetEntries(); @@ -190,9 +261,9 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { const Int_t kSQL_Y_strips = 16; const Int_t kSQL_20_strips = 16; - const Double_t kSQLY_energy_thr = 1.; - const Double_t kSQLX_energy_thr = 1.; - const Double_t kSQL20_energy_thr = 1.2; +// const Double_t kSQLY_energy_thr = 1.; +// const Double_t kSQLX_energy_thr = 1.; +// const Double_t kSQL20_energy_thr = 1.2; const Double_t kSQ20_norm_thickness = 20.; @@ -210,6 +281,82 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { 6, 4, -4, 13 }; + Int_t thinChThresholds[16] = { + 113, 121, 125, 125, + 125, 132, 130, 115, + 118, 128, 132, 130, + 130, 130, 130, 112 + }; + + Int_t thickLYChThresholds[16] = { + 4000, 130, 82, 110, + 87, 112, 82, 108, + 82, 112, 82, 108, + 80, 110, 90, 105, + }; + + Int_t thickLXChThresholds[32] = { + 60, 45, 45, 45, + 45, 4000, 4000, 45, + 45, 45, 45, 45, + 45, 45, 45, 45, + ///// + 45, 45, 45, 45, + 45, 45, 45, 45, + 45, 45, 45, 45, + 45, 45, 45, 60 + }; + + Int_t minTimeLthin[16] { + 410, 410, 420, 415, + 420, 410, 415, 415, + 420, 420, 420, 420, + 410, 425, 0, 0 + }; + + Int_t maxTimeLthin[16] { + 520, 525, 520, 530, + 520, 530, 525, 525, + 530, 530, 530, 520, + 510, 525, 0, 0 + }; + + Double_t minTimeLX[32] { + 0, 0, 325, 328, + 330, 0, 0, 330, + 326, 328, 326, 335, + 326, 328, 328, 330, + // + 420, 425, 420, 420, + 420, 420, 420, 420, + 430, 425, 420, 420, + 415, 415, 425, 415 + }; + + Double_t maxTimeLX[32] { + 0, 0, 345, 345, + 350, 0, 0, 350, + 345, 345, 345, 350, + 340, 342, 345, 345, + // + 460, 460, 460, 460, + 460, 460, 460, 460, + 460, 460, 460, 460, + 460, 460, 460, 453 + }; + + const UShort_t CsIleftThr = 180; + + const Double_t thinXoffset = 1.; + const Double_t thinYoffset = -1.8; + + const Float_t MWPC1_X_offset = -1.0; + const Float_t MWPC1_Y_offset = -2.1375; + const Float_t MWPC2_X_offset = 0.2; + const Float_t MWPC2_Y_offset = -1.125; + + + fw->cd(); ////////////////////////////////// @@ -222,6 +369,10 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { tr->GetEvent(eventNo); + trigger = revent->trigger; + + if (trigger != 3) continue; + // cout << eventNo << endl; // cout << revent->SQX_L[0] << endl; // SQLXE[0] = revent->SQX_L[0]; @@ -239,16 +390,143 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { SQLYmult = 0; SQ20EcorrMult = 0; SQ20EcorrHitMult = 0; + SQ20EcorrHitCMult = 0; SQ20timeMult = 0; SQLXtimeMult = 0; // SQLXtimeSum = 0; SQLYtimeMult = 0; + CsI_L_veto = kFALSE; + expectedThinStrip = -1; + + + /////////////////////////////////////////////// + //beam diagnostics + /////////////////////////////////////////////// + + /////////////////////////////////////////////// + //TOF + /////////////////////////////////////////////// + + dEbeam = 0.; + TOF = 0.; + if (revent->tF5[0]>0 && revent->tF5[1]>0 && revent->tF5[2]>0 && revent->tF5[3]>0 + && revent->tF3[0]>0 && revent->tF3[1]>0 && revent->tF3[2]>0 && revent->tF3[3]>0 + && revent->F5[0]>0 && revent->F5[1]>0 && revent->F5[2]>0 && revent->F5[3]>0 +// && ( (revent->tF5[0]+revent->tF5[1]+revent->tF5[2]+revent->tF5[3]) - (revent->tF3[0]+revent->tF3[1]+revent->tF3[2]+revent->tF3[3]) )/4.*0.125+89.165<200 +// && ( (tF5[0]+tF5[1]+tF5[2]+tF5[3]) - (tF3[0]+tF3[1]+tF3[2]+tF3[3]) )/4.*0.125+89.165>100 +// && (F5[0]+F5[1]+F5[2]+F5[3])/4. < 2500 + ) { + dEbeam = (revent->F5[0]+revent->F5[1]+revent->F5[2]+revent->F5[3])/4.; + TOF = ( (revent->tF5[0]+revent->tF5[1]+revent->tF5[2]+revent->tF5[3]) - (revent->tF3[0]+revent->tF3[1]+revent->tF3[2]+revent->tF3[3]) )/4.*0.125+89.165; + } + + const Float_t l12 = 546.; + const Float_t lt = 270.; + + /////////////////////////////////////////////// + //MWPC's + /////////////////////////////////////////////// + + + x1p = -100.; + y1p = -100.; + x2p = -100.; + y2p = -100.; + xt = -100.; + yt = -100.; + + x1c = -100.; + y1c = -100.; + x2c = -100.; + y2c = -100.; + xtc = -100.; + ytc = -100.; + + x1MultC = 0; + y1MultC = 0; + x2MultC = 0; + y2MultC = 0; + + + Bool_t flagMWPC = 1; + + tF5 = 0.125*(revent->tF5[0]+revent->tF5[1]+revent->tF5[2]+revent->tF5[3])/4.; +// tF5 = 0.125*revent->tF5[0]; + + if(revent->nx1>10 || revent->ny1>10 || revent->nx2>10 || revent->ny2>10) flagMWPC=0; + if((0.125*revent->tMWPC[0]-tF5)<60 || (0.125*revent->tMWPC[0]-tF5)>77) flagMWPC=0; + if((0.125*revent->tMWPC[1]-tF5)<60 || (0.125*revent->tMWPC[1]-tF5)>80) flagMWPC=0; + if((0.125*revent->tMWPC[2]-tF5)<70 || (0.125*revent->tMWPC[2]-tF5)>90) flagMWPC=0; + if((0.125*revent->tMWPC[3]-tF5)<60 || (0.125*revent->tMWPC[3]-tF5)>80) flagMWPC=0; + + //MWPC: work with wire multiplicity equal to 1 + if (flagMWPC!=0) { + if (revent->x1[0]<1000 && revent->y1[0]<1000 && revent->nx1==1 && revent->ny1==1) { + x1p = (revent->x1[0]+0.5)*1.25-20. + MWPC1_X_offset; + y1p = (revent->y1[0]+0.5)*1.25-20. + MWPC1_Y_offset; + +// (x[0]+n/2.+0.5)*1.25-20.; + } + + if (revent->x2[0]<1000 && revent->y2[0]<1000 && revent->nx2==1 && revent->ny2==1) { + x2p = (revent->x2[0]+0.5)*1.25-20. + MWPC2_X_offset; + y2p = (revent->y2[0]+0.5)*1.25-20. + MWPC2_Y_offset; + } + + if (x1p > -50. && y1p > -50. && x2p > -50. && y2p > -50.) { + xt = x1p - (x1p - x2p)*((l12+lt)/l12); + yt = y1p - (y1p - y2p)*((l12+lt)/l12); +// xt = (l12*x1p + (l12+lt)*(x2p-x1p))/(l12 - (x2p-x1p)*TMath::Tan(TMath::DegToRad()*12.)); +// yt = (y2p-y1p)*(xt-x1p)/(x2p-x1p) + y1p; + } + }//if wire multiplicity == 1 + +// cout << flagMWPC << endl; + + //MWPC: work with cluster multiplicity equal to 2 + if (flagMWPC!=0) { + x1MultC = GetClustersMWPC(revent->nx1, revent->x1); + x2MultC = GetClustersMWPC(revent->nx2, revent->x2); + y1MultC = GetClustersMWPC(revent->ny1, revent->y1); + y2MultC = GetClustersMWPC(revent->ny2, revent->y2); + + if (x1MultC==1 && y1MultC==1) { + x1c = GetClusterPositionMWPC(revent->nx1, revent->x1, MWPC1_X_offset); + y1c = GetClusterPositionMWPC(revent->ny1, revent->y1, MWPC1_Y_offset); + } + + if (x2MultC==1 && y2MultC==1) { + x2c = GetClusterPositionMWPC(revent->nx2, revent->x2, MWPC2_X_offset); + y2c = GetClusterPositionMWPC(revent->ny2, revent->y2, MWPC2_Y_offset); + } + + if (x1c > -50. && y1c > -50. && x2c > -50. && y2c > -50.) { + xtc = x1c - (x1c - x2c)*((l12+lt)/l12); + ytc = y1c - (y1c - y2c)*((l12+lt)/l12); +// xtc = (l12*x1c + (l12+lt)*(x2c-x1c))/(l12 - (x2c-x1c)*TMath::Tan(TMath::DegToRad()*12.)); +// ytc = (y2c-y1c)*(xtc-x1c)/(x2c-x1c) + y1c; + } + } + + if (x1MultC!=1 || y1MultC!=1 || x2MultC!=1 || y2MultC!=1) continue; + + + + /////////////////////////////////////////////// + //work with Si detectors + /////////////////////////////////////////////// //cout << endl; + //CdI left -veto: + for(Int_t i = 0; i<16; i++) { + if (revent->CsI_L[i]>CsIleftThr) CsI_L_veto = kTRUE; + } + + for (Int_t i = 0; i < 32; i++) { // cout << pSQX_L_EC.Energy(1, i) << endl; SQLXE[i] = 0; @@ -256,19 +534,16 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { SQLXtime[i] = 0; SQLXEtimeFiltered[i] = 0.; - energy = pSQX_L_EC.Energy(revent->SQX_L[i], i); - -// if (i == 5 && i == 6) continue; - - if (energy > kSQLX_energy_thr) { + if (revent->SQX_L[i] > thickLXChThresholds[i]) { + energy = pSQX_L_EC.Energy(revent->SQX_L[i], i); SQLXE[i] = energy; SQLXEsum += SQLXE[i]; SQLXmult++; } - if (i<16 && i!=5 && i!=6 && i!=0 && i!=1 - && SQLXE[i]>kSQLX_energy_thr - && revent->tSQX_L[i]*0.3>320 && revent->tSQX_L[i]*0.3<380) { + if (i<16 && i!=0 && i!=1 + && SQLXE[i]>0 + && revent->tSQX_L[i]*0.3>minTimeLX[i] && revent->tSQX_L[i]*0.3tSQX_L[i]*0.3; SQLXtimeMult++; // SQLXtimeSum = SQLXtimeSum + SQLXtime[i]; @@ -278,8 +553,8 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { // cout << i << endl; } if (i>15 - && SQLXE[i]>kSQLX_energy_thr - && revent->tSQX_L[i]*0.3>410 && revent->tSQX_L[i]*0.3<470) { + && SQLXE[i]>0 + && revent->tSQX_L[i]*0.3>minTimeLX[i] && revent->tSQX_L[i]*0.3tSQX_L[i]*0.3; SQLXtimeMult++; // SQLXtimeSum = SQLXtimeSum + SQLXtime[i]; @@ -295,7 +570,7 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { SQRXEsum += SQRXE[i]; SQRXmult++; } - } + }// for 32 for (Int_t i = 0; i < 16; i++) { // cout << pSQX_L_EC.Energy(1, i) << endl; @@ -305,32 +580,33 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { SQ20time[i] = 0.; SQ20Ecorr[i] = 0.; SQ20EcorrHit[i] = 0.; + SQ20EcorrHitC[i] = 0.; // SQRXE[i] = 0; SQLYtime[i] = 0.; SQLYEtimeFiltered[i] = 0.; - energy = pSQ20_EC.Energy(revent->SQ20[i], i); - -// if (energy > 1.) { + if (revent->SQ20[i] > thinChThresholds[i]) { + energy = pSQ20_EC.Energy(revent->SQ20[i], i); SQ20E[i] = energy; SQ20Esum += SQ20E[i]; // cout << energy << endl; // cout << i << "\t" << energy << "\t" << pSQ20_EC.a[i][0] << "\t" << pSQ20_EC.a[i][1] << "\t" << pSQ20_EC.a[i][2] << endl; -// } + } SQ20time[i] = revent->tSQ20[i]*0.3; - if (SQ20E[i]>kSQL20_energy_thr && SQ20time[i]>400 && SQ20time[i]<500) SQ20timeMult++; + if (SQ20E[i]>0 && SQ20time[i]>minTimeLthin[i] && SQ20time[i]SQY_L[i], i); // cout << energy << "\t" << pSQY_L_EC.a[i][0] << "\t" << pSQY_L_EC.a[i][1] << "\t" << pSQY_L_EC.a[i][2] << endl; // cout << energy << endl; - if (energy > kSQLY_energy_thr) { + if (revent->SQY_L[i] > thickLYChThresholds[i]) { + energy = pSQY_L_EC.Energy(revent->SQY_L[i], i); SQLYE[i] = energy; SQLYEsum += SQLYE[i]; SQLYmult++; @@ -338,23 +614,13 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { Double_t time = revent->tSQY_L[i]*0.3 + timeCorr[i]; - if (SQLYE[i]>kSQLY_energy_thr && time>325 && time<333) { -// cout << "!!!!!!!!!!!!" << endl; + if (SQLYE[i]>0 && time>325 && time<333) { SQLYtime[i] = time; SQLYtimeMult++; SQLYEtimeFiltered[i]=energy; SQLYEtimeFilteredSum = SQLYEtimeFilteredSum + SQLYEtimeFiltered[i]; } - - - // energy = pSQX_R_EC.Energy(revent->SQX_R[i], i); - // - // if (energy > 1.0) { - // SQRXE[i] = energy; - // SQRXEsum += SQRXE[i]; - // SQRXmult++; - // } }//for 16 // cout << SQLXtimeMult << endl; @@ -371,11 +637,11 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { for (Int_t yi = 0; yi < kSQL_Y_strips; yi++) { // cout << "y bin: " << yi+1 << "\t\t"; - if (SQLYE[yi]>kSQLY_energy_thr) { + if (SQLYE[yi]>0) { for (Int_t xi = 0; xi < kSQL_X_strips; xi++) { - if (SQLXE[xi] > kSQLX_energy_thr) { + if (SQLXE[xi] > 0) { currThickness = hThickness->GetBinContent(xi+1, yi+1); // if (currThickness<30.) { // Int_t probableThinStrip = xi-2/2; @@ -383,7 +649,7 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { // } for (Int_t x20 = 0; x20 < kSQL_20_strips; x20++) { - if (SQ20E[x20] > kSQL20_energy_thr/*SQ20timeMult==1*/ && currThickness<30.) { + if (SQ20E[x20] > 0/*SQ20timeMult==1*/ && currThickness<30.) { // SQ20Ecorr[x20] = SQ20E[x20] + 1.; // currThickness = kSQ20_norm_thickness + (currThickness - kSQ20_norm_thickness)*0.5; @@ -409,48 +675,6 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { }//if for correction according energy multiplicity - /////////////////////////////////////////////// - //beam diagnostics - /////////////////////////////////////////////// - - dEbeam = 0.; - TOF = 0.; - if (revent->tF5[0]>0 && revent->tF5[1]>0 && revent->tF5[2]>0 && revent->tF5[3]>0 - && revent->tF3[0]>0 && revent->tF3[1]>0 && revent->tF3[2]>0 && revent->tF3[3]>0 - && revent->F5[0]>0 && revent->F5[1]>0 && revent->F5[2]>0 && revent->F5[3]>0 -// && ( (revent->tF5[0]+revent->tF5[1]+revent->tF5[2]+revent->tF5[3]) - (revent->tF3[0]+revent->tF3[1]+revent->tF3[2]+revent->tF3[3]) )/4.*0.125+89.165<200 -// && ( (tF5[0]+tF5[1]+tF5[2]+tF5[3]) - (tF3[0]+tF3[1]+tF3[2]+tF3[3]) )/4.*0.125+89.165>100 -// && (F5[0]+F5[1]+F5[2]+F5[3])/4. < 2500 - ) { - dEbeam = (revent->F5[0]+revent->F5[1]+revent->F5[2]+revent->F5[3])/4.; - TOF = ( (revent->tF5[0]+revent->tF5[1]+revent->tF5[2]+revent->tF5[3]) - (revent->tF3[0]+revent->tF3[1]+revent->tF3[2]+revent->tF3[3]) )/4.*0.125+89.165; - } - - const Float_t l12 = 546.; - const Float_t lt = 270.; - - x1p = -100.; - y1p = -100.; - x2p = -100.; - y2p = -100.; - xt = -100.; - yt = -100.; - - if (revent->x1[0]<1000 && revent->y1[0]<1000 && revent->nx1==1 && revent->ny1==1) { - x1p = revent->x1[0]*1.25-20.; - y1p = revent->y1[0]*1.25-20.; - } - - if (revent->x2[0]<1000 && revent->y2[0]<1000 && revent->nx2==1 && revent->ny2==1) { - x2p = revent->x2[0]*1.25-20.; - y2p = revent->y2[0]*1.25-20.; - } - - if (x1p > -50. && y1p > -50. && x2p > -50. && y2p > -50.) { - xt = x1p - (x1p - x2p)*(l12/lt); - yt = y1p - (y1p - y2p)*(l12/lt); - } - /////////////////////////////////////////////// @@ -469,56 +693,42 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { if (SQLXtimeMult==1) { for (Int_t i = 0; i<32; i++) { - if (SQLXEtimeFiltered[i]>0) x1mm = -30. + (i+1/2)*60./32.; + if (SQLXEtimeFiltered[i]>0) x1mm = -30. + (i+1./2.)*60./32.; } } // cout << SQLYtimeMult << endl; if (SQLYtimeMult==1) { for (Int_t i = 0; i<16; i++) { - if (SQLYEtimeFiltered[i]>0) y1mm = -30. + (i+1/2)*60./16.; + if (SQLYEtimeFiltered[i]>0) y1mm = -30. + (i+1./2.)*60./16.; } } + //for MWPC wire multiplicity == 1 if (SQLXtimeMult==1 && SQLYtimeMult==1 && xt>-50. && yt>-50.) { vHit1mm.SetXYZ(x1mm-xt, y1mm-yt, z1mm); -// cout << vHit1mm.X() << "\t" << vHit1mm.Y() << endl; -// cout << "\t" << zThin/z1mm*vHit1mm.X() << "\t" << zThin/z1mm*vHit1mm.Y() << endl; vHit1mm *= zThin/z1mm; -// cout << "\t" << vHit1mm.X() << "\t" << vHit1mm.Y() << "\t" << vHit1mm.Z() << endl; -// cout << "\t" << xt << "\t" << yt << "\t" << 0 << endl; -// cout << zThin/z1mm << endl; angleLeft = vHit1mm.Angle(vNorm); -// xThin = vHit1mm.X()+xt; -// yThin = vHit1mm.Y()+yt; - -// cout << angleLeft << endl; } -// cout << angleLeft << endl; - -// cout << x1mm << "\t" << y1mm << endl; - -// cout << vHit1mm.X() << "\t" << vHit1mm.Y() << endl; - - if (SQLXtimeMult==1 && SQLYtimeMult==1 && vHit1mm.X()+xt>-30. && vHit1mm.X()+xt<30. && vHit1mm.Y()+yt>-30. && vHit1mm.Y()+yt<30.) { - mapYbin = (vHit1mm.Y()+yt+30.-0.5*60./16.)*16./60.; - mapXbin = (vHit1mm.X()+xt+30.-0.5*60./32.)*32./60.; + + + mapYbin = (vHit1mm.Y()+yt+thinYoffset+30.-0.5*60./16.)*16./60.; + mapXbin = (vHit1mm.X()+xt+thinXoffset+30.-0.5*60./32.)*32./60.; + //coordinate in the thin detector plane xThin = vHit1mm.X()+xt; yThin = vHit1mm.Y()+yt; -// cout << "\t\t" << xThin << "\t" << yThin << endl; - -// cout << mapYbin << "\t" << "\t" << y1mm << "\t" << xt << endl; //calculation of corrected energy currThickness = hThickness->GetBinContent(mapXbin+1, mapYbin+1); + for (Int_t x20 = 0; x20 < kSQL_20_strips; x20++) { - if (SQ20E[x20] > kSQL20_energy_thr/*SQ20timeMult==1*/ && currThickness<30.) { + if (SQ20E[x20] > 0/*SQ20timeMult==1*/ && currThickness<30.) { SQ20EcorrHit[x20] = (kSQ20_norm_thickness/currThickness)*SQ20E[x20]; SQ20EcorrHitMult++; // cout << mapXbin << "\t" << mapYbin << "\t" << x20 @@ -530,13 +740,58 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { }//if inside the map -// cout << mapYbin << "\t" << "\t" << y1mm << endl; -// cout << mapXbin << "\t" << xThin << "\t" << x1mm << endl; + //end of MWPC wire multiplicity == 1 + //for MWPC cluster multiplicity == 1 -// cout << mapXbin << "\t" << mapYbin << endl; + if (SQLXtimeMult==1 && SQLYtimeMult==1 && xtc>-50. && ytc>-50.) { + vHit1mm.SetXYZ(x1mm-xtc, y1mm-ytc, z1mm); + vHit1mm *= zThin/z1mm; + angleLeft = vHit1mm.Angle(vNorm); + } - tw->Fill(); + + if (SQLXtimeMult==1 && SQLYtimeMult==1 + && vHit1mm.X()+xtc>-30. && vHit1mm.X()+xtc<30. + && vHit1mm.Y()+ytc>-30. && vHit1mm.Y()+ytc<30.) { + + + mapYbin = (vHit1mm.Y()+ytc+thinYoffset+30.-0.5*60./16.)*16./60.; + mapXbin = (vHit1mm.X()+xtc+thinXoffset+30.-0.5*60./32.)*32./60.; + //coordinate in the thin detector plane + xThin = vHit1mm.X()+xtc; + yThin = vHit1mm.Y()+ytc; + + //calculation of corrected energy + currThickness = hThickness->GetBinContent(mapXbin+1, mapYbin+1); + + for (Int_t x20 = 0; x20 < kSQL_20_strips; x20++) { + if (SQ20E[x20] > 0/*SQ20timeMult==1*/ && currThickness<30.) { + SQ20EcorrHitC[x20] = (kSQ20_norm_thickness/currThickness)*SQ20E[x20]; + SQ20EcorrHitCMult++; + // cout << mapXbin << "\t" << mapYbin << "\t" << x20 + // << "\t" << currThickness + // << "\t" << SQ20E[x20] << "\t" << SQ20EcorrHit[x20] << "\t" << endl; + + } + }//for*/ + + if (SQ20EcorrHitCMult==1) { + expectedThinStrip = (xThin + 25. - thinXoffset)*(16./50.)-0.5; + } + + }//if inside the map + + //end of MWPC wire multiplicity == 1 + + if ( + SQLXtimeMult>0 && SQLXEtimeFilteredSum>0 + && trigger==3 + && x1MultC==1 && y1MultC==1 && x2MultC==1 && y2MultC==1 + ) { + + tw->Fill(); + }//if TTree::Fill if (eventNo%100000 == 0 && eventNo !=0) { cout << eventNo << " events processed." << endl; @@ -552,8 +807,15 @@ void fillTree(Int_t nofile = 0, const Int_t noevents = 0) { } -void fillChain(const Int_t from = 0, const Int_t to = 9, const Int_t noevents = 0) { +void fillChain(const TString beam = "he", const Int_t from = 0, const Int_t to = 9, const Int_t noevents = 0) { + + TStopwatch stopwatch; + stopwatch.Start(); + for (Int_t i = from; i <= to; i++) { - fillTree(i, noevents); + fillTree(beam, i, noevents); } + + cout << "Finished in "<< stopwatch.RealTime() << " seconds" << endl; + cout << "Finished in "<< stopwatch.RealTime()/60. << " minutes" << endl; } diff --git a/analysis/rootlogon.C b/analysis/rootlogon.C index 3f20d1793e190a1d1547f9998f28cbb35193f216..3863e8bf86ead29e66101bcc3f810151387be775 100644 --- a/analysis/rootlogon.C +++ b/analysis/rootlogon.C @@ -4,4 +4,6 @@ // gInterpreter->AddIncludePath("/home/vratik/workspace/AculUtils/TELoss"); gSystem->Load("../libUserAnalysis.so"); // gSystem->Load("../../AculUtils/libTELoss.so"); + + gStyle->SetOptStat(1); } diff --git a/analysis/searchOfThresholds.cxx b/analysis/searchOfThresholds.cxx new file mode 100644 index 0000000000000000000000000000000000000000..c2813096fef3a2fc3f320f08b7d3d0b2b375d445 --- /dev/null +++ b/analysis/searchOfThresholds.cxx @@ -0,0 +1,382 @@ +#include "TSystem.h" +#include "TFile.h" +#include "TTree.h" +#include "TChain.h" +#include "TCanvas.h" +#include "TBox.h" +#include "TCut.h" +#include "TCutG.h" +#include "TStopwatch.h" + +using std::cout; +using std::endl; + +void searchOfThresholds() { + + Long64_t drawEntries = 10000000; + +// TFile *fr = new TFile("~/data/exp1804/h5_11_00.root"); + TChain *tr = new TChain("AnalysisxTree"); + tr->Add("~/data/exp1804/h5_14_0?.root"); +// tr->Add("~/data/exp1804/be10_0?_?0.root"); + + TChain *trCal = new TChain("cal"); +// trCal->Add("~/data/exp1804/h5_11_0?_calib.root"); +// trCal->Add("~/data/exp1804/h5_14_0?_calib.root"); +// trCal->Add("~/data/exp1804/be10_0?_?0_calib.root"); + + + + //calibration file for thin detector + TFile *fcal = new TFile("~/data/exp1804/calib/si_20_03.root"); + TTree *tAlphaCal = (TTree*)fcal->Get("AnalysisxTree"); + + //calibration file for 1mm left detector + TFile *fcal1000 = new TFile("~/data/exp1804/calib/si_1000_LR_02.root"); + TTree *tAlphaCal1000 = (TTree*)fcal1000->Get("AnalysisxTree"); + + tr->GetListOfFiles()->Print(); + trCal->GetListOfFiles()->Print(); + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + + TString varName; + TString condition; + + Int_t firstThinStrip = 5; + + + TStopwatch stopwatch; + stopwatch.Start(); + + +//return; + + +///////////////////////////////////////////////////////////////////// +// c0_1: thresholds search for SQ20 in channels: + + + /*TCanvas *c0_1 = new TCanvas("c0_1", "amplitudes in ADC channel unit", 1600, 1200); + c0_1->Divide(4,4); + + firstThinStrip = 0; + + //according to calibration data +// Int_t thinChThresholds[16] = { +// 110, 112, 113, 115, +// 112, 118, 113, 110, +// 104, 115, 120, 113, +// 115, 119, 113, 102 +// }; + + //according to experimental data + Int_t thinChThresholds[16] = { + 113, 121, 125, 125, + 125, 132, 130, 115, + 118, 128, 132, 130, + 130, 130, 130, 112 + }; + + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 10000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 16; i++) { + c0_1->cd(i+1); + gPad->SetLogy(); + + varName.Form("SQ20[%d]", i+firstThinStrip); + condition.Form("SQ20[%d]>0 && SQ20[%d]<1200", i+firstThinStrip, i+firstThinStrip); + + tAlphaCal->SetLineColor(kBlack); + tAlphaCal->Draw(varName, condition, "", drawEntries); + c0_1->Update(); + + condition.Form("SQ20[%d]>%d && SQ20[%d]<1200", i+firstThinStrip, thinChThresholds[i], i+firstThinStrip); + cout << condition << endl; + + tAlphaCal->SetLineColor(kRed); + tAlphaCal->Draw(varName, condition, "same", drawEntries); + c0_1->Update(); + + + }//*/ + +// return; + + +///////////////////////////////////////////////////////////////////// +// c0_1_check: check of thresholds for SQ20 in channels +// on real experimental data: + + + /*TCanvas *c0_1_check = new TCanvas("c0_1_check", "amplitudes in ADC channel unit", 1600, 1200); + c0_1_check->Divide(4,4); + + firstThinStrip = 0; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 16; i++) { + c0_1_check->cd(i+1); + gPad->SetLogy(); + + varName.Form("SQ20[%d]", i+firstThinStrip); + condition.Form("SQ20[%d]>0 && SQ20[%d]<200", i+firstThinStrip, i+firstThinStrip); + + tr->SetLineColor(kBlack); + tr->Draw(varName, condition, "", drawEntries); + c0_1_check->Update(); + + condition.Form("SQ20[%d]>%d && SQ20[%d]<200", i+firstThinStrip, thinChThresholds[i], i+firstThinStrip); + cout << condition << endl; + + tr->SetLineColor(kRed); + tr->Draw(varName, condition, "same", drawEntries); + c0_1_check->Update(); + + + }//*/ + +// return; + + +///////////////////////////////////////////////////////////////////// +// c0_2: thresholds search for SQX_L in channels: + + + /*TCanvas *c0_2 = new TCanvas("c0_2", "amplitudes in ADC channel unit (0-15)", 1600, 1200); + c0_2->Divide(4,4); + TCanvas *c0_2_2 = new TCanvas("c0_2_2", "amplitudes in ADC channel unit (15-32)", 1600, 1200); + c0_2_2->Divide(4,4); + + firstThinStrip = 0; + + Int_t thickChThresholds[32] = { + 60, 45, 45, 45, + 45, 4000, 4000, 45, + 45, 45, 45, 45, + 45, 45, 45, 45, + ///// + 45, 45, 45, 45, + 45, 45, 45, 45, + 45, 45, 45, 45, + 45, 45, 45, 60 + }; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 10000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 32; i++) { + if (i<16) c0_2->cd(i+1); + if (i>15) c0_2_2->cd(i-16+1); + gPad->SetLogy(); + + varName.Form("SQX_L[%d]", i+firstThinStrip); + condition.Form("SQX_L[%d]>200 && SQX_L[%d]<2000", i+firstThinStrip, i+firstThinStrip); + + tAlphaCal1000->SetLineColor(kBlack); + tAlphaCal1000->Draw(varName, condition, "", drawEntries); + c0_2->Update(); + c0_2_2->Update(); + + condition.Form("SQX_L[%d]>%d && SQX_L[%d]<2000", i+firstThinStrip, thickChThresholds[i], i+firstThinStrip); + cout << condition << endl; + + tAlphaCal1000->SetLineColor(kRed); + tAlphaCal1000->Draw(varName, condition, "same", drawEntries); + c0_2->Update(); + c0_2_2->Update(); + + + }//*/ + +// return; + +///////////////////////////////////////////////////////////////////// +// c0_2_check: check of thresholds for SQX_L in channels +// on real experimental data: + + + /*TCanvas *c0_2_check = new TCanvas("c0_2_check", "amplitudes in ADC channel unit (0-15); exp", 1600, 1200); + c0_2_check->Divide(4,4); + TCanvas *c0_2_2_check = new TCanvas("c0_2_2_check", "amplitudes in ADC channel unit (15-32); exp", 1600, 1200); + c0_2_2_check->Divide(4,4); + + firstThinStrip = 0; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 32; i++) { + if (i<16) c0_2_check->cd(i+1); + if (i>15) c0_2_2_check->cd(i-16+1); + gPad->SetLogy(); + + varName.Form("SQX_L[%d]", i+firstThinStrip); + condition.Form("SQX_L[%d]>0 && SQX_L[%d]<400", i+firstThinStrip, i+firstThinStrip); + + tr->SetLineColor(kBlack); + tr->Draw(varName, condition, "", drawEntries); + c0_2_check->Update(); + c0_2_2_check->Update(); + + condition.Form("SQX_L[%d]>%d && SQX_L[%d]<400", i+firstThinStrip, thickChThresholds[i], i+firstThinStrip); + cout << condition << endl; + + tr->SetLineColor(kRed); + tr->Draw(varName, condition, "same", drawEntries); + c0_2_check->Update(); + c0_2_2_check->Update(); + + + }//*/ + +// return; + + +///////////////////////////////////////////////////////////////////// +// c0_3: thresholds search for SQY_L in channels: + + + /*TCanvas *c0_3 = new TCanvas("c0_3", "amplitudes in ADC channel unit (0-15); Y-strips", 1600, 1200); + c0_3->Divide(4,4); + + firstThinStrip = 0; + + //according to calibration data +// Int_t thickYChThresholds[16] = { +// 4000, 130, 80, 110, +// 90, 120, 80, 110, +// 80, 110, 80, 110, +// 80, 110, 90, 115, +// }; + + //according to experimental data + Int_t thickYChThresholds[16] = { + 4000, 130, 82, 110, + 87, 112, 82, 108, + 82, 112, 82, 108, + 80, 110, 90, 105, + }; + + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 10000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 16; i++) { + c0_3->cd(i+1); + + gPad->SetLogy(); + + varName.Form("SQY_L[%d]", i+firstThinStrip); + condition.Form("SQY_L[%d]>0 && SQY_L[%d]<2000", i+firstThinStrip, i+firstThinStrip); + + tAlphaCal1000->SetLineColor(kBlack); + tAlphaCal1000->Draw(varName, condition, "", drawEntries); + c0_3->Update(); + + condition.Form("SQY_L[%d]>%d && SQY_L[%d]<2000", i+firstThinStrip, thickYChThresholds[i], i+firstThinStrip); + cout << condition << endl; + + tAlphaCal1000->SetLineColor(kRed); + tAlphaCal1000->Draw(varName, condition, "same", drawEntries); + c0_3->Update(); + + + }//*/ + +// return; + +///////////////////////////////////////////////////////////////////// +// c0_3_check: check of thresholds for SQY_L in channels +// on real experimental data: + + + /*TCanvas *c0_3_check = new TCanvas("c0_3_check", "amplitudes in ADC channel unit (0-15); Y-strips; exp", 1600, 1200); + c0_3_check->Divide(4,4); + + firstThinStrip = 0; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 16; i++) { + c0_3_check->cd(i+1); + + gPad->SetLogy(); + + varName.Form("SQY_L[%d]", i+firstThinStrip); + condition.Form("SQY_L[%d]>0 && SQY_L[%d]<200", i+firstThinStrip, i+firstThinStrip); + + tr->SetLineColor(kBlack); + tr->Draw(varName, condition, "", drawEntries); + c0_3_check->Update(); + + condition.Form("SQY_L[%d]>%d && SQY_L[%d]<200", i+firstThinStrip, thickYChThresholds[i], i+firstThinStrip); + cout << condition << endl; + + tr->SetLineColor(kRed); + tr->Draw(varName, condition, "same", drawEntries); + c0_3_check->Update(); + + }//*/ + + + +///////////////////////////////////////////////////////////////////// +// c0_3_check: check of thresholds for SQY_L in channels +// on real experimental data: + + + TCanvas *c2_check = new TCanvas("c2_check", "amplitudes in ADC channel unit; CsI_left; exp", 1200, 1100); + c2_check->Divide(4,4); + + firstThinStrip = 0; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 16; i++) { + c2_check->cd(i+1); + +// gPad->SetLogy(); + + varName.Form("CsI_L[%d]", i+firstThinStrip); + condition.Form("CsI_L[%d]>0 && CsI_L[%d]<2000", i+firstThinStrip, i+firstThinStrip); + + tr->SetLineColor(kBlack); + tr->Draw(varName, condition, "", drawEntries); + c2_check->Update(); + + condition.Form("CsI_L[%d]>%d && CsI_L[%d]<2000", i+firstThinStrip, 180, i+firstThinStrip); + cout << condition << endl; + + tr->SetLineColor(kRed); + tr->Draw(varName, condition, "same", drawEntries); + c2_check->Update(); + + }//*/ + + + return; + + +} diff --git a/analysis/searchTimeConditions.cxx b/analysis/searchTimeConditions.cxx new file mode 100644 index 0000000000000000000000000000000000000000..7b888aab151f8810ef4db4c22cfc620db5539b9e --- /dev/null +++ b/analysis/searchTimeConditions.cxx @@ -0,0 +1,337 @@ +#include "TSystem.h" +#include "TFile.h" +#include "TTree.h" +#include "TChain.h" +#include "TCanvas.h" +#include "TBox.h" +#include "TCut.h" +#include "TCutG.h" +#include "TStopwatch.h" + +using std::cout; +using std::endl; + +void searchTimeConditions() { + + Long64_t drawEntries = 10000000; + +// TFile *fr = new TFile("~/data/exp1804/h5_11_00.root"); + TChain *tr = new TChain("AnalysisxTree"); + tr->Add("~/data/exp1804/h5_14_0?.root"); +// tr->Add("~/data/exp1804/be10_0?_?0.root"); + + tr->SetMarkerStyle(20); + tr->SetMarkerSize(0.4); + + + TChain *trCal = new TChain("cal"); +// trCal->Add("~/data/exp1804/h5_11_0?_calib.root"); + trCal->Add("~/data/exp1804/h5_14_0?_calib.root"); +// trCal->Add("~/data/exp1804/be10_0?_?0_calib.root"); + + tr->AddFriend(trCal); + + TFile *fcal = new TFile("~/data/exp1804/calib/si_20_03_calib.root"); + TTree *tAlphaCal = (TTree*)fcal->Get("cal"); + tAlphaCal->SetMarkerColor(kRed); + tAlphaCal->SetMarkerStyle(20); + tAlphaCal->SetMarkerSize(.6); + +// TFile *fc = new TFile("cutsIdentification.root", "READ"); +// TCutG *cUS = (TCutG*)fc->Get("c7UpperShadow"); +// TCutG *cLS = (TCutG*)fc->Get("c7LowerShadow"); +// cLS->SetLineColor(kMagenta); +// TCutG *cMA = (TCutG*)fc->Get("c7MainAlpha"); +// TCutG *cA = (TCutG*)fc->Get("c7All"); +// TCutG *cBL = (TCutG*)fc->Get("cBeamLeft"); +// TCutG *cBR = (TCutG*)fc->Get("cBeamRight"); + + +// const Int_t drawEntries = 6000000; +// const Int_t drawEntries = tr->GetEntries(); + + tr->GetListOfFiles()->Print(); + trCal->GetListOfFiles()->Print(); + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + + TString varName; + TString condition; + +// cout << drawEntries << " entries in input chain will be processed." << endl; + + TString timeSiLeftCondition; + timeSiLeftCondition.Form("(SQLXtime[0]>0"); + for (Int_t j = 1; j<32; j++) { + condition.Form("|| SQLXtime[%d]>0", j); + timeSiLeftCondition.Append(condition); + } + timeSiLeftCondition.Append(")"); + + Int_t firstThinStrip = 5; + + Int_t conTimeMult = 1; + + TStopwatch stopwatch; + stopwatch.Start(); + +///////////////////////////////////////////////////////////////////// +// c1 + + + TCanvas *c1 = new TCanvas("c1", "time in thin detector", 1600, 1200); + c1->Divide(4,4); + + firstThinStrip = 0; + + Int_t minTimeLthin[16] { + 410, 410, 420, 415, + 420, 410, 415, 415, + 420, 420, 420, 420, + 410, 425, 0, 0 + }; + + Int_t maxTimeLthin[16] { + 520, 525, 520, 530, + 520, 530, 525, 525, + 530, 530, 530, 520, + 510, 525, 0, 0 + }; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 10000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 15; i++) { + c1->cd(i+1); + + varName.Form("SQ20E[%d]:tSQ20[%d]*0.3", i+firstThinStrip, i+firstThinStrip); + condition.Form("tSQ20[%d]*0.3>300 && tSQ20[%d]*0.3<700 && SQ20E[%d]>0 && SQ20E[%d]<4", + i+firstThinStrip, i+firstThinStrip, i+firstThinStrip, i+firstThinStrip); + + tr->SetMarkerColor(kBlack); + tr->Draw(varName, condition, "col", drawEntries); + c1->Update(); + + condition.Form("tSQ20[%d]>0 && SQ20E[%d]>0 && SQ20E[%d]<40 && SQ20timeMult==%d", + i+firstThinStrip, i+firstThinStrip, i+firstThinStrip, conTimeMult); + tr->SetMarkerColor(kRed); + tr->Draw(varName, condition, "same", drawEntries); + c1->Update(); + + }//*/ + +// return; + +///////////////////////////////////////////////////////////////////// +// c2 + + TCanvas *c2 = new TCanvas("c2", "time issues in 32 X-strips: 0-15", 1600, 1200); + c2->Divide(4,4); + + firstThinStrip = 0; + + Double_t minTimeLX[32] { + 0, 0, 325, 328, + 330, 0, 0, 330, + 326, 328, 326, 335, + 326, 328, 328, 330, + // + 420, 425, 420, 420, + 420, 420, 420, 420, + 430, 425, 420, 420, + 415, 415, 425, 415 + }; + + Double_t maxTimeLX[32] { + 0, 0, 345, 345, + 350, 0, 0, 350, + 345, 345, 345, 350, + 340, 342, 345, 345, + // + 460, 460, 460, 460, + 460, 460, 460, 460, + 460, 460, 460, 460, + 460, 460, 460, 453 + }; + + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 16; i++) { + c2->cd(i+1); + + varName.Form("SQLXE[%d]:0.3*tSQX_L[%d]", i+firstThinStrip, i+firstThinStrip); + condition.Form("0.3*tSQX_L[%d]>310 && 0.3*tSQX_L[%d]<380 && SQLXE[%d]>0 && trigger==3", + i+firstThinStrip, i+firstThinStrip, i+firstThinStrip); + + tr->SetMarkerColor(kRed); + tr->SetMarkerStyle(20); + tr->SetMarkerSize(.4); + tr->Draw(varName, condition, "col", drawEntries); + + varName.Form("SQLXE[%d]:SQLXtime[%d]", i+firstThinStrip, i+firstThinStrip); + condition.Form("SQLXtime[%d]>0 && SQLXE[%d]>0 && SQLXtimeMult==%d", + i+firstThinStrip, i+firstThinStrip, conTimeMult); + + varName.Form("SQLXEtimeFiltered[%d]:0.3*tSQX_L[%d]", i+firstThinStrip, i+firstThinStrip); + condition.Form("0.3*tSQX_L[%d]>200 && 0.3*tSQX_L[%d]<500 && SQLXEtimeFilteredSum>0 && SQLXtimeMult==%d", + i+firstThinStrip, i+firstThinStrip, conTimeMult); + + cout << condition << endl; + + tr->Draw(varName, condition, "same", drawEntries); + + + c2->Update(); + }//for */ + + + +// return; + + +///////////////////////////////////////////////////////////////////// +// c2_1 + + TCanvas *c2_1 = new TCanvas("c2_1", "time issues in 32 X-strips: 16-31", 1600, 1200); + c2_1->Divide(4,4); + + firstThinStrip = 16; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + for (Int_t i = 0; i < 16; i++) { + c2_1->cd(i+1); + + varName.Form("SQLXE[%d]:0.3*tSQX_L[%d]", i+firstThinStrip, i+firstThinStrip); + condition.Form("0.3*tSQX_L[%d]>390 && 0.3*tSQX_L[%d]<480 && SQLXE[%d]>0", i+firstThinStrip, i+firstThinStrip, i+firstThinStrip); + + tr->Draw(varName, condition, "col", drawEntries); + + varName.Form("SQLXE[%d]:SQLXtime[%d]", i+firstThinStrip, i+firstThinStrip); + condition.Form("SQLXtime[%d]>0 && SQLXE[%d]>0 && SQLXtimeMult==%d", + i+firstThinStrip, i+firstThinStrip, conTimeMult); + tr->Draw(varName, condition, "same", drawEntries); + + c2_1->Update(); + }//for */ + +// return; + +///////////////////////////////////////////////////////////////////// +// c3 + + TCanvas *c3 = new TCanvas("c3", "time issues in 16 Y-strips", 1600, 1200); + c3->Divide(4,4); + + firstThinStrip = 0; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + Float_t timeCorr[16] = {0, -2, 0, 2, + -2.5, 5, -6, 4, + 7, 3, 3, 4, + 6, 4, -4, 13 + }; + + + for (Int_t i = 0; i < 16; i++) { + c3->cd(i+1); + + + varName.Form("SQLYE[%d]:0.3*tSQY_L[%d]+%f", i+firstThinStrip, i+firstThinStrip, timeCorr[i]); + condition.Form("0.3*tSQY_L[%d]>300 && 0.3*tSQY_L[%d]<360 && SQLYE[%d]>0", i+firstThinStrip, i+firstThinStrip, i+firstThinStrip); + + tr->SetMarkerColor(kRed); + tr->SetMarkerStyle(20); + tr->SetMarkerSize(.4); + tr->Draw(varName, condition, "col", drawEntries); +// tr->Draw(varName, condition, "", drawEntries); + + varName.Form("SQLYEtimeFiltered[%d]:SQLYtime[%d]", i+firstThinStrip, i+firstThinStrip); + condition.Form("0.3*tSQY_L[%d]>200 && 0.3*tSQY_L[%d]<500 && SQLYEtimeFilteredSum>0 && SQLYtimeMult==%d", + i+firstThinStrip, i+firstThinStrip, conTimeMult); + + cout << condition << endl; + + tr->Draw(varName, condition, "same", drawEntries); + c3->Update(); + + condition.Form("0.3*tSQY_L[%d]+%f>325 && 0.3*tSQY_L[%d]+%f<333 && SQLYE[%d]>0", + i+firstThinStrip, timeCorr[i], i+firstThinStrip, timeCorr[i], i+firstThinStrip); + + tr->SetMarkerColor(kRed); + tr->Draw(varName, condition, "same", drawEntries); + + // tAlphaCal->Draw(varName, "", "same"); + + c3->Update(); + }//for */ + +// return; + +///////////////////////////////////////////////////////////////////// +// c4 + + TCanvas *c4 = new TCanvas("c4", "Multiplicity", 1800, 900); + c4->Divide(3,2); + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + + c4->cd(1); + varName.Form("SQLXmult"); + condition.Form("trigger==3 && SQLXmult>0 && SQLXmult<8"); + tr->Draw(varName, condition, "", drawEntries); + c4->Update(); + + c4->cd(4); + varName.Form("SQLXtimeMult"); + condition.Form("trigger==3 && SQLXtimeMult>0 && SQLXtimeMult<8"); + tr->Draw(varName, condition, "", drawEntries); + c4->Update(); + + c4->cd(2); + varName.Form("SQLYmult"); + condition.Form("trigger==3 && SQLYmult>0 && SQLYmult<8"); + tr->Draw(varName, condition, "", drawEntries); + c4->Update(); + + c4->cd(5); + varName.Form("SQLYtimeMult"); + condition.Form("trigger==3 && SQLYtimeMult>0 && SQLYtimeMult<18"); + tr->Draw(varName, condition, "", drawEntries); + c4->Update(); + + c4->cd(3); + varName.Form("SQ20EcorrMult"); + condition.Form("trigger==3 && SQ20EcorrMult>0"); + tr->Draw(varName, condition, "", drawEntries); + c4->Update(); + + c4->cd(6); + varName.Form("SQ20timeMult"); + condition.Form("trigger==3 && SQ20timeMult>0"); + tr->Draw(varName, condition, "", drawEntries); + c4->Update(); + + //*/ + + return; + +} diff --git a/analysis/showBananas2.cxx b/analysis/showBananas2.cxx index 4d102bec303d3f933cc8ea292a4d9f82e6a85fe6..f543c5879d440fec6af305583922fa814995fb2f 100644 --- a/analysis/showBananas2.cxx +++ b/analysis/showBananas2.cxx @@ -76,7 +76,7 @@ void showBananas2() { // c1 - /*TCanvas *c1 = new TCanvas("c1", "time in thin detector", 1600, 1200); + TCanvas *c1 = new TCanvas("c1", "time in thin detector", 1600, 1200); c1->Divide(4,4); firstThinStrip = 0; @@ -104,7 +104,7 @@ void showBananas2() { ///////////////////////////////////////////////////////////////////// // c2 - /*TCanvas *c2 = new TCanvas("c2", "time issues in 32 X-strips: 0-15", 1600, 1200); + TCanvas *c2 = new TCanvas("c2", "time issues in 32 X-strips: 0-15", 1600, 1200); c2->Divide(4,4); firstThinStrip = 0; @@ -147,7 +147,7 @@ void showBananas2() { -// return; + return; ///////////////////////////////////////////////////////////////////// diff --git a/analysis/showBananasHe.cxx b/analysis/showBananasHe.cxx new file mode 100644 index 0000000000000000000000000000000000000000..ff2b55291daed8daf5a31f58cbe2e1f5c8f9ac3c --- /dev/null +++ b/analysis/showBananasHe.cxx @@ -0,0 +1,566 @@ +#include "TSystem.h" +#include "TFile.h" +#include "TTree.h" +#include "TChain.h" +#include "TCanvas.h" +#include "TBox.h" +#include "TCut.h" +#include "TCutG.h" +#include "TStopwatch.h" +#include "TStyle.h" + +using std::cout; +using std::endl; + +void showBananasHe(Long64_t drawEntries = 100000, const TString beam = "he", TString opt = "") { + +// Long64_t drawEntries = 10000000; + + TChain *trCal = new TChain("cal"); + TChain *trCalBe = new TChain("cal"); + if (beam.EqualTo("he")) trCal->Add("~/data/exp1804/h5_14_0?_calib.root"); + if (beam.EqualTo("be")) trCal->Add("~/data/exp1804/be10_03_?0_calib.root"); + if (beam.Contains("he") && beam.Contains("be")) { + trCalBe->Add("~/data/exp1804/be10_03_?0_calib.root"); + trCal->Add("~/data/exp1804/h5_14_0?_calib.root"); + } + + trCal->GetListOfFiles()->Print(); + + cout << trCal->GetEntries() << " calibrated events." << endl; + + TString varName; + TString condition; + + Int_t firstThinStrip = 5; + Double_t max1mmDeposit = 31.; +// Double_t max1mmDeposit = 65.; + +// Double_t max20deposit = 5.; + Double_t max20deposit = 4.; + + Long64_t drawnEventsCorr2nd[16]; + + + TStopwatch stopwatch; + stopwatch.Start(); + +//return; + + +///////////////////////////////////////////////////////////////////// +// amplitudes and their cuts for left telescope +// may be found in searchOfThresholds.cxx +// +// time conditions for better choice of events +// may be found in searchTimeConditions.cxx +// +///////////////////////////////////////////////////////////////////// + + +///////////////////////////////////////////////////////////////////// +// c1 + + /*TCanvas *c1 = new TCanvas("c1", "1D energy deposits in left Si detectors: different corrections", 1400, 700); + c1->Divide(3,2); + + + + c1->cd(1); + trCal->Draw("SQLXEsum", "SQLXEsum>0 && trigger==3", "", 1000000); + c1->Update(); + + c1->cd(2); + trCal->Draw("SQLXEtimeFilteredSum", "SQLXEtimeFilteredSum>0 && trigger==3", "", 1000000); + c1->Update(); + + c1->cd(4); + trCal->Draw("SQ20E", "SQ20E>0 && trigger==3", "", 1000000); + c1->Update(); + + c1->cd(5); + trCal->Draw("SQ20Ecorr", "SQ20Ecorr>0 && trigger==3", "", 1000000); + c1->Update(); + + c1->cd(6); + trCal->Draw("SQ20EcorrHit", "SQ20EcorrHit>0 && trigger==3", "", 1000000); + c1->Update(); + + //*/ + + +// return; + +///////////////////////////////////////////////////////////////////// +// c2 + + /*TCanvas *c2 = new TCanvas("c5", "Hits of reaction products in left telescope", 1400, 900); + c2->Divide(3,3); + +// cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; +// drawEntries = 10000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + c2->Update(); + + c2->cd(1); + varName.Form("yt:xt>>hTarget(32, -40., 40., 32, -40., 40.)"); + // condition.Form("yThin>-80. && xThin>-80."); + condition.Form("trigger==3"); + trCal->Draw(varName, condition, "col", drawEntries); + c2->Update(); + + c2->cd(2); + varName.Form("yThin:xThin>>hThin(32, -40., 40., 16, -40., 40.)"); + condition.Form("trigger==3 && yThin>-80. && xThin>-80."); + // condition.Form("yThin>-80. && xThin>-80."); + trCal->Draw(varName, condition, "col", drawEntries); + c2->Update(); + + + c2->cd(3); + varName.Form("y1mm:x1mm>>h(32, -30., 30., 16, -30., 30.)"); + condition.Form("trigger==3 && y1mm>-80. && x1mm>-80."); + trCal->Draw(varName, condition, "col", drawEntries); + + TBox *boxThinDet = new TBox(25., -25., -25., 25.); + boxThinDet->SetLineColor(kRed); + boxThinDet->SetLineWidth(3); + boxThinDet->SetFillStyle(0); + boxThinDet->Draw("same"); + + Double_t thinXoffset = 1.; + Double_t thinYoffset = -1.8; + + TBox *boxThinDetCorrPosition + = new TBox(25.+thinXoffset, -25.+thinYoffset, -25.+thinXoffset, 25.+thinYoffset); + boxThinDetCorrPosition->SetLineColor(kGreen); + boxThinDetCorrPosition->SetLineWidth(3); + boxThinDetCorrPosition->SetFillStyle(0); + boxThinDetCorrPosition->Draw("same"); + c2->Update(); + + c2->cd(4); + varName.Form("xThin"); + condition.Form("trigger==3 && yThin>-80. && xThin>10. && xThin<20."); + trCal->Draw(varName, condition, "col", drawEntries); + + TLine *lx1 = new TLine(16.67, 0., 16.67, 100000.); + lx1->SetLineColor(kRed); + lx1->SetLineWidth(3); + lx1->Draw("same"); + TLine *lx2 = new TLine(18.542, 0., 18.542, 100000.); + lx2->SetLineColor(kRed); + lx2->SetLineWidth(3); + lx2->Draw("same"); + c2->Update(); + + c2->cd(5); + varName.Form("yThin"); + condition.Form("trigger==3 && yThin>-10. && yThin>10. && xThin>-80."); + trCal->Draw(varName, condition, "col", drawEntries); + + TLine *ly1 = new TLine(12.869, 0., 12.869, 15000.); + ly1->SetLineColor(kRed); + ly1->SetLineWidth(3); + ly1->Draw("same"); + TLine *ly2 = new TLine(16.305, 0., 16.305, 15000.); + ly2->SetLineColor(kRed); + ly2->SetLineWidth(3); + ly2->Draw("same"); + c2->Update(); + + c2->cd(6); + varName.Form("y1mm:x1mm"); + condition.Form("trigger==3 && y1mm>-80. && x1mm>-80."); + condition.Form("trigger==3 && y1mm>-80. && x1mm>-10. && x1mm<10."); + trCal->Draw(varName, condition, "", drawEntries); + + + + c2->cd(7); + varName.Form("mapXbin"); + condition.Form("trigger==3 && mapXbin>-1"); + trCal->Draw(varName, condition, "col", drawEntries); + + TLine *lxMap = new TLine(25.5, 0., 25.5, 745252.); + lxMap->SetLineColor(kRed); + lxMap->SetLineWidth(3); + lxMap->Draw("same"); + c2->Update(); + + c2->cd(8); + varName.Form("mapYbin"); + condition.Form("trigger==3 && mapYbin>-1"); + trCal->Draw(varName, condition, "col", drawEntries); + TLine *lyMap = new TLine(7.5, 0., 7.5, 268000.); + lyMap->SetLineColor(kRed); + lyMap->SetLineWidth(3); + lyMap->Draw("same"); + c2->Update(); + + c2->cd(9); + // varName.Form("SQ20timeMult"); + // condition.Form("trigger==3 && SQ20timeMult>0"); + trCal->Draw("angleLeft", "angleLeft<3", "", drawEntries); + c2->Update(); + + //*/ + +// return; + +///////////////////////////////////////////////////////////////////// +// c3 + + + /*TCanvas *c3 = new TCanvas("c3", "Multiplicity - clusters", 1400, 900); + c3->Divide(3,3); + +// cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + + c3->cd(1); + trCal->Draw("x1MultC", "x1MultC<4", "", drawEntries); + c3->Update(); + + c3->cd(2); + trCal->Draw("y1MultC", "y1MultC<4", "", drawEntries); + c3->Update(); + + c3->cd(4); + trCal->Draw("x2MultC", "x2MultC<4", "", drawEntries); + c3->Update(); + + c3->cd(5); + trCal->Draw("y2MultC", "y2MultC<4", "", drawEntries); + c3->Update(); + //*/ + +// return; + + +///////////////////////////////////////////////////////////////////// +// c11 + + /*TCanvas *c11 = new TCanvas("c11", "dE-E uncorrected", 1600, 800); + c11->Divide(4,2); + + firstThinStrip = 5; + +// Int_t drawnEvents1[8]; + +// cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; +// drawEntries = 1000000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + c11->cd(1); + trCal->Draw("SQ20E[0]:SQ20E[0]+SQLXEsum>>h11(150, 0, 35, 150, 0, 5)", "", "col", 1000); + + for (Int_t i = 0; i < 8; i++) { + c11->cd(i+1); + + if (opt.Contains("logz")) gPad->SetLogz(); + + varName.Form("SQ20E[%d]:SQ20E[%d]+SQLXEsum", i+firstThinStrip, i+firstThinStrip); + condition.Form("SQ20E[%d]>0 && SQ20E[%d]<%f " +// "&& SQLXEsum>1.1 " + "&& SQLXEsum<%f " + // "&& SQLXtimeMult==1" + "&& trigger==3", + i+firstThinStrip, i+firstThinStrip, max20deposit, max1mmDeposit); + + trCal->Draw(varName, condition, "col", drawEntries); + + c11->Update(); + + trCalBe->SetMarkerStyle(20); + trCalBe->SetMarkerColor(kRed); + trCalBe->SetMarkerSize(.4); + trCalBe->Draw(varName, condition, "same", drawEntries); + + TH2F *h = (TH2F*)gPad->GetPrimitive("htemp"); + h->SetStats(); + + c11->Update(); + }//*/ + +// return; + + + + +///////////////////////////////////////////////////////////////////// +// c12 + + /*TCanvas *c12 = new TCanvas("c12", "dE-E corrected (first approximation)", 1600, 800); + c12->Divide(4,2); + +// Int_t drawnEvents2[8]; + + firstThinStrip = 5; + + cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; + drawEntries = 1000000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + c12->cd(1); + trCal->Draw("SQ20Ecorr[0]:SQ20Ecorr[0]+SQLXEsum>>h12(100, 0, 35, 100, 0, 5)", "", "col", 1000); + + for (Int_t i = 0; i < 8; i++) { + c12->cd(i+1); + + if (opt.Contains("logz")) gPad->SetLogz(); + + varName.Form("SQ20Ecorr[%d]:SQ20Ecorr[%d]+SQLXEsum", + i+firstThinStrip, i+firstThinStrip); + condition.Form("" + "SQ20Ecorr[%d]>0 " + "&& SQ20Ecorr[%d]<5 " +// "&& SQLXEsum>1.1 " + "&& SQLXEsum<%f " + "&& SQLXtimeMult==1" + "&& trigger==3", + i+firstThinStrip, + i+firstThinStrip, + max1mmDeposit); + + trCal->Draw(varName, condition, "col", drawEntries); + c12->Update(); + }//*/ + +// return; + +///////////////////////////////////////////////////////////////////// +// c13 + + /*TCanvas *c13 = new TCanvas("c13", "dE-E corrected (second approximation)", 1600, 800); + c13->Divide(4,2); + + firstThinStrip = 5; + +// cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; +// drawEntries = 1000000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + c13->cd(1); +// trCal->Draw("SQ20Ecorr[0]:SQ20Ecorr[0]+SQLXEsum>>h13(150, 0, 35, 150, 0, 5)", "", "col", 1000); + trCal->Draw("SQ20Ecorr[0]:SQ20Ecorr[0]+SQLXEsum>>h13(150, 0, 45, 150, 0, 5)", "", "col", 1000); + + for (Int_t i = 0; i < 8; i++) { + c13->cd(i+1); + + if (opt.Contains("logz")) gPad->SetLogz(); + + drawnEventsCorr2nd[i] = 0; + + varName.Form("SQ20EcorrHit[%d]:SQ20EcorrHit[%d]+SQLXEtimeFilteredSum", i+firstThinStrip, i+firstThinStrip); + condition.Form("" + "SQ20EcorrHit[%d]>0" + " && SQ20EcorrHit[%d]<%f " +// "&& SQLXEsum>1.1" + " && SQLXEsum<%f " + "&& SQLXtimeMult==1" + "&& SQLXEtimeFilteredSum>0.5" + "&& trigger==3", + i+firstThinStrip, + i+firstThinStrip, + max20deposit, + max1mmDeposit); + + + drawnEventsCorr2nd[i] = trCal->Draw(varName, condition, "col", drawEntries); + + TH2F *h = (TH2F*)gPad->GetPrimitive("htemp"); + h->SetStats(); + + c13->Update(); + + }//*/ + +// return; + + +///////////////////////////////////////////////////////////////////// +// c14 + + /*TCanvas *c14 = new TCanvas("c14", "dE-E corrected (second approximation, clusters)", 1600, 800); + c14->Divide(4,2); + + firstThinStrip = 5; + +// cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; +// drawEntries = 1000000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + c14->cd(1); + trCal->Draw("SQ20Ecorr[0]:SQ20Ecorr[0]+SQLXEsum>>h13(150, 0, 35, 150, 0, 5)", "", "col", 1000); + + for (Int_t i = 0; i < 8; i++) { + c14->cd(i+1); + if (opt.Contains("logz")) gPad->SetLogz(); + drawnEventsCorr2nd[i] = 0; + + varName.Form("SQ20EcorrHitC[%d]:SQ20EcorrHitC[%d]+SQLXEtimeFilteredSum", i+firstThinStrip, i+firstThinStrip); + condition.Form("" + "SQ20EcorrHitC[%d]>0" + " && SQ20EcorrHitC[%d]<%f " +// "&& SQLXEsum>1.1" + " && SQLXEsum<%f " + "&& SQLXtimeMult==1" + "&& SQLXEtimeFilteredSum>0.5" + "&& trigger==3", + i+firstThinStrip, + i+firstThinStrip, + max20deposit, + max1mmDeposit); + + + drawnEventsCorr2nd[i] = trCal->Draw(varName, condition, "col", drawEntries); + trCalBe->SetMarkerStyle(20); + trCalBe->SetMarkerColor(kRed); + trCalBe->SetMarkerSize(.4); + trCalBe->Draw(varName, condition, "same", drawEntries); + + TH2F *h = (TH2F*)gPad->GetPrimitive("htemp"); + h->SetStats(); + + c14->Update(); + + }//*/ + +// return; + +///////////////////////////////////////////////////////////////////// +// c15 + + /*TCanvas *c15 = new TCanvas("c15", "dE-E corrected (second approximation, clusters, cleaned)", 1600, 800); + c15->Divide(4,2); + + firstThinStrip = 5; + +// cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; +// drawEntries = 1000000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + c15->cd(1); + trCal->Draw("SQ20Ecorr[0]:SQ20Ecorr[0]+SQLXEsum>>h13(150, 0, 35, 150, 0, 5)", "", "col", 1000); + + for (Int_t i = 0; i < 8; i++) { + c15->cd(i+1); + if (opt.Contains("logz")) gPad->SetLogz(); + drawnEventsCorr2nd[i] = 0; + + varName.Form("SQ20EcorrHitC[%d]:SQ20EcorrHitC[%d]+SQLXEtimeFilteredSum", i+firstThinStrip, i+firstThinStrip); + condition.Form("" + "SQ20EcorrHitC[%d]>0" + " && SQ20EcorrHitC[%d]<%f " +// "&& SQLXEsum>1.1" + " && SQLXEsum<%f " + "&& SQLXtimeMult==1" + "&& SQLXEtimeFilteredSum>0.5" + "&& SQ20EcorrHitCMult" + "&& CsI_L_veto==0" + "&& trigger==3", + i+firstThinStrip, + i+firstThinStrip, + max20deposit, + max1mmDeposit); + + + drawnEventsCorr2nd[i] = trCal->Draw(varName, condition, "col", drawEntries); + trCalBe->SetMarkerStyle(20); + trCalBe->SetMarkerColor(kRed); + trCalBe->SetMarkerSize(.4); + trCalBe->Draw(varName, condition, "same", drawEntries); + + TH2F *h = (TH2F*)gPad->GetPrimitive("htemp"); + h->SetStats(); + + c15->Update(); + + }//*/ + +// return; + + +///////////////////////////////////////////////////////////////////// +// c16 + + TCanvas *c16 = new TCanvas("c16", "dE-E corrected (second approximation, clusters, cleaned, check for thin strip)", 1600, 800); + c16->Divide(4,2); + + firstThinStrip = 5; + +// cout << tr->GetEntries() << " events." << endl; + cout << trCal->GetEntries() << " calibrated events." << endl; +// drawEntries = 1000000000; + cout << drawEntries << " entries in input chain will be processed." << endl; + + gStyle->SetOptStat(); + + c16->cd(1); + trCal->Draw("SQ20Ecorr[0]:SQ20Ecorr[0]+SQLXEsum>>h13(150, 0, 35, 150, 0, 5)", "", "col", 1000); + + for (Int_t i = 0; i < 8; i++) { + c16->cd(i+1); + if (opt.Contains("logz")) gPad->SetLogz(); + drawnEventsCorr2nd[i] = 0; + + varName.Form("SQ20EcorrHitC[%d]:SQ20EcorrHitC[%d]+SQLXEtimeFilteredSum", i+firstThinStrip, i+firstThinStrip); + condition.Form("" + "SQ20EcorrHitC[%d]>0" + " && SQ20EcorrHitC[%d]<%f " +// "&& SQLXEsum>1.1" + " && SQLXEsum<%f " + "&& SQLXtimeMult==1" + "&& SQLXEtimeFilteredSum>0.5" + "&& SQ20EcorrHitCMult" + "&& CsI_L_veto==0" + "&& abs(expectedThinStrip-%d)<=1" + "&& trigger==3", + i+firstThinStrip, + i+firstThinStrip, + max20deposit, + max1mmDeposit, + i+firstThinStrip); + + + drawnEventsCorr2nd[i] = trCal->Draw(varName, condition, "col", drawEntries); + trCalBe->SetMarkerStyle(20); + trCalBe->SetMarkerColor(kRed); + trCalBe->SetMarkerSize(.4); + trCalBe->Draw(varName, condition, "same", drawEntries); + + TH2F *h = (TH2F*)gPad->GetPrimitive("htemp"); + h->SetStats(); + + c16->Update(); + + }//*/ + +// return; + + + + + +// TCanvas *c = new TCanvas(); +// trCal->Draw("expectedThinStrip", "", "", drawEntries); + + + + for (Int_t i = 0; i<8; i++) { + cout << "Events after 2nd correction in strip No." << i << ": " << drawnEventsCorr2nd[i] << endl; + } + + cout << "Finished in "<< stopwatch.RealTime() << " seconds" << endl; + cout << "Finished in "<< stopwatch.RealTime()/60. << " minutes" << endl; + + return; + +} diff --git a/analysis/showBeam.cxx b/analysis/showBeam.cxx index a4d9d406bb108ba85dd430803a3d3d3119639f73..bca881213cc006620e36dacdb8dbfa89f2c11057 100644 --- a/analysis/showBeam.cxx +++ b/analysis/showBeam.cxx @@ -7,21 +7,17 @@ using std::cout; using std::endl; -void showBeam(const Long64_t drawEntries = 100000) +void showBeam(const Long64_t drawEntries = 100000, const TString beam = "he") { TChain *tr = new TChain("AnalysisxTree"); -// tr->Add("~/data/exp1804/h5_14_0?.root"); -// tr->Add("~/data/exp1804/output.root"); - tr->Add("~/data/exp1804/h5_12.root"); + if (beam.Contains("he")) tr->Add("~/data/exp1804/h5_14_0?.root"); + if (beam.Contains("be")) tr->Add("~/data/exp1804/be10_03_?0.root"); TChain *trCal = new TChain("cal"); -// trCal->Add("~/data/exp1804/h5_11_0?_calib.root"); -// trCal->Add("~/data/exp1804/h5_14_0?_calib.root"); -// trCal->Add("~/data/exp1804/h5_14_0?_calib.root"); - -// tr->AddFriend(trCal); + if (beam.Contains("he")) trCal->Add("~/data/exp1804/h5_14_0?_calib.root"); + if (beam.Contains("be")) trCal->Add("~/data/exp1804/be10_03_?0_calib.root"); // const Int_t drawEntries = 6000000; // const Int_t drawEntries = tr->GetEntries(); @@ -37,27 +33,17 @@ void showBeam(const Long64_t drawEntries = 100000) // cout << drawEntries << " entries in input chain will be processed." << endl; +///////////////////////////////////////////////////////////////////// +// c1 + -// TCanvas *c1 = new TCanvas("c1", "R: 0-7", 1600, 800); -// c1->Divide(4,2); TCanvas *c1 = new TCanvas("beam", "beam info", 1200, 1200); c1->Divide(2,2); c1->cd(1); -// tr->Draw("trigger", "", "", drawEntries); c1->cd(3); -// af5 = (NeEvent->F5[0]+NeEvent->F5[1]+NeEvent->F5[2]+NeEvent->F5[3]+4.*gRandom->Uniform())/4. -// tf5 = (NeEvent->tF5[0]+NeEvent->tF5[1]+NeEvent->tF5[2]+NeEvent->tF5[3]+4.*gRandom->Uniform())/4.; -// tf3 = (NeEvent->tF3[0]+NeEvent->tF3[1]+NeEvent->tF3[2]+NeEvent->tF3[3]+4.*gRandom->Uniform())/4.; -// ToF = (tf5 - tf3)*0.125+89.165; -// af5:TOF -// tr->Draw("(NeEvent->F5[0]+NeEvent->F5[1]+NeEvent->F5[2]+NeEvent->F5[3])/4.", "", "", drawEntries); -// tr->Draw("(F5[0]+F5[1]+F5[2]+F5[3])/4.", "", "", drawEntries); -// tr->Draw("(tF5[0]+tF5[1]+tF5[2]+tF5[3])/4.;", "", "", drawEntries); -// tr->Draw("( (tF5[0]+tF5[1]+tF5[2]+tF5[3]) - (tF3[0]+tF3[1]+tF3[2]+tF3[3]) )/4.*0.125+89.165;", -// "tF5[0]>0 && tF5[1]>0 && tF5[2]>0 && tF5[3]>0 && tF3[0]>0 && tF3[1]>0 && tF3[2]>0 && tF3[3]>0", "", drawEntries); tr->Draw("(F5[0]+F5[1]+F5[2]+F5[3])/4.:( (tF5[0]+tF5[1]+tF5[2]+tF5[3]) - (tF3[0]+tF3[1]+tF3[2]+tF3[3]) )/4.*0.125+89.165;", "tF5[0]>0 && tF5[1]>0 && tF5[2]>0 && tF5[3]>0 " "&& tF3[0]>0 && tF3[1]>0 && tF3[2]>0 && tF3[3]>0 " @@ -108,6 +94,9 @@ void showBeam(const Long64_t drawEntries = 100000) h2->GetYaxis()->SetTitleSize(0.06); h2->GetYaxis()->SetTitleOffset(0.75); +///////////////////////////////////////////////////////////////////// +// c2 + TCanvas *c2 = new TCanvas("c2", "Bugs in MWPC", 1200, 1200); c2->Divide(2,2); @@ -126,6 +115,128 @@ void showBeam(const Long64_t drawEntries = 100000) c2->cd(4); tr->Draw("y2[0]", "ny2==0", "", drawEntries); + c2->Update(); + +///////////////////////////////////////////////////////////////////// +// c3 + + + TCanvas *c3 = new TCanvas("c3", "Wire multiplicity", 1200, 1200); + c3->Divide(2,2); + + c3->cd(1); +// tr->Draw("x1[0]", "x2[0]<1000 && nx2==1", "", drawEntries); + tr->Draw("nx1", "nx1<30000", "", drawEntries); + + c3->cd(2); + tr->Draw("ny1", "ny1<30000", "", drawEntries); + + c3->cd(3); + tr->Draw("nx2", "nx2<30000", "", drawEntries); + + c3->cd(4); + tr->Draw("ny2", "ny2<30000", "", drawEntries); + + + +///////////////////////////////////////////////////////////////////// +// c4 + + + TCanvas *c4 = new TCanvas("c4", "MWPC time", 1200, 1200); + c4->Divide(2,2); + + c4->cd(1); +// tr->Draw("x1[0]", "x2[0]<1000 && nx2==1", "", drawEntries); + tr->Draw("tMWPC[0]*0.125", "nx1<30000", "", drawEntries); + tr->SetLineColor(kRed); + tr->Draw("tMWPC[0]*0.125", "nx1<30000 && trigger==1", "same", drawEntries); + tr->SetLineColor(kBlue); + tr->Draw("tMWPC[0]*0.125", "nx1<30000 && trigger==3", "same", drawEntries); + + + tr->SetLineColor(kBlack); + c4->cd(2); + tr->Draw("tMWPC[1]*0.125", "ny1<30000", "", drawEntries); + tr->Draw("tMWPC[1]*0.125-tF5[0]*0.125", "ny1<30000", "", drawEntries); + + c4->cd(3); + tr->Draw("tMWPC[2]*0.125", "nx2<30000", "", drawEntries); + + c4->cd(4); + tr->Draw("tMWPC[3]*0.125", "ny2<30000", "", drawEntries); + + c4->Update(); + +///////////////////////////////////////////////////////////////////// +// c5 + + + TCanvas *c5 = new TCanvas("c5", "Multiplicity filtered for time", 1200, 1200); + c5->Divide(2,2); + + c5->cd(1); + tr->Draw("nx1", "nx1<30000 && (tMWPC[0]*0.125-tF5[0]*0.125)>60 && (tMWPC[0]*0.125-tF5[0]*0.125)<77", "", drawEntries); + + c5->cd(2); + tr->Draw("ny1", "ny1<30000 && (tMWPC[1]*0.125-tF5[0]*0.125)>60 && (tMWPC[1]*0.125-tF5[0]*0.125)<80", "", drawEntries); + + c5->cd(3); + tr->Draw("nx2", "nx2<30000 && (tMWPC[2]*0.125-tF5[0]*0.125)>70 && (tMWPC[2]*0.125-tF5[0]*0.125)<90", "", drawEntries); + + c5->cd(4); + tr->Draw("ny2", "ny2<30000 && (tMWPC[3]*0.125-tF5[0]*0.125)>60 && (tMWPC[3]*0.125-tF5[0]*0.125)<80", "", drawEntries); + + c5->Update(); + +///////////////////////////////////////////////////////////////////// +// c6 + + + TCanvas *c6 = new TCanvas("c6", "Cluster multiplicity", 1500, 1000); + c6->Divide(3,2); + + c6->cd(1); + trCal->Draw("x1MultC", "x1MultC>0", "", drawEntries); + + c6->cd(2); + trCal->Draw("y1MultC", "y1MultC>0", "", drawEntries); + + c6->cd(4); + trCal->Draw("x2MultC", "x2MultC>0", "", drawEntries); + + c6->cd(5); + trCal->Draw("y2MultC", "y2MultC>0", "", drawEntries); + + c6->cd(3); + tr->Draw("nx1", "ny2==1 && nx2==1 && ny1==1 && nx1==1", "", drawEntries); + + c6->cd(6); + trCal->Draw("x1MultC", "x1MultC==1 && y1MultC==1 && x2MultC==1 && y2MultC==1", "", drawEntries); +// tr->Draw("x1", "nx1<1", "", drawEntries); + + c6->Update(); + + +///////////////////////////////////////////////////////////////////// +// c1 + + + TCanvas *c7 = new TCanvas("c7", "beam info - clusters", 1200, 1200); + c7->Divide(2,2); + + c7->cd(1); + trCal->Draw("y1p:x1p", "y1p>-50. && x1p>-50.", "col", drawEntries); + + c7->cd(3); + trCal->Draw("y2p:x2p", "y2p>-50. && x2p>-50.", "col", drawEntries); + + c7->cd(2); + trCal->Draw("y1c:x1c", "y1c>-50. && x1c>-50.", "col", drawEntries); + + c7->cd(4); + trCal->Draw("y2c:x2c", "y2c>-50. && x2c>-50.", "col", drawEntries); + return;