// This code does re-ordering and adds TDC information to trees, like in previous scripts. However, in this case it will take Hit and Event trees from separate files. The event tree will contain a "cluster" from every single strip. The script will write a new Hit tree which contains the hits from all 5 time bins for whichever strip had the biggest hit during the BAT trigger. // THIS SCRIPT SHOULD BE COMPILED BEFORE USE #include #include "TTree.h" #include "TFile.h" void sortPulse(char eventFile[],char hitFile[],int WindowDelay) { cout << "==== Event input file: " << eventFile << endl; cout << "==== Input file with all strip data: " << hitFile << endl; /// Code to sort our trees so that the 5 samples taken after each hit on the BAT is sorted into the right order. /// Want to modify so that it runs on all the trees in a file, which we give by name. TFile e(eventFile); TFile h(hitFile); //open new file to store the output Char_t outFile[100]; sprintf(outFile, "Pulse%s", eventFile); // creates an output filename TFile f2(outFile,"recreate"); cout << "==== Output file: " << outFile << endl; // Set up the information about the sensors int num_sensors=1; int num_samples=5; int num_skip=1100; TString sensorName[num_sensors]; // Specify the names of the sensors here sensorName[0] = ""; // Also, set up the bad strip array. This is initialised to all zeroes. int sensor=0; // initialise variable to keep track of sensor we're working on int strip=0; // keeps track of individual noisy strip entries int num_noisystrips=7; // max no. of noisy strips on a sensor - choose this yourself int noisyStrips[num_sensors][num_noisystrips]; // Loop to initialise entries in array to zero for(sensor = 0; sensor < num_sensors ; sensor++){ for(strip = 0; strip < num_noisystrips ; strip++){ noisyStrips[sensor][strip]=0; } } // Then, we can specify specific noisy strips. Need to only choose the relevant sensor! /* // Gla2 - no proper noisy strips, but due to an oddity in the pattern of unbonded strips I want to ignore strip 3. noisyStrips[0][0]=3; noisyStrips[0][1]=127; */ // Gla1 - avoid a single dubious strip near the start (bond pattern) and also last strip in each Alink (possibly shares signal with next strip) noisyStrips[0][0]=387; noisyStrips[0][1]=511; // Last strip in each Alink is dodgy too noisyStrips[0][2]=415; noisyStrips[0][3]=447; noisyStrips[0][4]=479; // On 3D sensor, finally want to avoid a strip with no proper signal, plus the adjacent one noisyStrips[0][5]=450; noisyStrips[0][6]=451; /* // Gla0 noisyStrips[0][0]=259; noisyStrips[0][1]=383; // Fre2 // The strips are around 1173 and 205, which are no use after re-ordering (fall between header strips) and also 1153, which is noisy. Note that 1152 is not used for finding hits anyway (edge of array). // So, block 1172-74, 1204-06, 1153 noisyStrips[0][0]=1172; noisyStrips[0][1]=1173; noisyStrips[0][2]=1174; noisyStrips[0][3]=1204; noisyStrips[0][4]=1205; noisyStrips[0][5]=1206; noisyStrips[0][6]=1153; // On Fre0, we have a few noisy strips - 1467, 1469, 1471. Also have strips falling between header strips as above - 1416, 1448. Lastly, 1460 is unbonded. // So, block 1415-7, 1447-49, 1459-61, 1467-71 noisyStrips[0][0]=1415; noisyStrips[0][1]=1416; noisyStrips[0][2]=1417; noisyStrips[0][3]=1447; noisyStrips[0][4]=1448; noisyStrips[0][5]=1449; noisyStrips[0][6]=1459; noisyStrips[0][7]=1460; noisyStrips[0][8]=1461; noisyStrips[0][9]=1467; noisyStrips[0][10]=1468; noisyStrips[0][11]=1469; noisyStrips[0][12]=1470; noisyStrips[0][13]=1471; */ // Define variables used when adding new data to output file int OrderTrue = 0; int OrderTrueEvent = 0; int SampleTime = 0; int BadStrip=0; int iHit=0; int iSens=0; Int_t nentriesEvent; Int_t nentriesHit; Int_t *indexEvent; int BCntLoopEvent=0; int BCntOldEvent=0; // Used to check for consecutive values int BATTrigEvent=0; int BCntLoopHit=0; int ADCMaxEvent=0; int PosXMaxEvent=0; int BATTrigHit=0; int PosXHit=0; // The trees: TTree *treeEvent; TTree *treeHit; TTree *tsortedEvent; TTree *tsortedHit; // Strings I'll use TString pathEvent; TString pathHit; // First loop will loop over all the sensors, saving us from having to use new code for each sensor. uu keeps track. for(Int_t uu = 0; uu < num_sensors ; uu++) { //******************************* // Setup for each sensor //******************************* // Getting the relevant tree within the file. Here, I want to get the Event AND Hit trees. Need to use pathEvent = "VeloLCMSMoni/LCMSEvent" + sensorName[uu]; // Put together the path treeEvent = (TTree*)e.Get(pathEvent); nentriesEvent = (Int_t)treeEvent->GetEntries(); pathHit = "VeloLCMSMoni/LCMSHit" + sensorName[uu]; treeHit = (TTree*)h.Get(pathHit); nentriesHit = (Int_t)treeHit->GetEntries(); //Drawing variable BCntLoop with no graphics option, for the event tree. //Array V1 (see TTree::Draw) is then an array of these BCntLoop entries treeEvent->Draw("BCntLoop","","goff"); indexEvent = new Int_t[nentriesEvent]; //Sorts the array containing BCntLoop values (V1) and puts the sorted indices into indexEvent //"false" options sets increasing order sorting TMath::Sort(nentriesEvent,treeEvent->GetV1(),indexEvent,false); //Create an empty clone of the original tree //I want to be able to add extra entries to this cloned tree // Copy the existing tree structure, call the tree *tsorted tsortedEvent = (TTree*)treeEvent->CloneTree(0); // Add new branches for the true ordering of events in each batch of 5 (OrderTrue), the raw timing info (TimeDiff) and the calculated hit time in ns HitTime // Note that I'm trying to add the same info to a couple of different trees here. Have to be careful. tsortedEvent->Branch("OrderTrue",&OrderTrue,"OrderTrue/I"); tsortedEvent->Branch("WindowDelay",&WindowDelay,"WindowDelay/I"); tsortedEvent->Branch("SampleTime",&SampleTime,"SampleTime/I"); //****** // Next, we need to do something similar for the Hit tree. // MAJOR POINT - Attempting to sort the Hit tree leads to failure, probably due to too many entries. However, since all the entries in the hit tree corresponding to a single trigger do appear a distinct block, and these blocks are in the correct order, we don't need to do the ordering. //Drawing variable BCntLoop with no graphics option, for the hit tree. //Array V1 (see TTree::Draw) is then an array of these BCntLoop entries //treeHit->Draw("BCntLoop","","goff"); //Then, we create an array, with the same no. of entries as the Event tree //indexHit = new Int_t[nentriesHit]; //Sorts the array containing BCntLoop values (V1) and puts the sorted indices into indexHit //"false" options sets increasing order sorting //TMath::Sort(nentriesHit,treeHit->GetV1(),indexHit,false); //Create an empty clone of the original tree //I want to be able to add extra entries to this cloned tree // Copy the existing tree structure, call the tree *tsorted tsortedHit = (TTree*)treeHit->CloneTree(0); // Add new branches for the true ordering of events in each batch of 5 (OrderTrue), the raw timing info (TimeDiff) and the calculated hit time in ns HitTime tsortedHit->Branch("OrderTrue",&OrderTrue,"OrderTrue/I"); tsortedHit->Branch("BadStrip",&BadStrip,"BadStrip/I"); tsortedHit->Branch("WindowDelay",&WindowDelay,"WindowDelay/I"); tsortedHit->Branch("SampleTime",&SampleTime,"SampleTime/I"); // Have to specify the following so that we can look at specific variables after using GetEntry. Also have to be careful when dealing with different detector trees - use different variables for different detectors. I get the BCntLoop to check that there are no errors in the file. BCntLoopEvent=0; BCntOldEvent=0; // Used to check for consecutive values BATTrigEvent=0; ADCMaxEvent=0; PosXMaxEvent=0; treeEvent->SetBranchAddress("BCntLoop",&BCntLoopEvent); treeEvent->SetBranchAddress("BATTrig",&BATTrigEvent); treeEvent->SetBranchAddress("ADCMax",&ADCMaxEvent); treeEvent->SetBranchAddress("PosXMax",&PosXMaxEvent); BCntLoopHit=0; BATTrigHit=0; PosXHit=0; treeHit->SetBranchAddress("BCntLoop",&BCntLoopHit); treeHit->SetBranchAddress("BATTrig",&BATTrigHit); treeHit->SetBranchAddress("PosX",&PosXHit); //*************** // // LOOPING OVER ALL THE ENTRIES IN THE EVENT TREE // Need to use iSens to loop over the sensors. // Within this loop, we use the array of indices to grab the events in the correct order, writing them to a new tree. cout << "About to start loop" << endl; iHit=0; // iterator to keep track of our position within the Hit tree for (iSens=0;iSensGetEntry(indexEvent[iSens]); // Read in an event. Previous statements mean BCntLoop is accessible. // Check to make sure that the trigger numbering is working correctly by comparing the BATTrig in the Vetra file, eventNum in the TDC if (BATTrigEvent!=((iSens+num_skip)/num_samples)){ cout << "Event numbering problem at BATTrig=" << BATTrigEvent << " entry number" << iSens << endl; } // Apply the numbering to the "num_samples" time bins OrderTrue = iSens%num_samples; OrderTrueEvent = iSens%num_samples; // Check the bunch counters are consecutive (if OrderTrue>0) if ((OrderTrue>0)&&((BCntLoopEvent-BCntOldEvent) != 1)){ cout << "Non-consecutive entries at BCntLoopEvent=" << BCntLoopEvent << " BAT Trigger no=" << BATTrigEvent << " OrderTrue=" << OrderTrue << endl; } // cout << "Value of BCntLoopEvent is " << BCntLoopEvent << endl; // Calculate the SampleTime, based on the correct ordering and the delay of the acceptance window SampleTime=25*OrderTrue+(20-WindowDelay); BCntOldEvent=BCntLoopEvent; // Update to the old bunch count value tsortedEvent->Fill(); // Fill the Event tree //****************** // Next, deal with the Hit tree. // We want to look at the Hit tree every time that the SampleTime falls within the range .... // The sorting by BCntLoop will effectively order the hit tree first by BATTrig then by bunch counter // So, when we look at the tree we will loop over all the events with the same BATTrig // Then, select those hits where the hit location matches PosXMax from the Event tree // Write these hits to the new hit tree. This should give us all 5 samples associated with each hit. // Have to be careful about messing up the SampleTime variable // // if((SampleTime>=45)&&(SampleTime<=65)){ while (iHitGetEntry(iHit); // Using iHit to keep track of indices, read in an event from the sorted Hit tree. Don't have to think about ordering // As a test, check to see if the BATTrig for the current hit is lower than that for event. If this is true, something has gone wrong. //if (BATTrigEvent>BATTrigHit) cout << "Error with Hit tree ordering - missed an entry somewhere " << endl; //cout << "iSens=" << iSens << " iHit=" << iHit << " BATTrigEvent=" << BATTrigEvent << " BATTrigHit=" << BATTrigHit << " BCntLoopEvent=" << BCntLoopEvent << " BCntLoopHit=" << BCntLoopHit << endl; // Check to see if the BATTrig number in this Hit match the current BATTrig. If not, break. if (BATTrigEvent!=BATTrigHit) break; //cout << "About to check if this is the correct strip at strip no " << PosXHit << " and BATTrig " << BATTrigHit << endl; // Since we have not broken, this Hit must belong to the current BAT trigger. So, want to select the hit if it is in the correct position // We also add in a requirement about ADCMaxEvent, so we can choose signals of a certain polarity //if ((PosXMaxEvent==PosXHit)&&(ADCMaxEvent>0)) { if ((PosXMaxEvent==PosXHit)&&(ADCMaxEvent<0)) { // Write the hit to the tree, once we've calculated the OrderTrue value for this hit using BCntLoopEvent OrderTrue=BCntLoopHit-BCntLoopEvent+OrderTrueEvent; // If this falls outwith the correct 0-4 range, set it to 5 if ((OrderTrue<0)||(OrderTrue>4)) OrderTrue=5; // Then, calculate the correct SampleTime SampleTime=25*OrderTrue+(20-WindowDelay); // Lastly, we check for bad strips BadStrip=0; for(strip = 0; strip < num_noisystrips ; strip++){ if (noisyStrips[sensor][strip]==PosXHit) BadStrip=1; } tsortedHit->Fill(); // Fill the Hit tree //cout << "Adding hit tree entry at BATTrig " << BATTrigEvent << " and order " << OrderTrue << endl; } ++iHit; // Move on to the next entry in the sorted hit tree. Note that if the BCntLoop comparison fails, then we'll stick with the same event till the next one. This should work fine, but be on the lookout for bugs! } } // SampleTime test loop } // Loop over all values of iSens // write the histos once this is finished tsortedEvent->Write(); tsortedHit->Write(); delete [] indexEvent; cout << sensorName[uu] <<" done" << endl; } // End loop over all sensors // Write some output info cout << "No of entries in detector file: " << nentriesEvent << endl; cout << "Equivalent no of telescope triggers: " << (nentriesEvent)/5+220 << endl; }