epic-lfhcal-tbana/NewStructure/CompareHGCROCCalib.cc

0001 #include "CompareHGCROCCalib.h"
0002 #include <vector>
0003 #include <algorithm>
0004 #include <sstream>
0005 #include <string>
0006 #include <cstdlib>
0007 #include "TROOT.h"
0008 #include "TPaletteAxis.h"
0009 #include "TProfile.h"
0010 #include "TChain.h"
0011 #include "TileSpectra.h"
0012 #include "PlotHelper.h"
0013 #include "CommonHelperFunctions.h"
0014
0015 #include "waveform_fitting/waveform_fit_base.h"
0016 #include "waveform_fitting/crystal_ball_fit.h"
0017 #include "waveform_fitting/max_sample_fit.h"
0018
0019
0020 // ****************************************************************************
0021 // Checking and opening input files
0022 // ****************************************************************************
0023 bool CompareHGCROCCalib::CheckAndOpenIO(void){
0024     int matchingbranch;
0025
0026     // reading files from a text file
0027     if( !inputFileList.IsNull() ){
0028         // text file with one .root file per line
0029         std::cout << "******************************" << std::endl;
0030         std::cout << "Reading from the file list: " << inputFileList.Data() << std::endl;
0031         std::cout << "******************************" << std::endl;
0032
0033         std::fstream dummyTxt;
0034         dummyTxt.open( inputFileList.Data(), std::ios::in );
0035         if( !dummyTxt.is_open() ){
0036             std::cout << "Error opening " << inputFileList.Data() << ", does the file exist?" << std::endl;
0037             return false;
0038         }
0039
0040         //initialize root file name
0041         std::string dummyRootName;
0042         //set first root file name
0043         dummyTxt>>dummyRootName;
0044
0045         int goodsetup; int gooddata;
0046         while( dummyTxt.good() ){
0047             std::cout << dummyRootName.data() << std::endl;
0048             //check that file exists and can be opened
0049             TFile dummyRoot = TFile( dummyRootName.data() , "READ");
0050             if( dummyRoot.IsZombie() ){
0051                 std::cout << "Error opening " << dummyRootName.data() << ", does the file exists?" << std::endl;
0052                 return false;
0053             }
0054             dummyRoot.Close();
0055             //add file name to setup and data chain
0056             goodsetup   = TsetupIn->Add( dummyRootName.c_str() );
0057             gooddata    = TdataIn->Add( dummyRootName.c_str() );
0058             if( goodsetup == 0 ){
0059                 std::cout << "Issues retrieving Setup tree from " << dummyRootName.data() <<", file is ignored" << std::endl;
0060             }
0061             if( gooddata == 0 ){
0062                 std::cout << "Issues retrieving Data tree from " << dummyRootName.data() <<", file is ignored" << std::endl;
0063             }
0064             //set next root file name
0065             dummyTxt>>dummyRootName;
0066         }
0067
0068         // Setup TChain of setups and data
0069         // initialize global variable setup
0070         setup = Setup::GetInstance();
0071         std::cout<<"Setup add " << setup << std::endl;
0072         matchingbranch=TsetupIn->SetBranchAddress("setup",&rswptr);
0073         if( matchingbranch < 0){
0074             std::cout << "Error retrieiving Setup info from the tree" << std::endl;
0075             return false;
0076         }
0077         std::cout << "Entries: " << TsetupIn->GetEntries() << std::endl;
0078         TsetupIn->GetEntry(0);
0079         setup->Initialize(*rswptr);
0080
0081         // initialize data with correct branch
0082         matchingbranch=TdataIn->SetBranchAddress("event",&eventptr);
0083         if( matchingbranch < 0){
0084             std::cout << "Error retrieiving data from the tree" << std::endl;
0085             return false;
0086         }
0087
0088     }
0089
0090     std::cout <<"=============================================================" << std::endl;
0091     std::cout <<" Basic setup complete" << std::endl;
0092     std::cout <<"=============================================================" << std::endl;
0093     return true;
0094 }
0095
0096
0097 // ****************************************************************************
0098 // Primary process function
0099 // ****************************************************************************
0100 bool CompareHGCROCCalib::Process(void){
0101     bool status;
0102
0103     status = CompareTOT();
0104
0105     return status;
0106 }
0107
0108 bool CompareHGCROCCalib::CompareTOT(void){
0109
0110     int nentries    = TsetupIn->GetEntries();
0111
0112     std::vector<int>    *values     = new std::vector<int>[5]; // 0 - tile nr (from 0 to n_channels-1), 1 - inj val (run number), 2 - TOT, 3 - ADC, 4 - cellID
0113     int     n_channels = 0; // will be read in the loop
0114     int pointCounter = 0;
0115
0116     // for plotting
0117     int max_inj = 0;
0118     int max_adc = 0;
0119     int max_tot = 0;
0120
0121 // ******************************************************************************************
0122 // Iterate over all entries (runs) in the calib tree
0123 // ******************************************************************************************
0124     for(int ientry = 0; ientry<nentries; ientry++){
0125         TsetupIn->GetEntry(ientry);
0126         TdataIn->GetEntry(ientry);
0127
0128         if( ientry == 0 ){
0129             n_channels = event.GetNTiles();
0130         }
0131
0132
0133         for(int i = 0; i<event.GetNTiles(); i++){
0134             Hgcroc*     aTile = (Hgcroc*)event.GetTile(i);
0135             int tmpTileID   = i;
0136             int tmpInjVal   = event.GetRunNumber();
0137             int tmpTOT      = aTile->GetTOT();
0138             int tmpADC      = aTile->GetIntegratedADC();
0139             int tmpCellID   = aTile->GetCellID();
0140
0141             if( tmpInjVal > max_inj )   max_inj = tmpInjVal;
0142             if( tmpTOT  > max_tot )     max_tot = tmpTOT;
0143             if( tmpADC  > max_adc )     max_adc = tmpADC;
0144
0145             if(debug > 1) std::cout << "Inj val: " << tmpInjVal << "\t cellID: " << aTile->GetCellID() << "\t i-th tile: " << i << "\t TOT: " << tmpTOT << "\t ADC: " << tmpADC << std::endl;
0146
0147             values[0].push_back( tmpTileID );
0148             values[1].push_back( tmpInjVal );
0149             values[2].push_back( tmpTOT );
0150             values[3].push_back( tmpADC );
0151             values[4].push_back( tmpCellID );
0152             pointCounter++;
0153         }
0154
0155         if( debug > 1 ) std::cout << "Point counter: " << pointCounter << "\t n_points in one graph: " << pointCounter / n_channels<< std::endl;
0156     }
0157
0158     // we have all the values "flattened", so we will need to detangle it into n_channels TGraphs - separated based on cellID
0159     TGraph* graphTOT[n_channels];
0160     TGraph* graphADC[n_channels];
0161     int    cellID[n_channels];
0162
0163     for(int i=0; i<n_channels; i++){
0164         graphTOT[i] = new TGraph( pointCounter / n_channels );
0165         graphADC[i] = new TGraph( pointCounter / n_channels );
0166         for(int j=0; j<pointCounter; j++){
0167             if( values[0].at(j) == i ){
0168                 graphTOT[i]->SetPoint(j / n_channels, values[1].at(j), values[2].at(j) );
0169                 graphADC[i]->SetPoint(j / n_channels, values[1].at(j), values[3].at(j) );
0170                 cellID[i]   = values[4].at(j);
0171             }
0172         }
0173     }
0174
0175     Int_t textSizePixel = 30;
0176
0177     for(int i=0; i<n_channels; i++){
0178         TCanvas*    canvas  = new TCanvas("canvas", "canvas", 800, 600);
0179         DefaultCanvasSettings( canvas, 0.13, 0.02, 0.045, 0.11);
0180
0181         SetMarkerDefaultsTGraph( graphTOT[i], 24, 0.8, kRed+1, kRed+1);
0182         SetMarkerDefaultsTGraph( graphADC[i], 24, 0.8, kBlue+2, kBlue+2);
0183
0184         TH2F*   hDummy  = new TH2F("hDummy", "hDummy", 100, 0, max_inj+50, 100, 0, (max_tot > max_adc ) ? max_tot * 1.2 : max_adc * 1.2 );
0185         SetStyleHistoTH2ForGraphs(hDummy, "inj val", "TOT/ADC", 0.85*textSizePixel, textSizePixel, 0.85*textSizePixel, textSizePixel,0.9, 1.5, 510, 510, 43, 63);
0186
0187         hDummy->Draw();
0188         graphTOT[i]->Draw("same,p");
0189         graphADC[i]->Draw("same,p");
0190
0191         TLatex  latexCellID     = TLatex(0.94, 0.89, Form("#bf{Cell %d}", cellID[i]));
0192         latexCellID.SetNDC();
0193         latexCellID.SetTextFont(42);
0194         latexCellID.SetTextSize(0.045);
0195         latexCellID.SetTextAlign(32);
0196         latexCellID.Draw("same");
0197
0198         TLegend* legend = GetAndSetLegend(0.15, 0.82, 2);
0199         legend->AddEntry( graphTOT[i], "TOT", "p");
0200         legend->AddEntry( graphADC[i], "ADC", "p");
0201         legend->Draw("same");
0202
0203
0204         canvas->SaveAs( Form("%s/TOT_ADCvsInjVal_Cell%d.pdf", outputDirPlots.Data(), cellID[i] ) );
0205         delete canvas; delete hDummy;
0206     }
0207
0208     return true;
0209 }