diff --git a/armory/evt2hist.cpp b/armory/evt2hist.cpp index 5273212..c6764b4 100644 --- a/armory/evt2hist.cpp +++ b/armory/evt2hist.cpp @@ -23,6 +23,12 @@ #define MAX_CHANNELS_PER_BOARD 16 #define BOARD_START 2 +//#############################TODO +// 1) Get ADC data +// 2) Change to GUI +// 3) calibration gamma + +int rawEnergyThreshold = 100; class measurment{ @@ -39,14 +45,29 @@ public: UShort_t trace_length; Bool_t trace_out_of_range; - UShort_t id; + Int_t trailing; + Int_t leading; + Int_t gap; + Int_t baseline; + Int_t ADC[8]; - measurment(){}; + UShort_t id; + Int_t detID; + ULong64_t eventID; + + UShort_t trace[1024]; + + measurment(){ + Clear(); + }; + + ~measurment(){}; void Clear(){ ch = 0; slot = 0; crate = 0; + headerLength = 0; eventLength = 0; pileup = false; time = 0; @@ -55,12 +76,35 @@ public: trace_length = 0; trace_out_of_range = 0; id = 0; + detID = -1; + eventID = 0; + trailing = 0; + leading = 0; + gap = 0; + baseline = 0; + for( int i = 0; i < 8; i++) ADC[i] = -1; + for( int i = 0; i < 1024; i++) trace[i] = 0; } void Print(){ - printf("Crate: %d, Slot: %d, Ch: %d | id: %d \n", crate, slot, ch, id); + printf("============== eventID : %llu\n", eventID); + printf("Crate: %d, Slot: %d, Ch: %d | id: %d = detID : %d \n", crate, slot, ch, id, detID); printf("HeaderLength: %d, Event Length: %d, energy: %d, timeStamp: %llu\n", headerLength, eventLength, energy, time); printf("trace_length: %d, pile-up:%d\n", trace_length, pileup); + if( headerLength > 4 ){ + if( headerLength > 12 ){ + printf(" trailing : %d\n", trailing); + printf(" leading : %d\n", leading); + printf(" gap : %d\n", gap); + printf(" baseLine : %d\n", baseline); + } + printf(" ADC : \n"); + for( int i = 0; i < 8; i++) printf(" %-10d\n", ADC[i]); + } + if( eventLength > headerLength ){ + printf(" trace:\n"); + for( int i = 0 ; i < trace_length ; i++)printf("%3d| %-10d\n",i, trace[i]); + } } }; @@ -80,11 +124,11 @@ int main(int argn, char **argv) { TString inFileName = argv[1]; TString corrFile = ""; - std::vector> corr; + std::vector> eCorr; if( argn >= 3 ){ corrFile = argv[2]; - corr.clear(); - corr = LoadCorrectionParameters(corrFile); + eCorr.clear(); + eCorr = LoadCorrectionParameters(corrFile); } long int inFilePos; @@ -92,7 +136,7 @@ int main(int argn, char **argv) { gClock.Reset(); gClock.Start("timer"); - ULong64_t measureID = 0; + ULong64_t measureID = -1; measurment data; @@ -104,19 +148,21 @@ int main(int argn, char **argv) { return -404; } - //get file size - //fseek(inFile, 0L, SEEK_END); - //long int inFileSize = ftell(inFile); - //rewind(inFile); ///back to the File begining - - printf(" in file: %s\n", inFileName.Data()); + printf(" in file: \003[1;31m%s\033[m\n", inFileName.Data()); printf(" Gamma energy correction file : %s\n", corrFile == "" ? "Not provided." : corrFile.Data()); printf("--------------------------------\n"); + + //================ get file size + fseek(inFile, 0L, SEEK_END); + long int inFileSize = ftell(inFile); + rewind(inFile); ///back to the File begining + unsigned long long fpos = 0; + //================ Historgrams TH1F * he[NCRYSTAL]; for( int i = 0 ; i < NCRYSTAL; i++){ - he[i] = new TH1F(Form("he%02d", i), Form("e-%2d", i), 1000, 0, 6000); + he[i] = new TH1F(Form("he%02d", i), Form("e-%2d", i), 1000, 0, 2000); switch (i % 4){ case 0 : he[i]->SetLineColor(2); break; case 1 : he[i]->SetLineColor(4); break; @@ -127,26 +173,29 @@ int main(int argn, char **argv) { //================ Set Canvas TApplication * app = new TApplication ("app", &argn, argv); + + TCanvas * canvas = new TCanvas("fCanvas", "Online Spectrum", 1800, 1400); canvas->Divide(1, 9, 0); canvas->SetCrosshair(1); for( int i = 0; i < 9 ; i++){ - canvas->cd(i+1)->SetBottomMargin(0.06); + canvas->cd(i+1)->SetBottomMargin(0.1); canvas->cd(i+1)->SetRightMargin(0.002); } - - unsigned int header[4]; //read 4 header, unsigned int = 4 byte = 32 bits. - unsigned long long nWords = 0; - unsigned long long fpos = 0; + + //=============== Read File + unsigned int header[4]; //read 4 header, unsigned int = 4 byte = 32 bits. + while ( ! feof(inFile) ){ fread(header, sizeof(header), 1, inFile); measureID ++; /// see the Pixie-16 user manual, Table4-2 - data.ch = header[0] & 0xF ; + data.eventID = measureID; + data.ch = header[0] & 0xF ; data.slot = (header[0] >> 4) & 0xF; data.crate = (header[0] >> 8) & 0xF; data.headerLength = (header[0] >> 12) & 0x1F; @@ -154,60 +203,95 @@ int main(int argn, char **argv) { data.pileup = header[0] >> 31 ; data.time = ((ULong64_t)(header[2] & 0xFFFF) << 32) + header[1]; data.cfd = header[2] >> 16 ; - data.energy = header[3] & 0xFFFF; + data.energy = (header[3] & 0xFFFF ) /2; // I don;t know why it has to "rebin" data.trace_length = (header[3] >> 16) & 0x7FFF; data.trace_out_of_range = header[3] >> 31; data.id = data.crate*MAX_BOARDS_PER_CRATE*MAX_CHANNELS_PER_BOARD + (data.slot-BOARD_START)*MAX_CHANNELS_PER_BOARD + data.ch; + data.detID = mapping[data.id]; - nWords += data.eventLength; - - ///printf("----------------------nWords: %llu, fpos: %llu\n", nWords, fpos); - ///for(int i = 0; i < 4; i++){ - /// printf(" %x\n", header[i]); - ///} - ///data.Print(); - - int detID = mapping[data.id]; - if( 0 <= detID && detID < 100 ){ - if( corrFile != ""){ - //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); + ///======== read ADC + if( data.headerLength >= 4 ){ + unsigned int extraHeader[data.headerLength-4]; + fread(extraHeader, sizeof(extraHeader), 1, inFile); + //if( measureID < 10 ) for(int i = 0; i < data.headerLength - 4; i++) printf(" %x\n", extraHeader[i]); + if( data.headerLength > 12){ + data.trailing = extraHeader[0]; + data.leading = extraHeader[1]; + data.gap = extraHeader[2]; + data.baseline = extraHeader[3]; + } + + for( int i = 0; i < 8; i++){ + int startID = 0; + if( data.headerLength > 12) startID = 4; ///the 1st 4 words + data.ADC[i] = extraHeader[i+startID]; + } + } + ///====== read trace + if( data.eventLength > data.headerLength ){ + unsigned int traceBlock[data.trace_length / 2]; + fread(traceBlock, sizeof(traceBlock), 1, inFile); + + for( int i = 0; i < data.trace_length/2 ; i++){ + data.trace[2*i+0] = traceBlock[i] & 0xFFFF ; + data.trace[2*i+1] = (traceBlock[i] >> 16 ) & 0xFFFF ; } } - //=== jump to next measurement - if( data.eventLength > 4 ){ - fseek(inFile, sizeof(int) * (data.eventLength-4), SEEK_CUR); - fpos += ftell(inFile); + if( measureID < 10 ) { + printf("----------------------event Length: %u, fpos: %llu byte (%lld words)\n", data.eventLength, fpos, fpos/4); + for(int i = 0; i < 4; i++) printf(" %x\n", header[i]); + data.Print(); } + //=== jump to next measurement. obsolete, we read the whole block + ///if( data.eventLength > 4 ) { + /// if( fseek(inFile, sizeof(int) * (data.eventLength-data.headerLength), SEEK_CUR) != 0 ) break;; + ///} + fpos = ftell(inFile); + + + //==== Fill Histogram + if( 0 <= data.detID && data.detID < 100 && data.energy > rawEnergyThreshold ){ + if( corrFile != ""){ + ///========= apply correction + int order = (int) eCorr[data.detID].size(); + double eCal = 0; + for( int k = 0; k < order ; k++){ + eCal += eCorr[data.detID][k] * TMath::Power(data.energy, k); + } + he[data.detID]->Fill(eCal); + }else{ + he[data.detID]->Fill(data.energy); + } + } + + + //==== event stats, print status every 10000 events if ( measureID % 10000 == 0 ) { - //fseek(inFile, 0L, SEEK_END); - //inFileSize = ftell(inFile); - //fseek(inFile, fpos, SEEK_SET); + /// update file size, slow down? + fseek(inFile, 0L, SEEK_END); + inFileSize = ftell(inFile); + fseek(inFile, fpos, SEEK_SET); inFilePos = ftell(inFile); - float tempf = (float)inFilePos/(1024.*1024.*1024.); gClock.Stop("timer"); double time = gClock.GetRealTime("timer"); gClock.Start("timer"); - printf("Total measurements: \x1B[32m%llu \x1B[0m\nData Read: \x1B[32m %.3f GB\x1B[0m\nTime used:%3.0f min %5.2f sec\033[A\033[A\r", - measureID, tempf, TMath::Floor(time/60.), time - TMath::Floor(time/60.)*60.); + printf("Total measurements: \x1B[32m%llu \x1B[0m\nReading Pos: \x1B[32m %.3f/%.3f GB\x1B[0m\nTime used:%3.0f min %5.2f sec\033[A\033[A\r", + measureID, inFilePos/(1024.*1024.*1024.), inFileSize/1024./1024./1024, TMath::Floor(time/60.), time - TMath::Floor(time/60.)*60.); } - + //==== Plot Canvas gClock.Stop("timer"); int time = TMath::Floor(gClock.GetRealTime("timer")*1000); // in millisec gClock.Start("timer"); if( time % 1000 == 0 ){ for( int i = 0; i < NCRYSTAL; i++){ canvas->cd(i/4 +1); - canvas->cd(i/4 +1)->SetLogy(); + //canvas->cd(i/4 +1)->SetLogy(); if( i % 4 == 0 ) { he[i]->Draw(); }else{ @@ -218,15 +302,29 @@ int main(int argn, char **argv) { canvas->Update(); gSystem->ProcessEvents(); } + }//---- end of file loop + + + for( int i = 0; i < NCRYSTAL; i++){ + canvas->cd(i/4 +1); + //canvas->cd(i/4 +1)->SetLogy(); + if( i % 4 == 0 ) { + he[i]->Draw(); + }else{ + he[i]->Draw("same"); + } } + canvas->Modified(); + canvas->Update(); + gSystem->ProcessEvents(); + inFilePos = ftell(inFile); gClock.Stop("timer"); double time = gClock.GetRealTime("timer"); gClock.Start("timer"); - float tempf = (float)inFilePos/(1024.*1024.*1024.); - printf("Total measurements: \x1B[32m%llu \x1B[0m\nData Read: \x1B[32m %.3f GB\x1B[0m\nTime used:%3.0f min %5.2f sec\033[A\033[A\r", - measureID, tempf, TMath::Floor(time/60.), time - TMath::Floor(time/60.)*60.); + printf("Total measurements: \x1B[32m%llu \x1B[0m\nReading Pos: \x1B[32m %.3f/%.3f GB\x1B[0m\nTime used:%3.0f min %5.2f sec\033[A\033[A\r", + measureID, inFilePos/(1024.*1024.*1024.), inFileSize/1024./1024./1024, TMath::Floor(time/60.), time - TMath::Floor(time/60.)*60.); fclose(inFile); diff --git a/armory/makefile b/armory/makefile index bba2edb..d3e03af 100644 --- a/armory/makefile +++ b/armory/makefile @@ -4,7 +4,8 @@ CC=g++ #all: xia2root xia2ev2_nopart pixie2root scan evt2root evt2hist #all: xia2root to2root MergeEVT evt2hist pxi-time-order ev22txt EventBuilder #all: xia2root to2root MergeEVT pxi-time-order ev22txt EventBuilder -all: to2root MergeEVT ev22txt EventBuilder pxi-time-order +#all: to2root MergeEVT ev22txt EventBuilder pxi-time-order +all: to2root evt2hist MergeEVT ev22txt EventBuilder pxi-time-order #this is FSU evt to root xia2root: ../armory/xia2root.cpp