Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-12-16 09:28:16

0001 #pragma once
0002 #ifndef COMMONHELPERFUNCTIONS_H
0003 #define COMMONHELPERFUNCTIONS_H
0004 
0005   #include <iostream>
0006   #include <fstream>
0007   #include "TString.h"
0008   #include "TObjString.h"
0009   #include <vector>
0010   #include <map>
0011   #include <utility>
0012 
0013   struct Layer{
0014     Layer(): nCells(0), energy(0.), avX(0.), avY(0.) {}
0015     double nCells;
0016     double energy;
0017     double avX;
0018     double avY;
0019   } ;
0020 
0021   inline int GetMaxLayer(std::map<int,Layer> layers){
0022     int maxLayer      = -1;
0023     double maxELayer  = 0;
0024     std::map<int, Layer>::iterator ithLayer;
0025     for(ithLayer=layers.begin(); ithLayer!=layers.end(); ++ithLayer){
0026       if (maxELayer < ithLayer->second.energy ){
0027         maxELayer = ithLayer->second.energy;
0028         maxLayer  = ithLayer->first;
0029       }
0030     }
0031     return maxLayer;
0032   }
0033   inline double GetAverageLayer(std::map<int,Layer> layers){
0034     double avLayer    = 0;
0035     double totE  = 0;
0036     std::map<int, Layer>::iterator ithLayer;
0037     for(ithLayer=layers.begin(); ithLayer!=layers.end(); ++ithLayer){
0038       avLayer +=   ithLayer->first*ithLayer->second.energy;
0039       totE    += ithLayer->second.energy;
0040     }
0041     avLayer = avLayer/totE;
0042     return avLayer;
0043   }
0044   
0045   inline double GetXAverage(std::map<int,Layer> layers, int layerMax = -100){
0046     double avLayer    = 0;
0047     double totE  = 0;
0048     std::map<int, Layer>::iterator ithLayer;
0049     for(ithLayer=layers.begin(); ithLayer!=layers.end(); ++ithLayer){
0050       if ((layerMax != -100) && (ithLayer->first > layerMax) )
0051         continue;
0052       avLayer += ithLayer->second.avX*ithLayer->second.energy;
0053       totE    += ithLayer->second.energy;
0054     }
0055     avLayer = avLayer/totE;
0056     return avLayer;
0057   }
0058 
0059   inline double GetYAverage(std::map<int,Layer> layers, int layerMax = -100){
0060     double avLayer    = 0;
0061     double totE  = 0;
0062     std::map<int, Layer>::iterator ithLayer;
0063     for(ithLayer=layers.begin(); ithLayer!=layers.end(); ++ithLayer){
0064       if ( (layerMax !=-100) && (ithLayer->first > layerMax) )
0065         continue;
0066       avLayer += ithLayer->second.avY*ithLayer->second.energy;
0067       totE    += ithLayer->second.energy;
0068     }
0069     avLayer = avLayer/totE;
0070     return avLayer;
0071   }
0072   
0073   struct RunInfo{
0074     RunInfo(): runNr(0), species(""), pdg(0), energy(0), vop(0), vbr(0), lgSet(0), hgSet(0), posX(0), posY(0), shapetime(0), assemblyNr(0), year(-1), month(-1), readout(""), facility(""), beamline(""), samples(0), trigDelay(0) {}
0075     int runNr;
0076     TString species;
0077     int pdg;
0078     float energy;
0079     float vop;
0080     float vbr;
0081     int lgSet;
0082     int hgSet;
0083     float posX;
0084     float posY;
0085     float shapetime;
0086     int assemblyNr;
0087     int year;
0088     int month;
0089     TString readout;
0090     TString facility;
0091     TString beamline;
0092     int samples;
0093     int trigDelay;
0094     int trigDead;
0095     int phase;
0096     int nFPGA;
0097     int nASIC;
0098   } ;
0099 
0100   TString GetStringFromRunInfo(RunInfo, Int_t);
0101 
0102   //__________________________________________________________________________________________________________
0103   //__________________ Read run infos from text file _________________________________________________________
0104   //__________________________________________________________________________________________________________    
0105   // specialData: 0 - std. TB, 1 - SPE data ORNL
0106   inline std::map<int,RunInfo> readRunInfosFromFile(TString runListFileName, int debug, int specialData = 0 ){
0107     std::map<int,RunInfo> runs;
0108     //std::cout << "INFO: You have given the following run list file: " << runListFileName.Data() << std::endl;
0109     std::ifstream runListFile;
0110     runListFile.open(runListFileName,std::ios_base::in);
0111     if (!runListFile) {
0112       std::cout << "ERROR: file " << runListFileName.Data() << " not found!" << std::endl;
0113       return runs;
0114     }
0115 
0116     TString facility="";
0117     TString beamline="";
0118     TString readout="";
0119     int year = -1;
0120     int month = -1;
0121     for( TString tempLine; tempLine.ReadLine(runListFile, kTRUE); ) {
0122       // check if line should be considered
0123       if (tempLine.BeginsWith("%") || tempLine.BeginsWith("#")){
0124         continue;
0125       }
0126       if (debug > 1) std::cout << tempLine.Data() << std::endl;
0127       
0128       TObjArray *tempArr2  = tempLine.Tokenize(" ");
0129       if(tempArr2->GetEntries()>0){
0130         if (tempLine.BeginsWith("year")){
0131           year=((TString)((TObjString*)tempArr2->At(1))->GetString()).Atoi();
0132           continue;
0133         } else if (tempLine.BeginsWith("month")){
0134           month=((TString)((TObjString*)tempArr2->At(1))->GetString()).Atoi();
0135           continue;
0136         } else if (tempLine.BeginsWith("readout")){
0137           readout=((TString)((TObjString*)tempArr2->At(1))->GetString());
0138           continue;
0139         } else if (tempLine.BeginsWith("facility")){
0140           facility=((TString)((TObjString*)tempArr2->At(1))->GetString());
0141           continue;
0142         } else if (tempLine.BeginsWith("beam-line")){ 
0143           beamline=((TString)((TObjString*)tempArr2->At(1))->GetString());
0144           continue;
0145         }
0146       }
0147       // Separate the string according to tabulators
0148       TObjArray *tempArr  = tempLine.Tokenize(",");
0149       if(tempArr->GetEntries()<1){
0150         if (debug > 1) std::cout << "nothing to be done" << std::endl;
0151         delete tempArr;
0152         continue;
0153       } 
0154 
0155       // Put them to the correct variables    
0156       RunInfo tempRun;
0157       tempRun.facility= facility;
0158       tempRun.beamline= beamline;
0159       tempRun.readout = readout;
0160       tempRun.year    = year; 
0161       tempRun.month    = month; 
0162       tempRun.runNr    = ((TString)((TObjString*)tempArr->At(0))->GetString()).Atoi();
0163       tempRun.species  =  (TString)((TObjString*)tempArr->At(1))->GetString();
0164       tempRun.pdg      = ((TString)((TObjString*)tempArr->At(2))->GetString()).Atoi();
0165       tempRun.energy   = ((TString)((TObjString*)tempArr->At(3))->GetString()).Atof();
0166       tempRun.vop      = ((TString)((TObjString*)tempArr->At(4))->GetString()).Atof();
0167       tempRun.vbr      = ((TString)((TObjString*)tempArr->At(5))->GetString()).Atof();
0168       
0169       if (readout.CompareTo("CAEN") == 0){
0170         tempRun.hgSet    = ((TString)((TObjString*)tempArr->At(6))->GetString()).Atoi();
0171         tempRun.lgSet    = ((TString)((TObjString*)tempArr->At(7))->GetString()).Atoi();
0172         if (tempArr->GetEntries() > 10){
0173          tempRun.shapetime = ((TString)((TObjString*)tempArr->At(10))->GetString()).Atof();
0174         }
0175       } else {
0176         tempRun.trigDelay = ((TString)((TObjString*)tempArr->At(6))->GetString()).Atoi();
0177         tempRun.samples   = ((TString)((TObjString*)tempArr->At(7))->GetString()).Atoi();
0178         tempRun.trigDead  = ((TString)((TObjString*)tempArr->At(10))->GetString()).Atoi();
0179         tempRun.phase     = ((TString)((TObjString*)tempArr->At(11))->GetString()).Atoi();
0180         tempRun.nFPGA     = ((TString)((TObjString*)tempArr->At(12))->GetString()).Atoi();
0181         tempRun.nASIC     = ((TString)((TObjString*)tempArr->At(13))->GetString()).Atoi();
0182       }
0183       tempRun.posX    = ((TString)((TObjString*)tempArr->At(8))->GetString()).Atoi();
0184       tempRun.posY    = ((TString)((TObjString*)tempArr->At(9))->GetString()).Atoi();
0185       if (specialData == 1) tempRun.assemblyNr = ((TString)((TObjString*)tempArr->At(10))->GetString()).Atoi();
0186                   
0187       if (debug > 1) std::cout << "Run " << tempRun.runNr << "\t species: " << tempRun.species << "\t energy: "  << tempRun.energy << "\t Vop: " << tempRun.vop << "\t Vov: " << tempRun.vop-tempRun.vbr << "\t Xbeam: " << tempRun.posX<< "\t Ybeam: " << tempRun.posY << "\t shaping time: " << tempRun.shapetime << std::endl;
0188       runs[tempRun.runNr]=tempRun;
0189     }
0190     std::cout << year << "-" << month << "\t:\t" << facility.Data() << "-" << beamline.Data() << "\t Readout: " << readout.Data() << std::endl;
0191     std::cout << "registered " << runs.size() << " runs from  "<< runListFileName.Data() << std::endl;
0192     
0193     return runs;
0194   };
0195 
0196   inline int GetSpeciesIntFromRunInfo(RunInfo currRunInfo){
0197       if (currRunInfo.species.Contains("cosmics")){
0198           return  0; // cosmics
0199       } else if (currRunInfo.species.CompareTo("g") == 0){
0200           return  1; // gamma
0201       } else if (currRunInfo.species.Contains("muon") || currRunInfo.species.Contains("Muon") || currRunInfo.species.CompareTo("mu-") == 0){
0202           return  2; // muon
0203       } else if (currRunInfo.species.Contains("Electron") || currRunInfo.species.Contains("electron") || currRunInfo.species.CompareTo("e-") == 0 ){
0204           return  3; // electron
0205       } else if (currRunInfo.species.Contains("Positron") || currRunInfo.species.Contains("positron") || currRunInfo.species.CompareTo("e+") == 0 ){
0206           return  6; // positron
0207       } else if (currRunInfo.species.Contains("Pion") || currRunInfo.species.Contains("pion") || currRunInfo.species.CompareTo("pi-") == 0 || currRunInfo.species.CompareTo("pi+") == 0 ){
0208           return  4; // pion
0209       } else if (currRunInfo.species.Contains("Hadron") || currRunInfo.species.Contains("hadron") || currRunInfo.species.CompareTo("h+") == 0 || currRunInfo.species.CompareTo("h-") == 0 ){
0210           return  5; // hadron/proton
0211       }
0212       
0213       return -1;
0214   }
0215   
0216   inline Double_t ReturnMipPlotRangeDepVov(double Vov, bool isHG, ReadOut::Type type){
0217     if (type == ReadOut::Type::Caen){
0218       if (isHG){
0219         if (Vov < 2)
0220           return 550.;
0221         else if (Vov < 3)
0222           return 750.;
0223         else if (Vov < 4)
0224           return 950.;
0225         else if (Vov < 5)
0226           return 1150.;
0227         else
0228           return 1350.;
0229       } else {
0230         if (Vov < 2)
0231           return 85.;
0232         else if (Vov < 3)
0233           return 105.;
0234         else if (Vov < 4)
0235           return 125.;
0236         else if (Vov < 5)
0237           return 145.;
0238         else
0239           return 165.;      
0240       }
0241     } else {
0242       return 410.;
0243     }
0244   }
0245   
0246 
0247   
0248   
0249 #endif