/**********************************************************/ /* */ /* Modified by Ryan From */ /* */ /* PXI SCAN CODE -- J.M. Allmond (ORNL) -- July 2016 */ /* */ /**********************************************************/ #include #include #include #include #include #include "TFile.h" #include "TTree.h" #include "TMath.h" #include "TBenchmark.h" #define RAND ((float) rand() / ((unsigned int) RAND_MAX + 1)) // random number in interval (0,1) #define MAX_CRATES 2 #define MAX_BOARDS_PER_CRATE 13 #define MAX_CHANNELS_PER_BOARD 16 #define BOARD_START 2 #define MAX_ID MAX_CRATES*MAX_BOARDS_PER_CRATE*MAX_CHANNELS_PER_BOARD #define HEADER_LENGTH 4 //unit = words with 4 bytes per word #define MAX_SUB_LENGTH 2016 //unit = words with 4 bytes per word ; 2004 --> 40 micro-second trace + 4 word header #define RAWE_REBIN_FACTOR 2.0 // Rebin 32k pixie16 spectra to something smaller to fit better into 8k. #include "../mapping.h" #include "../armory/evtReader.h" unsigned long long int dataCount=0; unsigned long long int pileUpCount=0; unsigned long long int evtCount=0; int tick2ns = 4; /////////////////////////////////// // START OF MAIN FUNCTION // /////////////////////////////////// int main(int argc, char **argv) { printf("=====================================\n"); printf("=== evt.to --> root ===\n"); printf("=====================================\n"); // Check that the corrent number of arguments were provided. if (argc <= 3) { printf("Incorrect number of arguments:\n"); printf("%s [outFile] [timeWindow] [to File1] [to File2] .... \n", argv[0]); printf(" outFile : output root file name\n"); printf(" timeWindow : number of tick, 1 tick = %d ns. default = 100 \n", tick2ns); return 1; } TString outFileName = argv[1]; int timeWindow = atoi(argv[2]); int nFile = argc - 3; TString inFileName[nFile]; for( int i = 0 ; i < nFile ; i++){ inFileName[i] = argv[i+3]; } printf("====================================\n"); evtReader * evt = new evtReader(); DataBlock * data = evt->data; printf(" Number of input file : %d \n", nFile); printf(" out file : \033[1;31m%s\033[m\n", outFileName.Data()); printf(" Event building time window : %d tics = %d nsec \n", timeWindow, timeWindow*tick2ns); TFile * outRootFile = new TFile(outFileName, "recreate"); outRootFile->cd(); TTree * tree = new TTree("tree", outFileName); unsigned long long evID = 0; int multi = 0; int id[MAX_ID] = {0}; double e[MAX_ID] = {TMath::QuietNaN()}; unsigned long long e_t[MAX_ID] = {0}; int cfd[MAX_ID] = {0}; bool pileup[MAX_ID] = {0}; int qdc[MAX_ID][8] = {0}; int multiClover = 0 ; /// this is total multiplicity for all crystal int runID = 0; // date-run-fileID int multiBeam = 0; //unsigned short pileup[MAXMULTI]; tree->Branch("evID", &evID, "event_ID/l"); tree->Branch("multi", &multi, "multi/I"); tree->Branch("detID", id, "detID[multi]/I"); tree->Branch("e", e, "e[multi]/D"); tree->Branch("e_t", e_t, "e_timestamp[multi]/l"); //tree->Branch("pileup", pileup, "pileup[multi]/O"); tree->Branch("cfd", cfd, "cfd[multi]/I"); tree->Branch("qdc", qdc, "qdc[multi][8]/I"); tree->Branch("multiClover", &multiClover, "multiplicity_Clover/I"); tree->Branch("multiBeam", &multiBeam, "multiplicity_Beam/I"); tree->Branch("runID", &runID, "runID/I"); int countGP = 0; //gamma-particle coincident double totalDataSize = 0; int outFileCount = 0; for( int i = 0; i < nFile; i++){ evt->OpenFile(inFileName[i], false); if( evt->IsOpen() == false ) continue; printf("==================================================== %d / %d\n", i+1, nFile); printf("\033[1;31m%s \033[m\n", inFileName[i].Data()); int pos = inFileName[i].Last('/'); TString temp = inFileName[i]; temp.Remove(0, pos+1); pos = temp.Last('-'); temp.Remove(0, pos+1); pos = temp.Last('.'); temp.Remove(pos); runID = atoi(temp.Data()); long long int etime = -1; long long int tdif = -1; while ( evt->IsEndOfFile() == false ) { //main while loop if ( evt->ReadBlock() == -1) break; if ( data->crate == 0 ) continue; //Set reference time for event building if (etime == -1) { etime = data->time; tdif = 0; multi = 0; multiClover = 0; multiBeam = 0; }else { tdif = data->time - etime; } //Check for end of event, rewind, and break out of while loop if (tdif > timeWindow) { //Gate //if( multiClover > 2 && multiBeam ==0 ) { // outRootFile->cd(); tree->Fill(); //printf("------------------\n"); // // countGP++; //} evID ++; //clear data etime = data->time; tdif = 0; multi = 0; multiClover = 0; multiBeam = 0; int haha = data->crate*MAX_BOARDS_PER_CRATE*MAX_CHANNELS_PER_BOARD + (data->slot-BOARD_START)*MAX_CHANNELS_PER_BOARD + data->ch; id[multi] = mapping[haha]; e[multi] = data->energy; e_t[multi] = data->time; cfd[multi] = data->cfd_forced == 0 ? (data->cfd - 16384*(data->cfd_source == true ? 1 : 0))/2 : 0; pileup[multi] = data->pileup; for( int i = 0; i < 8; i++) qdc[multi][i] = data->QDCsum[i]; if( 0 <= id[multi] && id[multi] < 100 ) multiClover += 1; if( id[multi] == 300 || id[multi] == 301 || id[multi] == 307 || id[multi] == 308 ) multiBeam += 1; //printf("id: %d, multiBeam: %d, %llu, tdif %llu\n", id[multi], multiBeam, e_t[multi], tdif); multi++ ; }else{ //if within time window, fill array; int haha = data->crate*MAX_BOARDS_PER_CRATE*MAX_CHANNELS_PER_BOARD + (data->slot-BOARD_START)*MAX_CHANNELS_PER_BOARD + data->ch; id[multi] = mapping[haha]; e[multi] = data->energy; e_t[multi] = data->time; cfd[multi] = data->cfd_forced == 0 ? (data->cfd - 16384*(data->cfd_source == true ? 1 : 0))/2 : 0; //printf(" %d , %d , %d, %d \n", data->cfd_forced, data->cfd_source, data->cfd, cfd[multi]); pileup[multi] = data->pileup; for( int i = 0; i < 8; i++) qdc[multi][i] = data->QDCsum[i]; if( 0 <= id[multi] && id[multi] < 100 ) multiClover += 1; if( id[multi] == 300 || id[multi] == 301 || id[multi] == 307 || id[multi] == 308 ) multiBeam += 1; //printf("id: %d, multiBeam: %d, %llu, tdif %llu\n", id[multi], multiBeam, e_t[multi], tdif); multi++ ; } // total pileups if (data->pileup == 1) { pileUpCount++; } evt->PrintStatus(10000); } // end main while loop evt->PrintStatus(1); printf("\n\n\n"); printf(" total number of event built : %llu\n", evID); //printf(" total number of Gamma - GAGG coincdient : %d (%.3f %%)\n", countGP, countGP*1.0/evID); outRootFile->cd(); tree->Write(); totalDataSize += (evt->GetFileSize())/1024./1024./1024.; double rootFileSize = outRootFile->GetSize()/1024./1024./1024. ; // in GB printf(" ----------- root file size : %.3f GB\n", rootFileSize); printf(" ---------- total read size : %.3f GB\n", totalDataSize); printf(" ----------- reduction rate : %.3f %%\n", rootFileSize*100./totalDataSize); evt->CloseFile(); /* if( rootFileSize > 3.0 ) { break; }*/ ///try to open a new root file when file size > 2 GB /*if( rootFileSize > 2.0 ) { outRootFile->Close(); delete outRootFile; delete tree; outFileCount += 1; if( outFileCount > 5 ) break; TString outFileName2 = outFileName; outFileName2.Insert(outFileName.Sizeof() - 6, Form("_%03d",outFileCount)); outRootFile = new TFile( outFileName2, "recreate"); tree = new TTree("tree", "tree"); tree->Branch("evID", &evID, "event_ID/l"); tree->Branch("multi", &multi, "multi/I"); tree->Branch("detID", id, "detID[multi]/I"); tree->Branch("e", e, "e[multi]/D"); tree->Branch("e_t", e_t, "e_timestamp[multi]/l"); tree->Branch("qdc", qdc, "qdc[multi][8]/I"); tree->Branch("multiClover", &multiClover, "multiplicity_crystal/I"); tree->Branch("multiBeam", &multiBeam, "multiplicity_GAGG/I"); tree->Branch("runID", &runID, "runID/I"); }*/ } outRootFile->Close(); printf("\n\n\n==================== finished.\r\n"); //printf(" number of Gamma - GAGG coincdient : %d\n", countGP); return 0; }