Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:15:46

0001 #include <TROOT.h>
0002 #include <TString.h>
0003 #include <TObject.h>
0004 #include <TObjString.h>
0005 #include <TSystem.h>
0006 #include <TChain.h>
0007 #include <TMath.h>
0008 #include <TVector3.h>
0009 #include <iostream>
0010 #include <fstream>
0011 #include <TParticlePDG.h>
0012 #include <TDatabasePDG.h>
0013 #include <TRandom3.h>
0014 
0015 #include <TCanvas.h>
0016 #include <TPad.h>
0017 #include <TH1.h>
0018 #include <TH1D.h>
0019 #include <TH1F.h>
0020 #include <TH2.h>
0021 #include <TH3.h>
0022 #include <TFile.h>
0023 #include <TH2D.h>
0024 #include <TH2F.h>
0025 #include <TString.h>
0026 #include <TDatime.h>
0027 #include <TF1.h>
0028 #include <TF2.h>
0029 #include <THStack.h>
0030 #include <TGraph.h>
0031 #include <TStyle.h>
0032 #include <TGraphAsymmErrors.h>
0033 #include <TLine.h>
0034 #include <TLatex.h>
0035 #include <TArrow.h>
0036 #include <TGraphErrors.h>
0037 #include <TGaxis.h>
0038 #include <TLegend.h>
0039 #include <TFrame.h>
0040 #include <TLorentzVector.h>
0041 
0042 #include "PlottingHeader.h"
0043 #include "FittingHeader.h"
0044 #include "CommonVariables.h"
0045 
0046 const int gMaxLayers    = 16;
0047 const int gMaxSpecies   = 2;
0048 const int gMaxVoltages  = 10;
0049 const int gMaxBEn       = 11;
0050 const int gMaxGainS     = 7;
0051 
0052 struct mipDataPoint{
0053     mipDataPoint(): layer(0), channel(0), mpvL(0), empvL(0), max(0), emax(0), vop(0), runnr(0), beame(0), isHad (true), lgSet(0), hgSet(0) {}
0054     int layer;
0055     int channel;
0056     float mpvL;
0057     float empvL;
0058     float max;
0059     float emax;
0060     float vop;
0061     int runnr;
0062     float beame;
0063     bool isHad;
0064     int lgSet;
0065     int hgSet;
0066 } ;
0067 
0068 
0069 //__________________________________________________________________________________________________________
0070 //_____________________MAIN function !!! ___________________________________________________________________
0071 //__________________________________________________________________________________________________________
0072 void CompareDifferentRuns( TString configFileName     = "", 
0073                            TString outputDir          = "Compare/",
0074                             Int_t verbosity           = 0,
0075                             Int_t tbdata              = 0,        // 1: SPS 2023, 2: PS 2023
0076                             TString mappingFile       = "",
0077                             TString runListFileName   = "configs/SPS_RunNumbers.txt"
0078                             
0079                                   ){
0080     StyleSettingsThesis("pdf");
0081     SetPlotStyle();
0082     Double_t textSizeRel = 0.04;
0083     
0084       // make output directory
0085     TString dateForOutput = ReturnDateStr();
0086     TString outputDirPlots = Form("%s/%s",outputDir.Data(), dateForOutput.Data());
0087     TString outputDirPlotsH = Form("%s/%s/Hadron",outputDir.Data(), dateForOutput.Data());
0088     TString outputDirPlotsE = Form("%s/%s/Electron",outputDir.Data(), dateForOutput.Data());
0089     gSystem->Exec("mkdir -p "+outputDir);
0090     gSystem->Exec("mkdir -p "+outputDirPlots);
0091     gSystem->Exec("mkdir -p "+outputDirPlotsH);
0092     std::cout << outputDirPlotsH.Data() << std::endl;
0093     gSystem->Exec("mkdir -p "+outputDirPlotsE);
0094     TString labelTB = "";
0095     Double_t rangeMaxMIP[2] = {0,1000};
0096     if (tbdata == 1){
0097       labelTB = "SPS, LFHCal-Prototype 1";
0098       rangeMaxMIP[1]  = 1200;
0099     }
0100     
0101     // ********************************************************************************************************
0102     // read run list and corresponding settings
0103     // ********************************************************************************************************
0104     std::vector<runInfo> allRuns = readRunInfosFromFile(runListFileName, 1 );
0105 
0106     // ********************************************************************************************************    
0107     // read folder and name from file
0108     // ********************************************************************************************************
0109     std::vector<Int_t> runnumbers;
0110     std::vector<TString> fileNames;
0111     ifstream in;
0112     in.open(configFileName,ios_base::in);
0113     if (!in) {
0114         std::cout << "ERROR: file " << configFileName.Data() << " not found!" << std::endl;
0115         return;
0116     }
0117 
0118     for( TString tempLine; tempLine.ReadLine(in, kTRUE); ) {
0119       // check if line should be considered
0120       if (tempLine.BeginsWith("%") || tempLine.BeginsWith("#")){
0121           continue;
0122       }
0123       if (verbosity > 0) std::cout << tempLine.Data() << std::endl;
0124 
0125       // Separate the string according to tabulators
0126       TObjArray *tempArr  = tempLine.Tokenize("\t");
0127       if(tempArr->GetEntries()<1){
0128           if (verbosity > 1) std::cout << "nothing to be done" << std::endl;
0129           delete tempArr;
0130           continue;
0131       } else if (tempArr->GetEntries() == 1 ){
0132           // Separate the string according to space
0133           tempArr       = tempLine.Tokenize(" ");
0134           if(tempArr->GetEntries()<1){
0135               if (verbosity > 1) std::cout << "nothing to be done" << std::endl;
0136               delete tempArr;
0137               continue;
0138           } else if (tempArr->GetEntries() == 1  ) {
0139               if (verbosity > 1) std::cout << ((TString)((TObjString*)tempArr->At(0))->GetString()).Data() << " has not been reset, no value given!" << std::endl;
0140               delete tempArr;
0141               continue;
0142           }
0143       }
0144 
0145       TString tempname  = ((TString)((TObjString*)tempArr->At(1))->GetString());
0146       Int_t temprun     = ((TString)((TObjString*)tempArr->At(0))->GetString()).Atoi();
0147       runnumbers.push_back(temprun);
0148       fileNames.push_back(tempname);
0149     } 
0150   
0151     std::cout<<"=============================================================="<<std::endl;
0152     std::vector<runInfo> compRuns;
0153     std::vector<TGraphErrors*>graphsMPV;
0154     std::vector<TGraphErrors*>graphsMax;
0155     
0156     Int_t runsToComp = (Int_t)runnumbers.size();
0157     
0158     for (Int_t r = 0; r < runsToComp; r++){
0159       
0160       Int_t indexCRun = findCurrentRun(allRuns, runnumbers.at(r));
0161       runInfo currentRunInfo;
0162       if (indexCRun < 0){
0163         std::cout << "run not in current list of runs, provided" << std::endl;
0164         return;
0165       } else {
0166         std::cout << "Run "<< runnumbers.at(r) << "\t found at index " << indexCRun << std::endl;
0167         currentRunInfo = GetRunInfoObject(allRuns,indexCRun);
0168       }
0169       compRuns.push_back(currentRunInfo);
0170       std::cout <<  runnumbers.at(r) << "\t" << fileNames.at(r).Data() << "\t" << compRuns.at(r).vop << std::endl;
0171       
0172       TFile* tempFile = new TFile(fileNames.at(r).Data(), "OPEN");
0173       if (tempFile->IsZombie()){
0174           std::cout << fileNames.at(r).Data() << " is broken, please remove from list or fix!" << std::endl;
0175       }
0176       TGraphErrors* tempGraph = (TGraphErrors*)tempFile->Get("graphMPV_HG_channels");
0177       graphsMPV.push_back(tempGraph);
0178       tempGraph = (TGraphErrors*)tempFile->Get("graphMax_HG_channels");
0179       graphsMax.push_back(tempGraph);
0180     }
0181     std::cout<<"=============================================================="<<std::endl;
0182 
0183     std::vector<mipDataPoint> mipPoints;
0184     for (Int_t r = 0; r < runsToComp; r++){
0185       std::cout <<  runnumbers.at(r) << "\t" << fileNames.at(r).Data() << "\t" << compRuns.at(r).species.Data() << "\t" << compRuns.at(r).vop << std::endl;
0186       if (graphsMPV.at(r)->GetN() !=graphsMax.at(r)->GetN() ){
0187         std::cout << "MPV and Max graph don't have same size ... something is wrong!" << std::endl;
0188         continue;
0189       }
0190       for (Int_t e = 0; e < (Int_t)graphsMPV.at(r)->GetN(); e++){
0191         mipDataPoint currMip;
0192         currMip.layer   = (int)(graphsMPV.at(r)->GetX()[e]/10);
0193         currMip.channel = (int)(graphsMPV.at(r)->GetX()[e]-10*currMip.layer);
0194         currMip.mpvL    = (float)(graphsMPV.at(r)->GetY()[e]);
0195         currMip.empvL   = (float)(graphsMPV.at(r)->GetEY()[e]);
0196         currMip.max     = (float)(graphsMax.at(r)->GetY()[e]);
0197         currMip.emax    = (float)(graphsMax.at(r)->GetEY()[e]);
0198         currMip.vop     = (float)(compRuns.at(r).vop);
0199         currMip.runnr   = (int)(compRuns.at(r).runNr);
0200         currMip.beame   = (float)(compRuns.at(r).energy);
0201         if (compRuns.at(r).species.CompareTo("h") == 0)
0202           currMip.isHad   = true;
0203         else 
0204           currMip.isHad   = false;
0205         currMip.lgSet   = (float)(compRuns.at(r).lgSet);
0206         currMip.hgSet   = (float)(compRuns.at(r).hgSet);
0207         
0208 //         std::cout << graphsMPV.at(r)->GetX()[e] << "\t l: " << currMip.layer << "\t c: " << currMip.channel << "\t vop: " << currMip.vop << std::endl;
0209         mipPoints.push_back(currMip);
0210       }      
0211     }
0212     
0213     Float_t tempYVal[5000];
0214     Float_t tempXVal[5000];
0215     Float_t tempEYVal[5000];
0216     Float_t tempEXVal[5000];
0217  
0218     TGraphErrors* graphVoVvsMipMax[gMaxSpecies][gMaxBEn][gMaxGainS];
0219     TGraphErrors* graphRunNrvsMipMax[gMaxSpecies][gMaxBEn][gMaxGainS][gMaxVoltages];
0220     for (Int_t s = 0; s < gMaxSpecies; s++){
0221       for (Int_t gs = 0; gs < gMaxGainS; gs++){
0222         for (Int_t be = 0; be < gMaxBEn; be++){
0223           graphVoVvsMipMax[s][be][gs]   = new TGraphErrors();
0224           for (Int_t v = 0; v < gMaxVoltages; v++){
0225             graphRunNrvsMipMax[s][be][gs][v] = new TGraphErrors();
0226           }
0227         }
0228       }
0229     }
0230     for (Int_t p = 0; p< (Int_t)mipPoints.size(); p++){
0231       Int_t s   = mipPoints.at(p).isHad; // will need to be adjusted if we don't have h & e
0232       Int_t be = GetBeamEnergyIndex(mipPoints.at(p).beame, tbdata);
0233       Int_t gs = GetGainSetIndex(mipPoints.at(p).hgSet, mipPoints.at(p).lgSet, tbdata);
0234       Int_t v  = GetOverVoltageIndex(mipPoints.at(p).vop, tbdata);
0235       
0236       if (be == -1 || gs == -1){
0237         std::cout << "beam E: "<< mipPoints.at(p).beame << "\t or HG set: " << mipPoints.at(p).hgSet << "\t or LG set: " << mipPoints.at(p).hgSet << std::endl;
0238         continue;
0239       }
0240       if (v == -1){
0241         std::cout << "Vov: "<< mipPoints.at(p).vop << std::endl;
0242         continue;
0243       }
0244       graphVoVvsMipMax[s][be][gs]->SetPoint(graphVoVvsMipMax[s][be][gs]->GetN(),mipPoints.at(p).vop, mipPoints.at(p).max );
0245       graphVoVvsMipMax[s][be][gs]->SetPointError(graphVoVvsMipMax[s][be][gs]->GetN()-1,0, mipPoints.at(p).emax );
0246       graphRunNrvsMipMax[s][be][gs][v]->SetPoint(graphRunNrvsMipMax[s][be][gs][v]->GetN(),mipPoints.at(p).runnr, mipPoints.at(p).max );
0247       graphRunNrvsMipMax[s][be][gs][v]->SetPointError(graphRunNrvsMipMax[s][be][gs][v]->GetN()-1,0, mipPoints.at(p).emax );
0248     }
0249 
0250     for (Int_t s = 0; s < gMaxSpecies; s++){
0251       for (Int_t gs = 0; gs < gMaxGainS; gs++){
0252         for (Int_t be = 0; be < gMaxBEn; be++){
0253           if (graphVoVvsMipMax[s][be][gs]->GetN() != 0){
0254             graphVoVvsMipMax[s][be][gs]->Sort();
0255           }
0256           for (Int_t v = 0; v < gMaxVoltages; v++){
0257             if (graphRunNrvsMipMax[s][be][gs][v]->GetN() != 0){
0258               graphRunNrvsMipMax[s][be][gs][v]->Sort();
0259             }
0260           }
0261         }
0262       }
0263     }
0264 
0265     
0266     
0267                                             // beam type, beam energy, gain setting,  layer,   channel
0268     TGraphErrors* graphPerChannelVoVvsMipMax[gMaxSpecies] [gMaxBEn]  [gMaxGainS]   [gMaxLayers][9];
0269                                             // beam type, voltage, gain setting,  layer,   channel
0270     TGraphErrors* graphPerChannelBeamEvsMipMax[gMaxSpecies] [gMaxVoltages]  [gMaxGainS]   [gMaxLayers][9];
0271     bool layerExists[gMaxLayers];
0272     for (Int_t l = 0; l < gMaxLayers; l++){
0273       layerExists[l] = false;
0274       for (Int_t cb = 0; cb < 9; cb++){
0275         for (Int_t s = 0; s < gMaxSpecies; s++){
0276           for (Int_t gs = 0; gs < gMaxGainS; gs++){
0277             for (Int_t be = 0; be < gMaxBEn; be++){
0278               graphPerChannelVoVvsMipMax[s][be][gs][l][cb] = new TGraphErrors();
0279             }
0280             for (Int_t v = 0; v < gMaxVoltages; v++){
0281               graphPerChannelBeamEvsMipMax[s][v][gs][l][cb] = new TGraphErrors();
0282             }
0283           }
0284         }
0285       }
0286     }
0287     
0288     // filter for different conditions
0289     for (Int_t p = 0; p< (Int_t)mipPoints.size(); p++){
0290       Int_t l = mipPoints.at(p).layer;
0291       Int_t c = mipPoints.at(p).channel;
0292       Int_t s = mipPoints.at(p).isHad; // will need to be adjusted if we don't have h & e
0293       Int_t be = GetBeamEnergyIndex(mipPoints.at(p).beame, tbdata);
0294       Int_t gs = GetGainSetIndex(mipPoints.at(p).hgSet, mipPoints.at(p).lgSet, tbdata);
0295       Int_t v  = GetOverVoltageIndex(mipPoints.at(p).vop, tbdata);
0296       if (be == -1 || gs == -1){
0297         std::cout << "beam E: "<< mipPoints.at(p).beame << "\t or HG set: " << mipPoints.at(p).hgSet << "\t or LG set: " << mipPoints.at(p).hgSet << std::endl;
0298         continue;
0299       }
0300       if (v == -1){
0301         std::cout << "Vov: "<< mipPoints.at(p).vop << std::endl;
0302         continue;
0303       }
0304       if (!layerExists[l]) layerExists[l] = true;        
0305       // add point to channel graph vs voltage
0306       graphPerChannelVoVvsMipMax[s][be][gs][l][c]->SetPoint(graphPerChannelVoVvsMipMax[s][be][gs][l][c]->GetN(),mipPoints.at(p).vop, mipPoints.at(p).max );
0307       graphPerChannelVoVvsMipMax[s][be][gs][l][c]->SetPointError(graphPerChannelVoVvsMipMax[s][be][gs][l][c]->GetN()-1,0, mipPoints.at(p).emax );
0308       // add point to channel graph vs beamE
0309       graphPerChannelBeamEvsMipMax[s][v][gs][l][c]->SetPoint(graphPerChannelBeamEvsMipMax[s][v][gs][l][c]->GetN(),mipPoints.at(p).beame, mipPoints.at(p).max );
0310       graphPerChannelBeamEvsMipMax[s][v][gs][l][c]->SetPointError(graphPerChannelBeamEvsMipMax[s][v][gs][l][c]->GetN()-1,0, mipPoints.at(p).emax );
0311     }
0312         
0313     for (Int_t l = 0; l < gMaxLayers; l++){
0314       for (Int_t cb = 0; cb < 9; cb++){
0315         for (Int_t s = 0; s < gMaxSpecies; s++){
0316           for (Int_t gs = 0; gs < gMaxGainS; gs++){        
0317             for (Int_t be = 0; be < gMaxBEn; be++){
0318               if (graphPerChannelVoVvsMipMax[s][be][gs][l][cb]->GetN() != 0){
0319                 graphPerChannelVoVvsMipMax[s][be][gs][l][cb]->Sort();
0320               }
0321             }
0322             for (Int_t v = 0; v < gMaxVoltages; v++){
0323               if (graphPerChannelBeamEvsMipMax[s][v][gs][l][cb]->GetN() != 0){
0324                 graphPerChannelBeamEvsMipMax[s][v][gs][l][cb]->Sort();
0325               }
0326             }
0327           }
0328         }
0329       }
0330     }
0331     
0332     
0333     TCanvas* canvas1DSimple = new TCanvas("canvas1DSimple","",0,0,700,500);  // gives the page size
0334     DefaultCancasSettings( canvas1DSimple, 0.08, 0.02, 0.02, 0.08);
0335     
0336     TH1D* histDummyVvsMPV = new TH1D("histDummyVvsMPV","", 70, 0, 7 );
0337     SetStyleHistoTH1ForGraphs( histDummyVvsMPV, "#it{V}_{ov} (V)", "MPV_{lan-gaus} (HG ADC)",0.85* textSizeRel, textSizeRel, 0.85* textSizeRel, textSizeRel, 0.9, 0.95);
0338     histDummyVvsMPV->GetYaxis()->SetRangeUser(rangeMaxMIP[0], rangeMaxMIP[1]);
0339     histDummyVvsMPV->Draw();
0340     Bool_t lenteredAll[2]  = {0, 0};
0341     Bool_t plottedA      = 0; 
0342     TLegend* legend = GetAndSetLegend2( 0.12, 0.95-2*textSizeRel, 0.25, 0.95,textSizeRel, 2, "Beam:", 42,0.4);
0343     for (Int_t gs = 0; gs < gMaxGainS; gs++){
0344       for (Int_t be = 0; be < gMaxBEn; be++){
0345         if (graphVoVvsMipMax[1][be][gs]->GetN() > 0){
0346           SetMarkerDefaultsTGraphErr( graphVoVvsMipMax[1][be][gs], 24, 0.8, kRed+1,kRed+1);
0347           graphVoVvsMipMax[1][be][gs]->Draw("p,e,same");
0348           if (lenteredAll[1] == 0){
0349             legend->AddEntry(graphVoVvsMipMax[1][be][gs], "h", "p");
0350             lenteredAll[1]  = 1;
0351             plottedA       = 1;
0352           }
0353         }
0354         if (graphVoVvsMipMax[0][be][gs]->GetN() > 0){
0355           SetMarkerDefaultsTGraphErr( graphVoVvsMipMax[0][be][gs], 25, 0.8, kBlue+1,kBlue+1);
0356           graphVoVvsMipMax[0][be][gs]->Draw("p,e,same");
0357           if (lenteredAll[0] == 0){
0358             legend->AddEntry(graphVoVvsMipMax[0][be][gs], "e", "p");
0359             lenteredAll[0]  = 1;
0360             plottedA       = 1;
0361           }
0362         }
0363       }
0364     }
0365     DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0366     legend->Draw();
0367     if (plottedA) canvas1DSimple->SaveAs(Form("%s/VopVsMPVlangaus.pdf",outputDirPlots.Data()));
0368     
0369     TH1D* histDummyRunNrvsMPV = new TH1D("histDummyRunNrvsMPV","", 201, -0.5,200.5  );
0370     SetStyleHistoTH1ForGraphs( histDummyRunNrvsMPV, "Run Nr.", "MPV_{lan-gaus} (HG ADC)",0.85* textSizeRel, textSizeRel, 0.85* textSizeRel, textSizeRel, 0.9, 0.95);
0371     histDummyRunNrvsMPV->GetYaxis()->SetRangeUser(rangeMaxMIP[0], rangeMaxMIP[1]);
0372     histDummyRunNrvsMPV->Draw();
0373     for (Int_t gs = 0; gs < gMaxGainS; gs++){
0374       for (Int_t be = 0; be < gMaxBEn; be++){
0375         for (Int_t v = 0; v < gMaxVoltages; v++){
0376           if (graphRunNrvsMipMax[1][be][gs][v]->GetN() > 0){
0377             SetMarkerDefaultsTGraphErr( graphRunNrvsMipMax[1][be][gs][v], 24, 0.8, kRed+1,kRed+1);
0378             graphRunNrvsMipMax[1][be][gs][v]->Draw("p,e,same");
0379           }
0380           if (graphRunNrvsMipMax[0][be][gs][v]->GetN() > 0){
0381             SetMarkerDefaultsTGraphErr( graphRunNrvsMipMax[0][be][gs][v], 25, 0.8, kBlue+1,kBlue+1);
0382             graphRunNrvsMipMax[0][be][gs][v]->Draw("p,e,same");
0383           }
0384         }
0385       }
0386     }
0387     DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0388     legend->Draw();
0389     if (plottedA) canvas1DSimple->SaveAs(Form("%s/RunNrVsMPVlangaus.pdf",outputDirPlots.Data()));
0390 
0391     
0392     
0393     histDummyVvsMPV->Draw();
0394     Bool_t lenteredBE[gMaxBEn][gMaxSpecies];
0395     for (Int_t be = 0; be < gMaxBEn; be++ ){
0396       for (Int_t s = 0; s < gMaxSpecies; s++ ){
0397         lenteredBE[be][s] = 0;
0398       }
0399     }
0400     
0401     plottedA      = 0; 
0402     TLegend* legendE = GetAndSetLegend2( 0.12, 0.945-5*textSizeRel, 0.7, 0.945,textSizeRel, 4, "Beam species, E (GeV):", 42,0.15);
0403     for (Int_t gs = 0; gs < gMaxGainS; gs++){
0404       for (Int_t be = 0; be < gMaxBEn; be++){
0405         Float_t beamEnergy = GetBeamEnergyFromIndex(be, tbdata);
0406         for (Int_t s = 0; s < gMaxSpecies; s++){
0407             if (graphVoVvsMipMax[s][be][gs]->GetN() > 0){
0408               SetMarkerDefaultsTGraphErr( graphVoVvsMipMax[s][be][gs], (s==1) ? 24 : 25, 0.8, colorBeamE[be],colorBeamE[be]);
0409               graphVoVvsMipMax[s][be][gs]->Draw("p,e,same");
0410               if (lenteredBE[be][s] == 0){
0411                 legendE->AddEntry(graphVoVvsMipMax[s][be][gs], Form("%s, %d",(s==1) ? "h" : "e", (Int_t)beamEnergy), "p");
0412                 lenteredBE[be][s]  = 1;
0413                 plottedA          = 1;
0414               }
0415             }
0416         }
0417       }
0418     }
0419     DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0420     legendE->Draw();
0421     if (plottedA) canvas1DSimple->SaveAs(Form("%s/VopVsMPVlangaus_withBeamEnergy.pdf",outputDirPlots.Data()));
0422     
0423     histDummyRunNrvsMPV->Draw();
0424     for (Int_t gs = 0; gs < gMaxGainS; gs++){
0425       for (Int_t be = 0; be < gMaxBEn; be++){
0426         Float_t beamEnergy = GetBeamEnergyFromIndex(be, tbdata);
0427         for (Int_t s = 0; s < gMaxSpecies; s++){
0428           for (Int_t v = 0; v < gMaxVoltages; v++){
0429             if (graphRunNrvsMipMax[s][be][gs][v]->GetN() > 0){
0430               SetMarkerDefaultsTGraphErr( graphRunNrvsMipMax[s][be][gs][v], (s==1) ? 24 : 25, 0.8, colorBeamE[be],colorBeamE[be]);
0431               graphRunNrvsMipMax[s][be][gs][v]->Draw("p,e,same");
0432             }
0433           }
0434         }
0435       }
0436     }
0437     DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0438     legendE->Draw();
0439     if (plottedA) canvas1DSimple->SaveAs(Form("%s/RunNrVsMPVlangaus_withBeamEnergy.pdf",outputDirPlots.Data()));
0440 
0441 
0442     Bool_t lenteredV[gMaxBEn][gMaxSpecies];
0443     for (Int_t v = 0; v < gMaxBEn; v++ ){
0444       for (Int_t s = 0; s < gMaxSpecies; s++ ){
0445         lenteredV[v][s] = 0;
0446       }
0447     }
0448 
0449     histDummyRunNrvsMPV->Draw();
0450     plottedA      = 0; 
0451     TLegend* legendV = GetAndSetLegend2( 0.12, 0.945-4.5*textSizeRel, 0.7, 0.945,textSizeRel, 4, "Beam species, V_{ov} (V):", 42,0.15);
0452     for (Int_t v = 0; v < gMaxVoltages; v++){
0453       for (Int_t gs = 0; gs < gMaxGainS; gs++){
0454         Float_t voltage = GetOverVoltageFromIndex(v, tbdata);
0455         for (Int_t s = 0; s < gMaxSpecies; s++){
0456           for (Int_t be = 0; be < gMaxBEn; be++){
0457             if (graphRunNrvsMipMax[s][be][gs][v]->GetN() > 0){
0458               SetMarkerDefaultsTGraphErr( graphRunNrvsMipMax[s][be][gs][v], (s==1) ? 24 : 25, 0.8, colorVoltage[v],colorVoltage[v]);
0459               graphRunNrvsMipMax[s][be][gs][v]->Draw("p,e,same");
0460               if (lenteredV[v][s] == 0){
0461                 legendV->AddEntry(graphRunNrvsMipMax[s][be][gs][v], Form("%s, %1.1f",(s==1) ? "h" : "e", voltage), "p");
0462                 lenteredV[v][s]  = 1;
0463                 plottedA          = 1;
0464               }
0465             }
0466           }
0467         }
0468       }
0469     }
0470     DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0471     legendV->Draw();
0472     if (plottedA) canvas1DSimple->SaveAs(Form("%s/RunNrVsMPVlangaus_withVoltage.pdf",outputDirPlots.Data()));
0473     
0474 
0475     
0476     TH1D* histDummyBeamEvsMPV = new TH1D("histDummyBeamEvsMPV","", 361, -0.5, 360.5 );
0477     SetStyleHistoTH1ForGraphs( histDummyBeamEvsMPV, "#it{E}_{beam} (GeV)", "MPV_{lan-gaus} (HG ADC)",0.85* textSizeRel, textSizeRel, 0.85* textSizeRel, textSizeRel, 0.9, 0.95);
0478     histDummyBeamEvsMPV->GetYaxis()->SetRangeUser(rangeMaxMIP[0], rangeMaxMIP[1]);
0479     
0480     // ****************************************************************************************
0481     // Create electron beam plots
0482     // ****************************************************************************************    
0483     for (Int_t l = 0; l < gMaxLayers; l++ ){
0484       if (layerExists[l]){
0485         for (Int_t gs = 0; gs < gMaxGainS; gs++){
0486             for (Int_t be = 0; be < gMaxBEn; be++){
0487             canvas1DSimple->cd();
0488             canvas1DSimple->Clear();
0489             histDummyVvsMPV->Draw();
0490             Bool_t lentered[9]  = {0, 0, 0, 0, 0, 0, 0, 0, 0};
0491             Bool_t plotted      = 0; 
0492             TLegend* legend = GetAndSetLegend2( 0.12, 0.95-3*textSizeRel, 0.4, 0.95,textSizeRel, 4, Form("Layer %d, channel:",l), 42,0.2);
0493           
0494             for (Int_t cb = 1; cb < 9; cb++){
0495               if (graphPerChannelVoVvsMipMax[1][be][gs][l][cb]->GetN() > 0){
0496                 SetMarkerDefaultsTGraphErr( graphPerChannelVoVvsMipMax[1][be][gs][l][cb], markerReadBoard[cb-1], 0.8, colorReadBoard[cb-1],colorReadBoard[cb-1]);
0497                 graphPerChannelVoVvsMipMax[1][be][gs][l][cb]->Draw("p,e,same");
0498                 if (lentered[cb] == 0){
0499                   legend->AddEntry(graphPerChannelVoVvsMipMax[1][be][gs][l][cb], Form("%d",cb), "p");
0500                   lentered[cb]  = 1;
0501                   plotted       = 1;
0502                 }
0503               }
0504             }
0505             
0506             Int_t lowGainSet  = GetLGGainSetIndex(gs, tbdata);
0507             Int_t highGainSet = GetHGGainSetIndex(gs, tbdata);
0508             Float_t beamEnergy = GetBeamEnergyFromIndex(be, tbdata);
0509             DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0510             DrawLatex(0.95, 0.92-textSizeRel, Form("h-beam E =%.0f GeV", beamEnergy), true, textSizeRel, 42);        
0511             DrawLatex(0.95, 0.92-2*textSizeRel, Form("HG amp: %d, LG amp: %d", highGainSet, lowGainSet), true, textSizeRel, 42);
0512             legend->Draw();
0513             
0514             if (plotted) canvas1DSimple->SaveAs(Form("%s/VopVsMPVlangaus_Ebeam%03.0fGeV_Layer%02d_HG%02d_LG%02d.pdf",outputDirPlotsH.Data(), beamEnergy,  l, highGainSet, lowGainSet));
0515           }
0516         }
0517       }
0518     } 
0519     for (Int_t l = 0; l < gMaxLayers; l++ ){
0520       if (layerExists[l]){
0521         for (Int_t gs = 0; gs < gMaxGainS; gs++){
0522           for (Int_t v = 0; v < gMaxVoltages; v++){
0523             canvas1DSimple->cd();
0524             canvas1DSimple->Clear();
0525             histDummyBeamEvsMPV->Draw();
0526             Bool_t lentered[9]  = {0, 0, 0, 0, 0, 0, 0, 0, 0};
0527             Bool_t plotted      = 0; 
0528             TLegend* legend = GetAndSetLegend2( 0.12, 0.95-3*textSizeRel, 0.4, 0.95,textSizeRel, 4, Form("Layer %d, channel:",l), 42,0.2);
0529           
0530             for (Int_t cb = 1; cb < 9; cb++){
0531               if (graphPerChannelBeamEvsMipMax[1][v][gs][l][cb]->GetN() > 0){
0532                 SetMarkerDefaultsTGraphErr( graphPerChannelBeamEvsMipMax[1][v][gs][l][cb], markerReadBoard[cb-1], 0.8, colorReadBoard[cb-1],colorReadBoard[cb-1]);
0533                 graphPerChannelBeamEvsMipMax[1][v][gs][l][cb]->Draw("p,e,same");
0534                 if (lentered[cb] == 0){
0535                   legend->AddEntry(graphPerChannelBeamEvsMipMax[1][v][gs][l][cb], Form("%d",cb), "p");
0536                   lentered[cb]  = 1;
0537                   plotted       = 1;
0538                 }
0539               }
0540             }
0541             
0542             Int_t lowGainSet  = GetLGGainSetIndex(gs, tbdata);
0543             Int_t highGainSet = GetHGGainSetIndex(gs, tbdata);
0544             Float_t overV     = GetOverVoltageFromIndex(v, tbdata);
0545             if (lowGainSet == -1 || highGainSet == -1 || overV == -1) continue;
0546             DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0547             DrawLatex(0.95, 0.92-textSizeRel, Form("V_{ov} =%.1f V", overV), true, textSizeRel, 42);        
0548             DrawLatex(0.95, 0.92-2*textSizeRel, Form("HG amp: %d, LG amp: %d", highGainSet, lowGainSet), true, textSizeRel, 42);
0549             legend->Draw();
0550             
0551             if (plotted) canvas1DSimple->SaveAs(Form("%s/BeamEvsMPVlangaus_V%03.0fmV_Layer%02d_HG%02d_LG%02d.pdf",outputDirPlotsH.Data(), overV*1000,  l, highGainSet, lowGainSet));
0552           }
0553         }
0554       }
0555     }
0556 
0557     
0558     
0559     // ****************************************************************************************
0560     // Create electron beam plots
0561     // ****************************************************************************************
0562     for (Int_t l = 0; l < gMaxLayers; l++ ){
0563       if (layerExists[l]){
0564         for (Int_t gs = 0; gs < gMaxGainS; gs++){
0565             for (Int_t be = 0; be < gMaxBEn; be++){
0566             canvas1DSimple->cd();
0567             canvas1DSimple->Clear();
0568             histDummyVvsMPV->Draw();
0569             Bool_t lentered[9]  = {0, 0, 0, 0, 0, 0, 0, 0, 0};
0570             Bool_t plotted      = 0; 
0571             TLegend* legend = GetAndSetLegend2( 0.12, 0.95-3*textSizeRel, 0.4, 0.95,textSizeRel, 4, Form("Layer %d, channel:",l), 42,0.2);
0572           
0573             for (Int_t cb = 1; cb < 9; cb++){
0574               if (graphPerChannelVoVvsMipMax[0][be][gs][l][cb]->GetN() > 0){
0575                 SetMarkerDefaultsTGraphErr( graphPerChannelVoVvsMipMax[0][be][gs][l][cb], markerReadBoard[cb-1], 0.8, colorReadBoard[cb-1],colorReadBoard[cb-1]);
0576                 graphPerChannelVoVvsMipMax[0][be][gs][l][cb]->Draw("p,e,same");
0577                 if (lentered[cb] == 0){
0578                   legend->AddEntry(graphPerChannelVoVvsMipMax[0][be][gs][l][cb], Form("%d",cb), "p");
0579                   lentered[cb]  = 1;
0580                   plotted       = 1;
0581                 }
0582               }
0583             }
0584             
0585             Int_t lowGainSet  = GetLGGainSetIndex(gs, tbdata);
0586             Int_t highGainSet = GetHGGainSetIndex(gs, tbdata);
0587             Float_t beamEnergy = GetBeamEnergyFromIndex(be, tbdata);
0588             DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0589             DrawLatex(0.95, 0.92-textSizeRel, Form("e-beam E =%.0f GeV", beamEnergy), true, textSizeRel, 42);        
0590             DrawLatex(0.95, 0.92-2*textSizeRel, Form("HG amp: %d, LG amp: %d", highGainSet, lowGainSet), true, textSizeRel, 42);
0591             legend->Draw();
0592             
0593             if (plotted) canvas1DSimple->SaveAs(Form("%s/VopVsMPVlangaus_Ebeam%03.0fGeV_Layer%02d_HG%02d_LG%02d.pdf",outputDirPlotsE.Data(), beamEnergy,  l, highGainSet, lowGainSet));
0594           }
0595         }
0596       }
0597     }
0598     for (Int_t l = 0; l < gMaxLayers; l++ ){
0599       if (layerExists[l]){
0600         for (Int_t gs = 0; gs < gMaxGainS; gs++){
0601           for (Int_t v = 0; v < gMaxVoltages; v++){
0602             canvas1DSimple->cd();
0603             canvas1DSimple->Clear();
0604             histDummyBeamEvsMPV->Draw();
0605             Bool_t lentered[9]  = {0, 0, 0, 0, 0, 0, 0, 0, 0};
0606             Bool_t plotted      = 0; 
0607             TLegend* legend = GetAndSetLegend2( 0.12, 0.95-3*textSizeRel, 0.4, 0.95,textSizeRel, 4, Form("Layer %d, channel:",l), 42,0.2);
0608           
0609             for (Int_t cb = 1; cb < 9; cb++){
0610               if (graphPerChannelBeamEvsMipMax[0][v][gs][l][cb]->GetN() > 0){
0611                 SetMarkerDefaultsTGraphErr( graphPerChannelBeamEvsMipMax[0][v][gs][l][cb], markerReadBoard[cb-1], 0.8, colorReadBoard[cb-1],colorReadBoard[cb-1]);
0612                 graphPerChannelBeamEvsMipMax[0][v][gs][l][cb]->Draw("p,e,same");
0613                 if (lentered[cb] == 0){
0614                   legend->AddEntry(graphPerChannelBeamEvsMipMax[0][v][gs][l][cb], Form("%d",cb), "p");
0615                   lentered[cb]  = 1;
0616                   plotted       = 1;
0617                 }
0618               }
0619             }
0620             
0621             Int_t lowGainSet  = GetLGGainSetIndex(gs, tbdata);
0622             Int_t highGainSet = GetHGGainSetIndex(gs, tbdata);
0623             Float_t overV     = GetOverVoltageFromIndex(v, tbdata);
0624             if (lowGainSet == -1 || highGainSet == -1 || overV == -1) continue;
0625             DrawLatex(0.95, 0.92, Form("%s", labelTB.Data()), true, textSizeRel, 42);        
0626             DrawLatex(0.95, 0.92-textSizeRel, Form("V_{ov} =%.1f V", overV), true, textSizeRel, 42);        
0627             DrawLatex(0.95, 0.92-2*textSizeRel, Form("HG amp: %d, LG amp: %d", highGainSet, lowGainSet), true, textSizeRel, 42);
0628             legend->Draw();
0629             
0630             if (plotted) canvas1DSimple->SaveAs(Form("%s/BeamEvsMPVlangaus_V%03.0fmV_Layer%02d_HG%02d_LG%02d.pdf",outputDirPlotsE.Data(), overV*1000,  l, highGainSet, lowGainSet));
0631           }
0632         }
0633       }
0634     }
0635 
0636     
0637 }