diff --git a/Analyzer.C b/Analyzer.C index 6e33d76..60df793 100644 --- a/Analyzer.C +++ b/Analyzer.C @@ -29,7 +29,7 @@ vector> eCorr; //############################################ histogram declaration -TH2F * heVID; +TH2F * hevID; TH1F * he[NCRYSTAL]; TH2F * hgg[NCRYSTAL][NCRYSTAL]; @@ -37,7 +37,7 @@ TH2F * hgg[NCRYSTAL][NCRYSTAL]; TH2F * hcoin; ///----- after calibration and BGO veto -TH2F * heCalVID; +TH2F * heCalvID; TH1F * heCal[NCRYSTAL]; TH2F * hcoinBGO; @@ -50,8 +50,12 @@ void Analyzer::Begin(TTree * tree){ printf("======================== Histograms declaration\n"); - heVID = new TH2F("heVID", "e vs ID; det ID; e [ch]", NCRYSTAL, 0, NCRYSTAL, rawEnergyRange[1] - rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); - heCalVID = new TH2F("heCalVID", Form("eCal vs ID (BGO veto > %.1f); det ID; Energy [keV]", BGO_threshold), NCRYSTAL, 0, NCRYSTAL, (energyRange[2] - energyRange[1])/energyRange[0], energyRange[1], energyRange[2]); + hevID = new TH2F("hevID", "e vs ID; det ID; e [ch]", NCRYSTAL, 0, NCRYSTAL, rawEnergyRange[1] - rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); + heCalvID = new TH2F("heCalvID", Form("eCal vs ID (BGO veto > %.1f); det ID; Energy [keV]", BGO_threshold), NCRYSTAL, 0, NCRYSTAL, (energyRange[2] - energyRange[1])/energyRange[0], energyRange[1], energyRange[2]); + + hevID->SetNdivisions(-409, "X"); + heCalvID->SetNdivisions(-409, "X"); + for( int i = 0; i < NCRYSTAL; i ++){ he[i] = new TH1F( Form("he%02d", i), Form("e -%02d", i), rawEnergyRange[1] - rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); heCal[i] = new TH1F(Form("heCal%02d", i), Form("eCal -%02d (BGO veto > %.1f); Energy [keV]; count / %d keV", i, BGO_threshold, energyRange[0]), (energyRange[2] - energyRange[1])/energyRange[0], energyRange[1], energyRange[2]); @@ -110,7 +114,7 @@ Bool_t Analyzer::Process(Long64_t entry){ //if( pileup[detID] == 1 ) continue; //======== Fill raw data - heVID->Fill( detID, e[detID]); + hevID->Fill( detID, e[detID]); he[detID]->Fill(e[detID]); @@ -137,7 +141,7 @@ Bool_t Analyzer::Process(Long64_t entry){ } - heCalVID->Fill( detID, eCal); + heCalvID->Fill( detID, eCal); heCal[detID]->Fill(eCal); for( int detJ = detID +1; detJ < NCRYSTAL; detJ++) { @@ -165,11 +169,11 @@ void Analyzer::Terminate(){ cCanvas->cd(1); cCanvas->cd(1)->SetLogz(1); - heVID->Draw("colz"); + hevID->Draw("colz"); cCanvas->cd(2); cCanvas->cd(2)->SetLogz(1); - heCalVID->Draw("colz"); + heCalvID->Draw("colz"); cCanvas->cd(3); cCanvas->cd(3)->SetLogz(1); diff --git a/Monitor_raw.C b/Monitor_raw.C index 0a7437d..cede9ba 100644 --- a/Monitor_raw.C +++ b/Monitor_raw.C @@ -11,7 +11,8 @@ //############################################ User setting -int rawEnergyRange[2] = {500, 6000}; // in ch +int rawEnergyRange[2] = {500, 6000}; // in ch, {min, max} +int energyRange[3] = {2, 40, 2000}; // in keV, {resol, min, max} TString e_corr = "correction_e.dat"; @@ -28,8 +29,10 @@ vector> eCorr; //############################################ histogram declaration -TH1F * hTDiff; -TH2F * hTDiffvEventID; +TH1F * he[NCRYSTAL]; +TH1F * heCal[NCRYSTAL]; + +TH2F * heCalvID; //############################################ BEGIN void Monitor_raw::Begin(TTree * tree){ @@ -38,9 +41,27 @@ void Monitor_raw::Begin(TTree * tree){ NumEntries = tree->GetEntries(); printf("======================== Creating histograms\n"); + + for( int i = 0 ; i < NCRYSTAL ; i++){ + he[i] = new TH1F (Form("he%02d", i), Form("e%02d", i), rawEnergyRange[1]-rawEnergyRange[0], rawEnergyRange[0], rawEnergyRange[1]); + heCal[i] = new TH1F (Form("heCal%02d", i), Form("e%02d (Cali.)", i), (energyRange[2]-energyRange[1])/energyRange[0], energyRange[1], energyRange[2]); + + switch (i%4){ + case 0: he[i]->SetLineColor(2);break; + case 1: he[i]->SetLineColor(kYellow+3);break; + case 2: he[i]->SetLineColor(kGreen+2);break; + case 3: he[i]->SetLineColor(4);break; + } + switch (i%4){ + case 0: heCal[i]->SetLineColor(2);break; + case 1: heCal[i]->SetLineColor(kYellow+3);break; + case 2: heCal[i]->SetLineColor(kGreen+2);break; + case 3: heCal[i]->SetLineColor(4);break; + } + } - hTDiff = new TH1F("hTDiff", "time different between this and next event", 2000, -1000, 1000); - hTDiffvEventID = new TH2F("hTDiffvEventID", "time different between this and next event; eventID; TDiff", 50, 0, maxEvent, 2000, -1000, 1000); + heCalvID = new TH2F("heCalvID", "ID vs Energy (Cali.); ID; Energy", NCRYSTAL, 0, NCRYSTAL, (energyRange[2]-energyRange[1])/energyRange[0], energyRange[1], energyRange[2]); + printf("======================== end of histograms creation.\n"); @@ -61,7 +82,6 @@ Bool_t Monitor_raw::Process(Long64_t entry){ if (ProcessedEntries>NumEntries*Frac-1) { TString msg; msg.Form("%llu", NumEntries/1000); int len = msg.Sizeof(); - //printf(" %3.0f%% (%*llu/%llu k) processed in %6.1f sec | expect %6.1f sec\033[A\n", printf(" %3.0f%% (%*llu/%llu k) processed in %6.1f sec | expect %6.1f sec\n", Frac*100, len, ProcessedEntries/1000,NumEntries/1000,StpWatch.RealTime(), StpWatch.RealTime()/Frac); StpWatch.Start(kFALSE); @@ -73,19 +93,22 @@ Bool_t Monitor_raw::Process(Long64_t entry){ Terminate(); } + b_ID->GetEntry(entry); b_energy->GetEntry(entry); - b_time_stamp->GetEntry(entry); + b_timestamp->GetEntry(entry); - ULong64_t t0 = t; + int detID = mapping[ID]; - if( entry < (Long64_t) NumEntries ) b_time_stamp->GetEntry(entry+1);; + if( 0 <= detID && detID < 100 ){ + he[detID]->Fill(e); + + //double eCal = ApplyCorrection(eCorr, detID, e); + double eCal = eCorr[detID][0]+eCorr[detID][1]*e; + heCal[detID]->Fill(eCal); + + heCalvID->Fill(detID, eCal); + } - ULong64_t t1 = t; - - int tDiff = (int) t1 - t0; - - hTDiff->Fill( tDiff ); - hTDiffvEventID->Fill( entry, tDiff); return kTRUE; } @@ -96,18 +119,35 @@ void Monitor_raw::Terminate(){ printf("============================== finishing.\n"); gROOT->cd(); - TCanvas * cc = new TCanvas("cc", "cc", 2000, 1000); + int nCrystalPerClover = 4; + int nClover = NCRYSTAL / nCrystalPerClover; + TCanvas * cc = new TCanvas("cc", "cc", 2000, 2000); if( cc->GetShowEventStatus() == 0 ) cc->ToggleEventStatus(); - - cc->Divide(2,1); + cc->Divide(1, 9, 0); - cc->cd(1); - hTDiff->Draw(); + for (Int_t i = 0; i < nClover; i++) { + int canvasID = i + 1; + cc->cd(canvasID); + cc->cd(canvasID)->SetGrid(); + cc->cd(canvasID)->SetTickx(2); + cc->cd(canvasID)->SetTicky(2); + cc->cd(canvasID)->SetBottomMargin(0.06); + cc->cd(canvasID)->SetLogy(); + + for( Int_t j = 0; j < nCrystalPerClover; j++){ + int hID = nCrystalPerClover*i+ j ; + heCal[hID]->Draw("same"); + //he[hID]->Draw("same"); + } + } + cc->SetCrosshair(1); - cc->cd(2); - cc->cd(2)->SetGrid(); - hTDiffvEventID->SetMarkerStyle(3); - hTDiffvEventID->Draw(""); + TCanvas * c1 = new TCanvas("c1", "c1", 1000, 1000); + if( c1->GetShowEventStatus() == 0 ) c1->ToggleEventStatus(); + c1->SetLogz(); + c1->SetGridx(); + heCalvID->SetNdivisions(-409, "X"); + heCalvID->Draw("colz"); } diff --git a/Monitor_raw.h b/Monitor_raw.h index 734f612..9bac8a4 100644 --- a/Monitor_raw.h +++ b/Monitor_raw.h @@ -19,15 +19,15 @@ public : // Declaration of leaf types Long64_t evID; - UShort_t detID; + UShort_t ID; UShort_t e; ULong64_t t; // List of branches TBranch *b_data_ID; //! - TBranch *b_det_ID; //! + TBranch *b_ID; //! TBranch *b_energy; //! - TBranch *b_time_stamp; //! + TBranch *b_timestamp; //! Monitor_raw(TTree * /*tree*/ =0) : fChain(0) { } virtual ~Monitor_raw() { } @@ -59,9 +59,9 @@ void Monitor_raw::Init(TTree *tree) fChain->SetMakeClass(1); fChain->SetBranchAddress("evID", &evID, &b_data_ID); - fChain->SetBranchAddress("detID", &detID, &b_det_ID); - fChain->SetBranchAddress("e", &e, &b_energy); - fChain->SetBranchAddress("t", &t, &b_time_stamp); + fChain->SetBranchAddress("id", &ID, &b_ID); + fChain->SetBranchAddress("e", &e, &b_energy); + fChain->SetBranchAddress("t", &t, &b_timestamp); } Bool_t Monitor_raw::Notify() diff --git a/armory/evt2hist.cpp b/armory/evt2hist.cpp index 06a3069..5273212 100644 --- a/armory/evt2hist.cpp +++ b/armory/evt2hist.cpp @@ -168,10 +168,11 @@ int main(int argn, char **argv) { ///} ///data.Print(); - int detID = map[data.id]; + int detID = mapping[data.id]; if( 0 <= detID && detID < 100 ){ if( corrFile != ""){ - double eCal = ApplyCorrection(corr, detID, data.energy); + //double eCal = ApplyCorrection(corr, detID, data.energy); + double eCal = corr[detID][0] + corr[detID][1]*data.energy; he[detID]->Fill(eCal); }else{ he[detID]->Fill(data.energy); @@ -230,7 +231,7 @@ int main(int argn, char **argv) { fclose(inFile); printf("\n============= reasched end of file\n"); - printf("Crtl+C to end program.\n"); + printf("\nCrtl+C to end program.\n"); app->Run(); diff --git a/armory/evt2root.cpp b/armory/evt2root.cpp index e761ef7..e7b3126 100644 --- a/armory/evt2root.cpp +++ b/armory/evt2root.cpp @@ -63,10 +63,9 @@ public: //############################################# int main(int argn, char **argv) { - if (argn != 2 && argn != 3 ) { + if (argn != 2 ) { printf("Usage :\n"); - printf("%s [evt File] [timeWindow] \n", argv[0]); - printf(" timeWindow : number of tick, 1 tick = 10 ns. default = 100 \n"); + printf("%s [evt File] \n", argv[0]); return 1; } @@ -75,7 +74,6 @@ int main(int argn, char **argv) { outFileName.Remove(inFileName.First('.')); outFileName.Append("_raw.root"); - long int inFilePos = 0; TBenchmark gClock; gClock.Reset(); gClock.Start("timer"); @@ -96,6 +94,7 @@ int main(int argn, char **argv) { fseek(inFile, 0L, SEEK_END); long int inFileSize = ftell(inFile); rewind(inFile); ///back to the File begining + long int inFilePos = 0; printf(" in file: %s\n", inFileName.Data()); printf("out file: %s\n", outFileName.Data()); @@ -106,21 +105,20 @@ int main(int argn, char **argv) { TTree * tree = new TTree("tree", "tree"); tree->Branch("evID", &measureID, "data_ID/L"); - tree->Branch("detID", &data.id, "det_ID/s"); + tree->Branch("id", &data.id, "ID/s"); tree->Branch("e", &data.energy, "energy/s"); - tree->Branch("t", &data.time, "time_stamp/l"); + tree->Branch("t", &data.time, "timestamp/l"); //tree->Branch("tdiff", &data.timeDiff, "time_Diff/L"); -//=======TODO online event building unsigned int header[4]; //read 4 header, unsigned int = 4 byte = 32 bits. unsigned long long nWords = 0; - ULong64_t timeLast = 0; + //ULong64_t timeLast = 0; //=============== Read File /// while ( ! feof(inFile) ){ - while ( inFilePos <= inFileSize ){ + while( inFilePos < inFileSize || feof(inFile) ){ fread(header, sizeof(header), 1, inFile); inFilePos = ftell(inFile); @@ -166,13 +164,12 @@ int main(int argn, char **argv) { //event stats, print status every 10000 events if ( measureID % 10000 == 0 ) { - inFilePos = ftell(inFile); float tempf = (float)inFileSize/(1024.*1024.*1024.); gClock.Stop("timer"); double time = gClock.GetRealTime("timer"); gClock.Start("timer"); printf("Total measurements: \x1B[32m%lld \x1B[0m\nPercent Complete: \x1B[32m%ld%% of %.3f GB\x1B[0m\nTime used:%3.0f min %5.2f sec\033[A\033[A\r", - measureID, (100*inFilePos/inFileSize), tempf, TMath::Floor(time/60.), time - TMath::Floor(time/60.)*60.); + measureID +1 , (100*inFilePos/inFileSize), tempf, TMath::Floor(time/60.), time - TMath::Floor(time/60.)*60.); } //cern fill tree @@ -187,7 +184,7 @@ int main(int argn, char **argv) { gClock.Start("timer"); float tempf = (float)inFileSize/(1024.*1024.*1024.); printf("Total measurements: \x1B[32m%lld \x1B[0m\nPercent Complete: \x1B[32m%ld%% of %.3f GB\x1B[0m\nTime used:%3.0f min %5.2f sec\033[A\r", - measureID, (100*inFilePos/inFileSize), tempf, TMath::Floor(time/60.), time - TMath::Floor(time/60.)*60.); + measureID+1, (100*inFilePos/inFileSize), tempf, TMath::Floor(time/60.), time - TMath::Floor(time/60.)*60.); fclose(inFile); diff --git a/armory/pixie2root.cpp b/armory/pixie2root.cpp index bfd3d17..1f6c20f 100644 --- a/armory/pixie2root.cpp +++ b/armory/pixie2root.cpp @@ -307,7 +307,7 @@ int main(int argc, char **argv) { ///========== need a mapping, can reduce the array size, speed up. - int ch = map[subevt[sevtmult].id]; + int ch = mapping[subevt[sevtmult].id]; if ( 0 <= ch && ch < NCRYSTAL ){ energy[ch] = subevt[sevtmult].energy; etimestamp[ch] = subevt[sevtmult].time; diff --git a/mapping.h b/mapping.h index 8cae591..d7f4199 100644 --- a/mapping.h +++ b/mapping.h @@ -12,7 +12,7 @@ Other : 200 - 299 #define NOTHER 52 // 0 1 2 3 4 5 6 7 8 9 -int map[130] = { 0, 1, 2, 3, 100, 4, 5, 6, 7, 101, // 0 +int mapping[130] = { 0, 1, 2, 3, 100, 4, 5, 6, 7, 101, // 0 8, 9, 10, 11, 102, -1, 12, 13, 14, 15, // 10 103, 16, 17, 18, 19, 104, -1, -1, -1, -1, // 20 -1, -1, 20, 21, 22, 23, 105, 24, 25, 26, // 30