quicklook/QLflightS4_ND/S4_QL.cpp


/*
* S4 Quick Look
* Author Marcelli-Malvezzi
* Version 1.00 - March 2006
* 
* Description - The aim of S4 QL software is to monitor the bahaviour of this detector.
*               It creates three canvases: the first one contains the histogram of S4 data, the second one 
*               is relative to the time behaviour of collected data while the last shows S4 rate 
*               (information from trigger Packet).
*               See documentation for a more detailed description of the output.
*
*
* Parameters:
*               TSTring base    - the path to the root directory for the specific Pamela unpack session
*                                 There is no default value, without this input the program will not run
*               TString outDir  - the path where to save the output image (Default = ./) 
*               TString format  - the format which will be used fo rsave the produced images (Default = "jpg") 
*               Float_t DeltaT  - the time interval in minute for calculation of average S4 data for minute, 
*                                 see S4_QL_2 plot (Default = 1 minute)
*
*
* Version 1.1  - June 2006
* Fixed bugs: the vector "trcss" was inizilized to have dimension 10, but for files with large number of events this
*             is not true; now this vector is inizialized at 100.
*
*             the threshold at which S4 is set and the trigger configuration can change in the file; all these changes are reported
*             in a pad
*
*             for a large namber of events is not possible to have vectors, so all graphs have been converted in histograms
*******/


#include <iostream>
#include <fstream>
#include <sstream>
#include <math.h>
#include "TStyle.h"
#include "TFile.h"
#include "TList.h"
#include "TTree.h"
#include "TLatex.h"
#include "TObjString.h"
#include "TCanvas.h"
#include "TGraph.h"
#include "TH1F.h"
#include "TF1.h"
#include "TGaxis.h"
#include "TString.h"
#include "TPaveText.h"
#include "EventHeader.h"
#include "PscuHeader.h"
#include "TMultiGraph.h"
#include "physics/S4/S4Event.h"
#include "varDump/VarDumpEvent.h"
#include "varDump/VarDumpRecord.h"
#include "physics/trigger/TriggerEvent.h"

using namespace std;

void S4_QL(TString base, TString outDir, TString format, ULong_t DeltaT){  //DeltaT in minute

  //------ Variables initialization ---------/
  Int_t tmpSize;
  ULong_t mintime, maxtime;
  Int_t   adcmax;
  Int_t   j=0;
  Int_t   S4_TRHadc;
  Int_t   S4_TRHmip;
  TString filename = ((TObjString*)base.Tokenize('/')->Last())->GetString();
  filename = ((TObjString*)filename.Tokenize('.')->First())->GetString();
  char *trc;
  TString trcs;
  TString str;
  TString trcsstot[31];
  TString trcsstot2;
  TString trcss[100]=""; 
  TString trgconf[31]={"TOF1","TOF2","TOF3","TOF4","TOF5","TOF6","TOF7","S4","CALO","CALO-S4","TOF1-S4","TOF2-S4","TOF3-S4","TOF4-S4","TOF5-S4","TOF6-S4","TOF7-S4","TOF1-CALO","TOF2-CALO","TOF3-CALO","TOF4-CALO","TOF5-CALO","TOF6-CALO","TOF7-CALO","TOF1-CALO-S4","TOF2-CALO-S4","TOF3-CALO-S4","TOF4-CALO-S4","TOF5-CALO-S4","TOF6-CALO-S4","TOF7-CALO-S4",};
  stringstream oss, oss1, oss2, oss3, s4soglia, buffer, conftrig, noentries;
  ULong_t lastime, firstime;
  Int_t vardumpentries = 0;
  char S4_TRH[10];
  char S4_TRH2[10];
  Int_t trigconf = 0;
  Long64_t nevents;
  string title;
  double obmin=0.;
  double obmax=0.; 
  double obt; 
  double s4rate;
  int nbin=0;
  //------to open headerfile, s4file, Trigfile and vardumpfile---------------------------/
  TFile  *file  =new TFile(base.Data()) ;
  if (!file){
    printf("No such file in the directory has been found");
    return;
  }
  
  TTree *VarDumpTr = (TTree*)file->Get("VarDump");
  TTree *PhysicsTr = (TTree*)file->Get("Physics");

  TBranch *S4Br = PhysicsTr->GetBranch("S4");
  TBranch *TriggerBr = PhysicsTr->GetBranch("Trigger");
  TBranch *headBr = PhysicsTr->GetBranch("Header");
  TBranch *VarDumpBr = VarDumpTr->GetBranch("VarDump");
  TBranch *headVarDumpBr = VarDumpTr->GetBranch("Header");

  pamela::S4::S4Event  *s4e  = 0;
  pamela::EventHeader  *eh   = 0;
  pamela::PscuHeader   *ph   = 0;
  pamela::trigger::TriggerEvent *trige = 0;
  pamela::VarDumpEvent  *vde = 0; 
  pamela::VarDumpRecord *vdr = 0;
 
  PhysicsTr->SetBranchAddress("S4", &s4e);
  PhysicsTr->SetBranchAddress("Header", &eh);
  PhysicsTr->SetBranchAddress("Trigger", &trige);
  VarDumpTr->SetBranchAddress("VarDump", &vde);

  nevents = headBr->GetEntries();
  const Int_t sizetot = nevents;

  //----------- If nevents < = 0 ---------------------------------/
  if (nevents<=0) {
    printf("nevents = %i \n", nevents);
    printf(" \n");

    TCanvas *canvas4 = new TCanvas("No entries", "No entries ", 400, 200);
    canvas4->SetFillColor(10);
    canvas4->cd();

    TLatex *l = new TLatex();
    l->SetTextAlign(12);
    l->SetTextSize(0.15);
    l->SetTextColor(2);
    noentries.str("");
    noentries<< "S4_QL:";
    l->DrawLatex(0.05, 0.7, noentries.str().c_str()); 
    noentries.str("");
    noentries<< "No entries for this files";
    l->DrawLatex(0.05, 0.5, noentries.str().c_str()); 

    if (outDir == "./") {
      oss.str("");
      oss << filename.Data() << "_S4_QL." << format.Data();
    } else {
      oss.str("");
      oss << outDir.Data() << filename.Data() << "_S4_QL." << format.Data();
    }
  
    canvas4->Update();
    canvas4->SaveAs(oss.str().c_str()); 

    return;
  }

  for (Int_t i = 0; i < nevents; i++){
    headBr->GetEntry(i);
    S4Br->GetEntry(i);
    TriggerBr->GetEntry(i);
    if (s4e->unpackError == 1 || s4e->S4_DATA==0)  continue;      
    //-------to set max adc value------------------------------------------//
    if (i==0){
      adcmax=s4e->S4_DATA; 
    }
    //-------to set trigger configuration---------------------------------//
    if ((s4e->S4_DATA)> adcmax) adcmax=(s4e->S4_DATA);
    trigconf = trige->trigconf;
    trc = 0;
    if ( trigconf & (1<<0) ) trc = "TOF1"; 
    if ( trigconf & (1<<1) ) { 
      if (trc==0) trc= "TOF2";
      else trc = Form("%s-TOF2",trc);
    }
    if ( trigconf & (1<<2) ) { 
      if (trc==0) trc= "TOF3";
      else trc = Form("%s-TOF3",trc);
    } 
    if ( trigconf & (1<<3) ) { 
      if (trc==0) trc= "TOF4";
      else trc = Form("%s-TOF4",trc);
    } 
    if ( trigconf & (1<<4) ) { 
      if (trc==0) trc= "TOF5";
      else trc = Form("%s-TOF5",trc);
    } 
    if ( trigconf & (1<<5) ) { 
      if (trc==0) trc= "TOF6";
      else trc = Form("%s-TOF6",trc);
    } 
    if ( trigconf & (1<<6) ) { 
      if (trc==0) trc= "TOF7";
      else trc = Form("%s-TOF7",trc);
    } 
    if ( trigconf & (1<<7) ) { 
      if (trc==0) trc= "S4";
      else trc = Form("%s-S4",trc);
    } 
    if ( trigconf & (1<<8) ) { 
      if (trc==0) trc= "CALO";
      else trc = Form("%s-CALO",trc);
    } 
    if ( trigconf & (1<<9) ) { 
      if (trc==0) trc= "CALIB_ON";
      else trc = Form("%s-CALIB_ON",trc);
    } 
    trcs = "";
    trcs = trc;
    if (i==0){ 
      trcss[0]=trcs; 
      j=j+1;
    }
    if (trcs!=trcss[j-1] && i>0){ 
      trcss[j]=trcs;
      j=j+1; 
    }
  }
  for(Int_t k=0; k<31; k++){
    for(Int_t p=0; p<j ; p++){
      if(trgconf[k]==trcss[p]){
        trcsstot[k]=trgconf[k];
      }
    }
  }
  for(Int_t k=0; k<31; k++){
    if(trcsstot[k]!= "")
      trcsstot2=trcsstot2+"/"+ trcsstot[k];
  }
  //---------to search S4 threshold and convert it into char (hexadecimal value)--------//
  vardumpentries =  VarDumpBr->GetEntries();
  if(vardumpentries==0){
    cout<<"TTree VarDump: Entries = 0"<<"\n"<<"         S4 Threshold not defined"<<"\n";
    str= " Not defined";
  }
  else{
    for(Int_t k=0; k< vardumpentries; k++){
      VarDumpBr->GetEntry(k);
      vdr = (pamela::VarDumpRecord*)vde->Records->At(118);
      S4_TRHadc=((vdr->VAR_VALUE)/256);      //trh (number of ADC channels) 
      S4_TRHmip=S4_TRHadc*2;  

      //-----to convert decimal value to TString of S4 threshold value-------//
      if(k==0){
        sprintf(S4_TRH, "%d" , S4_TRHmip);
        str= S4_TRH;
      }else{
        sprintf(S4_TRH2, "%d" , S4_TRHmip);
        if(!strcmp(S4_TRH, S4_TRH2)) continue;
        str= str+"/"+S4_TRH2;   
      }
      sprintf(S4_TRH, "%d" , S4_TRHmip);
    }
  } 

  //********************** First Histogram ************************************//
  TH1F *h1 = new TH1F("All events", "S4 distribution for file: "+filename, adcmax, 10, adcmax+(adcmax/10));
  TH1F *h2 = new TH1F("S4 Triggered events", "S4 distribution for file: "+filename, adcmax, 10, adcmax+(adcmax/10));

  //********************** Second Histogram ************************************//
  headBr->GetEntry(0);
  ph = eh->GetPscuHeader();
  firstime = ph->GetOrbitalTime();
  headBr->GetEntry(nevents-1);
  ph = eh->GetPscuHeader();
  lastime = ph->GetOrbitalTime();
  obmin=firstime;
  obmax=lastime;   
  const ULong_t nint=(((lastime-firstime)/(DeltaT*60000)));
  const Int_t size = (Int_t)(nint+1); 

  TH1F *Allev = new TH1F("Mean signal from S4-all events", filename+": Mean signal from S4 (all triggered events)", size, obmin, obmax);
  TH1F *Alltime = new TH1F("Mean signal from S4-all events", filename+": Mean signal from S4 (all triggered events)", size, obmin, obmax);
  TH1F *S4ev = new TH1F("Mean signal from S4-triggered events", filename+": Mean signal from S4 (only s4 triggered events)", size, obmin, obmax);
  TH1F *S4time = new TH1F("Mean signal from S4-triggered events", filename+": Mean signal from S4 (only s4 triggered events)", size, obmin, obmax);  

  //********************** Third Histogram ************************************//
  nbin=sizetot;
  title="";
  title=filename+": S4 rate from Trigger Packet";

  const ULong_t nint3=(lastime-firstime);
  const Int_t size3 = (Int_t)((nint3)/3); 

  TH1F *rate= new TH1F(title.c_str(), title.c_str(), nint3, obmin, obmax);
  TH1F *rateline= new TH1F(filename+": smoothed S4 rate", filename+": smoothed S4 rate", size3, obmin, obmax);

  //------------------------------------------------------------------------------------------------------------------
  //------- fill histograms ---------//
  Int_t n=0, p=0;
  for (Int_t i = 0; i < nevents; i++){
    TriggerBr->GetEntry(i);
    headBr->GetEntry(i);
    S4Br->GetEntry(i);
    ph = eh->GetPscuHeader();
    obt = ph->GetOrbitalTime();  
    if (s4e->unpackError == 1 && (s4e->S4_DATA)==0)  continue;   
    s4rate= trige->s4calcount[0];
    rate->Fill(obt, s4rate);
    h1->Fill(s4e->S4_DATA);  
    Allev->Fill(obt,s4e->S4_DATA);
    Alltime->Fill(obt);
    if((i==n) && (s4rate >0)){
      rateline->Fill(obt,s4rate);
    n=n+5;
      }
    if ((trige->patterntrig[0] == 0) && (trige->patterntrig[2] == 0) && (trige->patterntrig[3] == 0) && (trige->patterntrig[4] == 0) &&     (trige->patterntrig[5] == 0)){
      h2->Fill(s4e->S4_DATA);
      S4ev->Fill(obt,s4e->S4_DATA);
      S4time->Fill(obt);
      p=p+1;
    }
  }


  //****************************** Canvases *******************************//
  //------------------- First Canvas --------------------------------//
  TCanvas *canvas1 = new TCanvas("S4_QL_1", "S4 HISTO ", 1280, 1024);
  canvas1->SetFillColor(10);
  canvas1->Divide(1,2);

  TPad *all= new TPad ("","", 0, 0, 1, 1);
  all->SetFillColor(10);
  TPad *s4 = new TPad ("s4","s4", 0, 1, 1, 0); 
  s4->SetFillColor(10);
  TPad *pad = new TPad("pad","pad", .30,.80,.78,.90);
  pad->SetFillColor(10);
  TPad *pad2 = new TPad("pad2","pad2", .30,.80,.78,.90); 
  pad2->SetFillColor(10);
  TLatex *l = new TLatex();
  l->SetTextAlign(12);
  l->SetTextSize(0.35);
  l->SetTextColor(8);
  s4soglia.str("");
  s4soglia << "S4_THRESHOLD: "<<str.Data()<<" MIP";
  TLatex *ll = new TLatex();
  ll->SetTextAlign(12);
  ll->SetTextSize(0.35);
  conftrig.str("");
  conftrig<<"TRIGGER CONF: "<<trcsstot2.Data();
  if(adcmax<=100) adcmax=120;

  canvas1->cd(1);
  all->Draw();
  all->cd();
  h1->SetLineColor(kBlack);
  h1->SetFillColor(kRed);
  h1->GetXaxis()->SetTitle("ADC");
  h1->GetXaxis()->CenterTitle();
  h1->GetYaxis()->SetTitle("Number of Events");
  h1->GetYaxis()->CenterTitle(); 
  h1->Draw();
  pad->Draw();
  pad->cd();
  l->DrawLatex(0.05, 0.65, s4soglia.str().c_str());
  ll->DrawLatex(0.05, 0.15, conftrig.str().c_str());


  canvas1->cd(2);
  s4->Draw();
  s4->cd();
  h2->SetLineColor(kBlack);
  h2->SetFillColor(2);
  h2->GetXaxis()->SetTitle("ADC");
  h2->GetXaxis()->CenterTitle();
  h2->GetYaxis()->SetTitle("Number of Events");
  h2->GetYaxis()->CenterTitle(); 
  h2->Draw();
  pad2->Draw();
  pad2->cd();
  l->DrawLatex(0.05, 0.65, s4soglia.str().c_str()); 
  ll->DrawLatex(0.05, 0.15, conftrig.str().c_str());

  if (h1->GetMaximum()>0){
    all->SetLogy();
    all->SetLogx();
  } 
  if (h2->GetMaximum()>0){
    s4->SetLogy();
    s4->SetLogx();
  }
  canvas1->Update();
 

  if (outDir == "./") {
    oss1.str("");
    oss1 << filename.Data() << "_S4_QL_1." << format.Data();
  } else {
    oss1.str("");
    oss1 << outDir.Data() << filename.Data() << "_S4_QL_1." << format.Data(); 
  }
  canvas1->SaveAs(oss1.str().c_str()); 
 
  //------------------- Second Canvas --------------------------------//
  if (nint==0){
    cout<<"Number of Time Interval = 0"<<"\n"<<"      Set another value for input parameter DeltaTevtime "<<"\n"<<"RETURN"<<"\n";
    
    TCanvas *canvas2 = new TCanvas("Number of time interval=0", "Number of time interval=0", 900, 300);
    canvas2->SetFillColor(10);
    canvas2->cd();

    TLatex *l = new TLatex();
    l->SetTextAlign(12);
    l->SetTextSize(0.10);
    l->SetTextColor(2);
    noentries.str("");
    noentries<< "S4_QL: Time evolution";
    l->DrawLatex(0.05, 0.7, noentries.str().c_str()); 
    noentries.str("");
    noentries<< "Number of Time Interval = 0";
    l->SetTextColor(1);
    l->DrawLatex(0.05, 0.5, noentries.str().c_str()); 
    noentries.str("      Set another value for input parameter DeltaTevtime");
    noentries<< "";
    l->DrawLatex(0.05, 0.3, noentries.str().c_str()); 

    if (outDir == "./") {
      oss.str("");
      oss << filename.Data() << "_S4_QL_2." << format.Data();
    } else {
      oss.str("");
      oss << outDir.Data() << filename.Data() << "_S4_QL_2." << format.Data();
    }
    canvas2->Update();
    canvas2->SaveAs(oss.str().c_str()); 

    return;
  }
  TCanvas *canvas2 = new TCanvas("S4_QL_2", "S4 - Time evolution", 1280, 1024);
  canvas2->SetFillColor(10);
  canvas2->Divide(1,2);
  
  TPad *pad3 = new TPad("pad3","pad3", .7,.1,.88,.30);
  pad3->SetFillColor(10);
  TPad *pad4 = new TPad("pad4","pad4", .7,.1,.88,.30); 
  pad4->SetFillColor(10);

  TLatex *l3 = new TLatex();
  l3->SetTextAlign(12);
  l3->SetTextSize(.3);
  l3->SetTextColor(2);
  buffer.str("");
  buffer<<"DT :  "<<DeltaT<<" min";
  
  canvas2->cd(1);
  Allev->SetStats(kFALSE);
  Allev->SetMarkerColor(2);
  Allev->SetMarkerSize(.5);
  Allev->SetMarkerStyle(21);
  Allev->GetXaxis()->SetTitle("Time ( OBT in ms )");
  Allev->GetXaxis()->CenterTitle();
  Allev->GetYaxis()->SetTitle("Mean signal ( ADC )");
  Allev->GetYaxis()->CenterTitle();
  Allev->SetMinimum(20);
  Allev->Divide(Alltime);
  Allev->Draw("p");
  pad3->Draw();
  pad3->cd();
  l3->DrawLatex(0.05, 0.65, buffer.str().c_str());

  
  canvas2->cd(2);
  if(p != 0){
  S4ev->SetStats(kFALSE);
  S4ev->SetMarkerColor(2);
  S4ev->SetMarkerSize(.5);
  S4ev->SetMarkerStyle(21);
  S4ev->GetXaxis()->SetTitle("Time ( OBT in ms )");
  S4ev->GetXaxis()->CenterTitle();
  S4ev->GetYaxis()->SetTitle("Mean signal ( ADC )");
  S4ev->GetYaxis()->CenterTitle();
  S4ev->SetMinimum(20);
  S4ev->Divide(S4time);
  S4ev->Draw("p");
  pad4->Draw();
  pad4->cd();
  l3->DrawLatex(0.05, 0.65, buffer.str().c_str());
  }else{
    TLatex *l = new TLatex();
    l->SetTextAlign(12);
    l->SetTextSize(0.07);
    l->SetTextColor(2);
    noentries.str("");
    noentries<< "S4_QL: Time evolution";
    l->DrawLatex(0.06, 0.6, noentries.str().c_str()); 
    noentries.str("");
    noentries<< "No S4 triggered events.";
    l->SetTextColor(1);
    l->DrawLatex(0.06, 0.5, noentries.str().c_str()); 
  }

  if (outDir == "./") {
    oss2.str("");
    oss2 << filename.Data() << "_S4_QL_2." << format.Data();
  } else {
    oss2.str("");
    oss2 << outDir.Data() << filename.Data() << "_S4_QL_2." << format.Data();
  }
  canvas2->SaveAs(oss2.str().c_str());

  //------------------- Third Canvas --------------------------------//
  TCanvas *canvas3 = new TCanvas("S4 - Rate", "S4 - Rate", 1280, 1024);
  canvas3->SetFillColor(10);
  canvas3->Divide(1,2);

  canvas3->cd(1);
  rate->SetStats(kFALSE);
  rate->SetMarkerColor(2);
  rate->SetMarkerSize(.5);
  rate->SetMarkerStyle(21);
  rate->GetXaxis()->SetTitle("Time ( OBT in ms )");
  rate->GetXaxis()->CenterTitle();
  rate->GetYaxis()->SetTitle("S4 rate (Hz)");
  rate->GetYaxis()->CenterTitle();
  if(rate->GetMaximum() > 1000) gPad->SetLogy();
  rate->Draw("9p");

  canvas3->cd(2);
  rateline->SetStats(kFALSE);
  rateline->SetMarkerColor(4);
  rateline->SetMarkerSize(.5);
  rateline->SetMarkerStyle(21);
  rateline->GetXaxis()->SetTitle("Time ( OBT in ms )");
  rateline->GetXaxis()->CenterTitle();
  rateline->GetYaxis()->SetTitle("S4 rate (Hz)");
  rateline->GetYaxis()->CenterTitle();
  rateline->SetMaximum(rate->GetMaximum());
  if(rateline->GetMaximum() > 1000) gPad->SetLogy();
  rateline->Draw("9p");
 

  if (outDir == "./") {
    oss3.str("");
    oss3 << filename.Data() << "_S4_QL_3." << format.Data();
  } else {
    oss3.str("");
    oss3 << outDir.Data() << filename.Data() << "_S4_QL_3." << format.Data();
  }
  canvas3->SaveAs(oss3.str().c_str());
  file->Close();
}


int main(int argc, char* argv[]){
  TString path;
  TString outDir  ="./";
  TString format  ="jpg";
  ULong_t DeltaT  =1;
   
    
 if (argc < 2){
    printf("You have to insert at least the file to analyze \n");
    printf("Try '--help' for more information. \n");
    exit(1);
  }  

  if (!strcmp(argv[1], "--help")){
        printf( "Usage: ND_QL FILE [OPTION] \n");
        printf( "\t --help                 Print this help and exit \n");
        printf( "\t -outDir[path]          Path where to put the output [default ./] \n");
        printf( "\t -format[jpg|gif|ps]    Format for output files [default 'jpg'] \n");
        printf( "\t -DeltaT                Time interval to control time evolution of acquired data in minute [default = 1 min] \n");
        exit(1);
  }

  path=argv[1];

  for (int i = 2; i < argc; i++){
    
    if (!strcmp(argv[i], "-outDir")){
      if (++i >= argc){
        printf( "-outDir needs arguments. \n");
        printf( "Try '--help' for more information. \n");
        exit(1);
      } 
      else{
        outDir = argv[i];
        continue;
      }
    }

    if (!strcmp(argv[i], "-format")){
      if (++i >= argc){
        printf( "-format needs arguments. \n");
        printf( "Try '--help' for more information. \n");
        exit(1);
      }
      else{
        format = argv[i];
        continue;
      }
    }

    if (!strcmp(argv[i], "-DeltaT")){
      if (++i >= argc){
        printf( "-DeltaT needs arguments. \n");
        printf( "Try '--help' for more information. \n");
        exit(1);
      } 
      else{
        DeltaT = atol(argv[i]);
        continue;
      }
    }


  }

   S4_QL(argv[1], outDir, format, DeltaT);

}

1	pam-rm2	1.1
2			/*
3			* S4 Quick Look
4			* Author Marcelli-Malvezzi
5			* Version 1.00 - March 2006
6			*
7			* Description - The aim of S4 QL software is to monitor the bahaviour of this detector.
8			* It creates three canvases: the first one contains the histogram of S4 data, the second one
9			* is relative to the time behaviour of collected data while the last shows S4 rate
10			* (information from trigger Packet).
11			* See documentation for a more detailed description of the output.
12			*
13			*
14			* Parameters:
15			* TSTring base - the path to the root directory for the specific Pamela unpack session
16			* There is no default value, without this input the program will not run
17			* TString outDir - the path where to save the output image (Default = ./)
18			* TString format - the format which will be used fo rsave the produced images (Default = "jpg")
19			* Float_t DeltaT - the time interval in minute for calculation of average S4 data for minute,
20			* see S4_QL_2 plot (Default = 1 minute)
21			*
22			*
23			* Version 1.1 - June 2006
24			* Fixed bugs: the vector "trcss" was inizilized to have dimension 10, but for files with large number of events this
25			* is not true; now this vector is inizialized at 100.
26			*
27			* the threshold at which S4 is set and the trigger configuration can change in the file; all these changes are reported
28			* in a pad
29			*
30			* for a large namber of events is not possible to have vectors, so all graphs have been converted in histograms
31			*******/
32
33
34			#include <iostream>
35			#include <fstream>
36			#include <sstream>
37			#include <math.h>
38			#include "TStyle.h"
39			#include "TFile.h"
40			#include "TList.h"
41			#include "TTree.h"
42			#include "TLatex.h"
43			#include "TObjString.h"
44			#include "TCanvas.h"
45			#include "TGraph.h"
46			#include "TH1F.h"
47			#include "TF1.h"
48			#include "TGaxis.h"
49			#include "TString.h"
50			#include "TPaveText.h"
51			#include "EventHeader.h"
52			#include "PscuHeader.h"
53			#include "TMultiGraph.h"
54			#include "physics/S4/S4Event.h"
55			#include "varDump/VarDumpEvent.h"
56			#include "varDump/VarDumpRecord.h"
57			#include "physics/trigger/TriggerEvent.h"
58
59			using namespace std;
60
61			void S4_QL(TString base, TString outDir, TString format, ULong_t DeltaT){ //DeltaT in minute
62
63			//------ Variables initialization ---------/
64			Int_t tmpSize;
65			ULong_t mintime, maxtime;
66			Int_t adcmax;
67			Int_t j=0;
68			Int_t S4_TRHadc;
69			Int_t S4_TRHmip;
70			TString filename = ((TObjString*)base.Tokenize('/')->Last())->GetString();
71			filename = ((TObjString*)filename.Tokenize('.')->First())->GetString();
72			char *trc;
73			TString trcs;
74			TString str;
75			TString trcsstot[31];
76			TString trcsstot2;
77			TString trcss[100]="";
78			TString trgconf[31]={"TOF1","TOF2","TOF3","TOF4","TOF5","TOF6","TOF7","S4","CALO","CALO-S4","TOF1-S4","TOF2-S4","TOF3-S4","TOF4-S4","TOF5-S4","TOF6-S4","TOF7-S4","TOF1-CALO","TOF2-CALO","TOF3-CALO","TOF4-CALO","TOF5-CALO","TOF6-CALO","TOF7-CALO","TOF1-CALO-S4","TOF2-CALO-S4","TOF3-CALO-S4","TOF4-CALO-S4","TOF5-CALO-S4","TOF6-CALO-S4","TOF7-CALO-S4",};
79			stringstream oss, oss1, oss2, oss3, s4soglia, buffer, conftrig, noentries;
80			ULong_t lastime, firstime;
81			Int_t vardumpentries = 0;
82			char S4_TRH[10];
83			char S4_TRH2[10];
84			Int_t trigconf = 0;
85			Long64_t nevents;
86			string title;
87			double obmin=0.;
88			double obmax=0.;
89			double obt;
90			double s4rate;
91			int nbin=0;
92			//------to open headerfile, s4file, Trigfile and vardumpfile---------------------------/
93			TFile *file =new TFile(base.Data()) ;
94			if (!file){
95			printf("No such file in the directory has been found");
96			return;
97			}
98
99			TTree VarDumpTr = (TTree)file->Get("VarDump");
100			TTree PhysicsTr = (TTree)file->Get("Physics");
101
102			TBranch *S4Br = PhysicsTr->GetBranch("S4");
103			TBranch *TriggerBr = PhysicsTr->GetBranch("Trigger");
104			TBranch *headBr = PhysicsTr->GetBranch("Header");
105			TBranch *VarDumpBr = VarDumpTr->GetBranch("VarDump");
106			TBranch *headVarDumpBr = VarDumpTr->GetBranch("Header");
107
108			pamela::S4::S4Event *s4e = 0;
109			pamela::EventHeader *eh = 0;
110			pamela::PscuHeader *ph = 0;
111			pamela::trigger::TriggerEvent *trige = 0;
112			pamela::VarDumpEvent *vde = 0;
113			pamela::VarDumpRecord *vdr = 0;
114
115			PhysicsTr->SetBranchAddress("S4", &s4e);
116			PhysicsTr->SetBranchAddress("Header", &eh);
117			PhysicsTr->SetBranchAddress("Trigger", &trige);
118			VarDumpTr->SetBranchAddress("VarDump", &vde);
119
120			nevents = headBr->GetEntries();
121			const Int_t sizetot = nevents;
122
123			//----------- If nevents < = 0 ---------------------------------/
124			if (nevents<=0) {
125			printf("nevents = %i \n", nevents);
126			printf(" \n");
127
128			TCanvas *canvas4 = new TCanvas("No entries", "No entries ", 400, 200);
129			canvas4->SetFillColor(10);
130			canvas4->cd();
131
132			TLatex *l = new TLatex();
133			l->SetTextAlign(12);
134			l->SetTextSize(0.15);
135			l->SetTextColor(2);
136			noentries.str("");
137			noentries<< "S4_QL:";
138			l->DrawLatex(0.05, 0.7, noentries.str().c_str());
139			noentries.str("");
140			noentries<< "No entries for this files";
141			l->DrawLatex(0.05, 0.5, noentries.str().c_str());
142
143			if (outDir == "./") {
144			oss.str("");
145			oss << filename.Data() << "_S4_QL." << format.Data();
146			} else {
147			oss.str("");
148			oss << outDir.Data() << filename.Data() << "_S4_QL." << format.Data();
149			}
150
151			canvas4->Update();
152			canvas4->SaveAs(oss.str().c_str());
153
154			return;
155			}
156
157			for (Int_t i = 0; i < nevents; i++){
158			headBr->GetEntry(i);
159			S4Br->GetEntry(i);
160			TriggerBr->GetEntry(i);
161			if (s4e->unpackError == 1 \|\| s4e->S4_DATA==0) continue;
162			//-------to set max adc value------------------------------------------//
163			if (i==0){
164			adcmax=s4e->S4_DATA;
165			}
166			//-------to set trigger configuration---------------------------------//
167			if ((s4e->S4_DATA)> adcmax) adcmax=(s4e->S4_DATA);
168			trigconf = trige->trigconf;
169			trc = 0;
170			if ( trigconf & (1<<0) ) trc = "TOF1";
171			if ( trigconf & (1<<1) ) {
172			if (trc==0) trc= "TOF2";
173			else trc = Form("%s-TOF2",trc);
174			}
175			if ( trigconf & (1<<2) ) {
176			if (trc==0) trc= "TOF3";
177			else trc = Form("%s-TOF3",trc);
178			}
179			if ( trigconf & (1<<3) ) {
180			if (trc==0) trc= "TOF4";
181			else trc = Form("%s-TOF4",trc);
182			}
183			if ( trigconf & (1<<4) ) {
184			if (trc==0) trc= "TOF5";
185			else trc = Form("%s-TOF5",trc);
186			}
187			if ( trigconf & (1<<5) ) {
188			if (trc==0) trc= "TOF6";
189			else trc = Form("%s-TOF6",trc);
190			}
191			if ( trigconf & (1<<6) ) {
192			if (trc==0) trc= "TOF7";
193			else trc = Form("%s-TOF7",trc);
194			}
195			if ( trigconf & (1<<7) ) {
196			if (trc==0) trc= "S4";
197			else trc = Form("%s-S4",trc);
198			}
199			if ( trigconf & (1<<8) ) {
200			if (trc==0) trc= "CALO";
201			else trc = Form("%s-CALO",trc);
202			}
203			if ( trigconf & (1<<9) ) {
204			if (trc==0) trc= "CALIB_ON";
205			else trc = Form("%s-CALIB_ON",trc);
206			}
207			trcs = "";
208			trcs = trc;
209			if (i==0){
210			trcss[0]=trcs;
211			j=j+1;
212			}
213			if (trcs!=trcss[j-1] && i>0){
214			trcss[j]=trcs;
215			j=j+1;
216			}
217			}
218			for(Int_t k=0; k<31; k++){
219			for(Int_t p=0; p<j ; p++){
220			if(trgconf[k]==trcss[p]){
221			trcsstot[k]=trgconf[k];
222			}
223			}
224			}
225			for(Int_t k=0; k<31; k++){
226			if(trcsstot[k]!= "")
227			trcsstot2=trcsstot2+"/"+ trcsstot[k];
228			}
229			//---------to search S4 threshold and convert it into char (hexadecimal value)--------//
230			vardumpentries = VarDumpBr->GetEntries();
231			if(vardumpentries==0){
232			cout<<"TTree VarDump: Entries = 0"<<"\n"<<" S4 Threshold not defined"<<"\n";
233			str= " Not defined";
234			}
235			else{
236			for(Int_t k=0; k< vardumpentries; k++){
237			VarDumpBr->GetEntry(k);
238			vdr = (pamela::VarDumpRecord*)vde->Records->At(118);
239			S4_TRHadc=((vdr->VAR_VALUE)/256); //trh (number of ADC channels)
240			S4_TRHmip=S4_TRHadc*2;
241
242			//-----to convert decimal value to TString of S4 threshold value-------//
243			if(k==0){
244			sprintf(S4_TRH, "%d" , S4_TRHmip);
245			str= S4_TRH;
246			}else{
247			sprintf(S4_TRH2, "%d" , S4_TRHmip);
248			if(!strcmp(S4_TRH, S4_TRH2)) continue;
249			str= str+"/"+S4_TRH2;
250			}
251			sprintf(S4_TRH, "%d" , S4_TRHmip);
252			}
253			}
254
255			//******************** First Histogram **********************************//
256			TH1F *h1 = new TH1F("All events", "S4 distribution for file: "+filename, adcmax, 10, adcmax+(adcmax/10));
257			TH1F *h2 = new TH1F("S4 Triggered events", "S4 distribution for file: "+filename, adcmax, 10, adcmax+(adcmax/10));
258
259			//******************** Second Histogram **********************************//
260			headBr->GetEntry(0);
261			ph = eh->GetPscuHeader();
262			firstime = ph->GetOrbitalTime();
263			headBr->GetEntry(nevents-1);
264			ph = eh->GetPscuHeader();
265			lastime = ph->GetOrbitalTime();
266			obmin=firstime;
267			obmax=lastime;
268			const ULong_t nint=(((lastime-firstime)/(DeltaT*60000)));
269			const Int_t size = (Int_t)(nint+1);
270
271			TH1F *Allev = new TH1F("Mean signal from S4-all events", filename+": Mean signal from S4 (all triggered events)", size, obmin, obmax);
272			TH1F *Alltime = new TH1F("Mean signal from S4-all events", filename+": Mean signal from S4 (all triggered events)", size, obmin, obmax);
273			TH1F *S4ev = new TH1F("Mean signal from S4-triggered events", filename+": Mean signal from S4 (only s4 triggered events)", size, obmin, obmax);
274			TH1F *S4time = new TH1F("Mean signal from S4-triggered events", filename+": Mean signal from S4 (only s4 triggered events)", size, obmin, obmax);
275
276			//******************** Third Histogram **********************************//
277			nbin=sizetot;
278			title="";
279			title=filename+": S4 rate from Trigger Packet";
280
281			const ULong_t nint3=(lastime-firstime);
282			const Int_t size3 = (Int_t)((nint3)/3);
283
284			TH1F *rate= new TH1F(title.c_str(), title.c_str(), nint3, obmin, obmax);
285			TH1F *rateline= new TH1F(filename+": smoothed S4 rate", filename+": smoothed S4 rate", size3, obmin, obmax);
286
287			//------------------------------------------------------------------------------------------------------------------
288			//------- fill histograms ---------//
289			Int_t n=0, p=0;
290			for (Int_t i = 0; i < nevents; i++){
291			TriggerBr->GetEntry(i);
292			headBr->GetEntry(i);
293			S4Br->GetEntry(i);
294			ph = eh->GetPscuHeader();
295			obt = ph->GetOrbitalTime();
296			if (s4e->unpackError == 1 && (s4e->S4_DATA)==0) continue;
297			s4rate= trige->s4calcount[0];
298			rate->Fill(obt, s4rate);
299			h1->Fill(s4e->S4_DATA);
300			Allev->Fill(obt,s4e->S4_DATA);
301			Alltime->Fill(obt);
302			if((i==n) && (s4rate >0)){
303			rateline->Fill(obt,s4rate);
304			n=n+5;
305			}
306			if ((trige->patterntrig[0] == 0) && (trige->patterntrig[2] == 0) && (trige->patterntrig[3] == 0) && (trige->patterntrig[4] == 0) && (trige->patterntrig[5] == 0)){
307			h2->Fill(s4e->S4_DATA);
308			S4ev->Fill(obt,s4e->S4_DATA);
309			S4time->Fill(obt);
310			p=p+1;
311			}
312			}
313
314
315			//**************************** Canvases *****************************//
316			//------------------- First Canvas --------------------------------//
317			TCanvas *canvas1 = new TCanvas("S4_QL_1", "S4 HISTO ", 1280, 1024);
318			canvas1->SetFillColor(10);
319			canvas1->Divide(1,2);
320
321			TPad *all= new TPad ("","", 0, 0, 1, 1);
322			all->SetFillColor(10);
323			TPad *s4 = new TPad ("s4","s4", 0, 1, 1, 0);
324			s4->SetFillColor(10);
325			TPad *pad = new TPad("pad","pad", .30,.80,.78,.90);
326			pad->SetFillColor(10);
327			TPad *pad2 = new TPad("pad2","pad2", .30,.80,.78,.90);
328			pad2->SetFillColor(10);
329			TLatex *l = new TLatex();
330			l->SetTextAlign(12);
331			l->SetTextSize(0.35);
332			l->SetTextColor(8);
333			s4soglia.str("");
334			s4soglia << "S4_THRESHOLD: "<<str.Data()<<" MIP";
335			TLatex *ll = new TLatex();
336			ll->SetTextAlign(12);
337			ll->SetTextSize(0.35);
338			conftrig.str("");
339			conftrig<<"TRIGGER CONF: "<<trcsstot2.Data();
340			if(adcmax<=100) adcmax=120;
341
342			canvas1->cd(1);
343			all->Draw();
344			all->cd();
345			h1->SetLineColor(kBlack);
346			h1->SetFillColor(kRed);
347			h1->GetXaxis()->SetTitle("ADC");
348			h1->GetXaxis()->CenterTitle();
349			h1->GetYaxis()->SetTitle("Number of Events");
350			h1->GetYaxis()->CenterTitle();
351			h1->Draw();
352			pad->Draw();
353			pad->cd();
354			l->DrawLatex(0.05, 0.65, s4soglia.str().c_str());
355			ll->DrawLatex(0.05, 0.15, conftrig.str().c_str());
356
357
358			canvas1->cd(2);
359			s4->Draw();
360			s4->cd();
361			h2->SetLineColor(kBlack);
362			h2->SetFillColor(2);
363			h2->GetXaxis()->SetTitle("ADC");
364			h2->GetXaxis()->CenterTitle();
365			h2->GetYaxis()->SetTitle("Number of Events");
366			h2->GetYaxis()->CenterTitle();
367			h2->Draw();
368			pad2->Draw();
369			pad2->cd();
370			l->DrawLatex(0.05, 0.65, s4soglia.str().c_str());
371			ll->DrawLatex(0.05, 0.15, conftrig.str().c_str());
372
373			if (h1->GetMaximum()>0){
374			all->SetLogy();
375			all->SetLogx();
376			}
377			if (h2->GetMaximum()>0){
378			s4->SetLogy();
379			s4->SetLogx();
380			}
381			canvas1->Update();
382
383
384			if (outDir == "./") {
385			oss1.str("");
386			oss1 << filename.Data() << "_S4_QL_1." << format.Data();
387			} else {
388			oss1.str("");
389			oss1 << outDir.Data() << filename.Data() << "_S4_QL_1." << format.Data();
390			}
391			canvas1->SaveAs(oss1.str().c_str());
392
393			//------------------- Second Canvas --------------------------------//
394			if (nint==0){
395			cout<<"Number of Time Interval = 0"<<"\n"<<" Set another value for input parameter DeltaTevtime "<<"\n"<<"RETURN"<<"\n";
396
397			TCanvas *canvas2 = new TCanvas("Number of time interval=0", "Number of time interval=0", 900, 300);
398			canvas2->SetFillColor(10);
399			canvas2->cd();
400
401			TLatex *l = new TLatex();
402			l->SetTextAlign(12);
403			l->SetTextSize(0.10);
404			l->SetTextColor(2);
405			noentries.str("");
406			noentries<< "S4_QL: Time evolution";
407			l->DrawLatex(0.05, 0.7, noentries.str().c_str());
408			noentries.str("");
409			noentries<< "Number of Time Interval = 0";
410			l->SetTextColor(1);
411			l->DrawLatex(0.05, 0.5, noentries.str().c_str());
412			noentries.str(" Set another value for input parameter DeltaTevtime");
413			noentries<< "";
414			l->DrawLatex(0.05, 0.3, noentries.str().c_str());
415
416			if (outDir == "./") {
417			oss.str("");
418			oss << filename.Data() << "_S4_QL_2." << format.Data();
419			} else {
420			oss.str("");
421			oss << outDir.Data() << filename.Data() << "_S4_QL_2." << format.Data();
422			}
423			canvas2->Update();
424			canvas2->SaveAs(oss.str().c_str());
425
426			return;
427			}
428			TCanvas *canvas2 = new TCanvas("S4_QL_2", "S4 - Time evolution", 1280, 1024);
429			canvas2->SetFillColor(10);
430			canvas2->Divide(1,2);
431
432			TPad *pad3 = new TPad("pad3","pad3", .7,.1,.88,.30);
433			pad3->SetFillColor(10);
434			TPad *pad4 = new TPad("pad4","pad4", .7,.1,.88,.30);
435			pad4->SetFillColor(10);
436
437			TLatex *l3 = new TLatex();
438			l3->SetTextAlign(12);
439			l3->SetTextSize(.3);
440			l3->SetTextColor(2);
441			buffer.str("");
442			buffer<<"DT : "<<DeltaT<<" min";
443
444			canvas2->cd(1);
445			Allev->SetStats(kFALSE);
446			Allev->SetMarkerColor(2);
447			Allev->SetMarkerSize(.5);
448			Allev->SetMarkerStyle(21);
449			Allev->GetXaxis()->SetTitle("Time ( OBT in ms )");
450			Allev->GetXaxis()->CenterTitle();
451			Allev->GetYaxis()->SetTitle("Mean signal ( ADC )");
452			Allev->GetYaxis()->CenterTitle();
453			Allev->SetMinimum(20);
454			Allev->Divide(Alltime);
455			Allev->Draw("p");
456			pad3->Draw();
457			pad3->cd();
458			l3->DrawLatex(0.05, 0.65, buffer.str().c_str());
459
460
461			canvas2->cd(2);
462			if(p != 0){
463			S4ev->SetStats(kFALSE);
464			S4ev->SetMarkerColor(2);
465			S4ev->SetMarkerSize(.5);
466			S4ev->SetMarkerStyle(21);
467			S4ev->GetXaxis()->SetTitle("Time ( OBT in ms )");
468			S4ev->GetXaxis()->CenterTitle();
469			S4ev->GetYaxis()->SetTitle("Mean signal ( ADC )");
470			S4ev->GetYaxis()->CenterTitle();
471			S4ev->SetMinimum(20);
472			S4ev->Divide(S4time);
473			S4ev->Draw("p");
474			pad4->Draw();
475			pad4->cd();
476			l3->DrawLatex(0.05, 0.65, buffer.str().c_str());
477			}else{
478			TLatex *l = new TLatex();
479			l->SetTextAlign(12);
480			l->SetTextSize(0.07);
481			l->SetTextColor(2);
482			noentries.str("");
483			noentries<< "S4_QL: Time evolution";
484			l->DrawLatex(0.06, 0.6, noentries.str().c_str());
485			noentries.str("");
486			noentries<< "No S4 triggered events.";
487			l->SetTextColor(1);
488			l->DrawLatex(0.06, 0.5, noentries.str().c_str());
489			}
490
491			if (outDir == "./") {
492			oss2.str("");
493			oss2 << filename.Data() << "_S4_QL_2." << format.Data();
494			} else {
495			oss2.str("");
496			oss2 << outDir.Data() << filename.Data() << "_S4_QL_2." << format.Data();
497			}
498			canvas2->SaveAs(oss2.str().c_str());
499
500			//------------------- Third Canvas --------------------------------//
501			TCanvas *canvas3 = new TCanvas("S4 - Rate", "S4 - Rate", 1280, 1024);
502			canvas3->SetFillColor(10);
503			canvas3->Divide(1,2);
504
505			canvas3->cd(1);
506			rate->SetStats(kFALSE);
507			rate->SetMarkerColor(2);
508			rate->SetMarkerSize(.5);
509			rate->SetMarkerStyle(21);
510			rate->GetXaxis()->SetTitle("Time ( OBT in ms )");
511			rate->GetXaxis()->CenterTitle();
512			rate->GetYaxis()->SetTitle("S4 rate (Hz)");
513			rate->GetYaxis()->CenterTitle();
514			if(rate->GetMaximum() > 1000) gPad->SetLogy();
515			rate->Draw("9p");
516
517			canvas3->cd(2);
518			rateline->SetStats(kFALSE);
519			rateline->SetMarkerColor(4);
520			rateline->SetMarkerSize(.5);
521			rateline->SetMarkerStyle(21);
522			rateline->GetXaxis()->SetTitle("Time ( OBT in ms )");
523			rateline->GetXaxis()->CenterTitle();
524			rateline->GetYaxis()->SetTitle("S4 rate (Hz)");
525			rateline->GetYaxis()->CenterTitle();
526			rateline->SetMaximum(rate->GetMaximum());
527			if(rateline->GetMaximum() > 1000) gPad->SetLogy();
528			rateline->Draw("9p");
529
530
531			if (outDir == "./") {
532			oss3.str("");
533			oss3 << filename.Data() << "_S4_QL_3." << format.Data();
534			} else {
535			oss3.str("");
536			oss3 << outDir.Data() << filename.Data() << "_S4_QL_3." << format.Data();
537			}
538			canvas3->SaveAs(oss3.str().c_str());
539			file->Close();
540			}
541
542
543			int main(int argc, char* argv[]){
544			TString path;
545			TString outDir ="./";
546			TString format ="jpg";
547			ULong_t DeltaT =1;
548
549
550			if (argc < 2){
551			printf("You have to insert at least the file to analyze \n");
552			printf("Try '--help' for more information. \n");
553			exit(1);
554			}
555
556			if (!strcmp(argv[1], "--help")){
557			printf( "Usage: ND_QL FILE [OPTION] \n");
558			printf( "\t --help Print this help and exit \n");
559			printf( "\t -outDir[path] Path where to put the output [default ./] \n");
560			printf( "\t -format[jpg\|gif\|ps] Format for output files [default 'jpg'] \n");
561			printf( "\t -DeltaT Time interval to control time evolution of acquired data in minute [default = 1 min] \n");
562			exit(1);
563			}
564
565			path=argv[1];
566
567			for (int i = 2; i < argc; i++){
568
569			if (!strcmp(argv[i], "-outDir")){
570			if (++i >= argc){
571			printf( "-outDir needs arguments. \n");
572			printf( "Try '--help' for more information. \n");
573			exit(1);
574			}
575			else{
576			outDir = argv[i];
577			continue;
578			}
579			}
580
581			if (!strcmp(argv[i], "-format")){
582			if (++i >= argc){
583			printf( "-format needs arguments. \n");
584			printf( "Try '--help' for more information. \n");
585			exit(1);
586			}
587			else{
588			format = argv[i];
589			continue;
590			}
591			}
592
593			if (!strcmp(argv[i], "-DeltaT")){
594			if (++i >= argc){
595			printf( "-DeltaT needs arguments. \n");
596			printf( "Try '--help' for more information. \n");
597			exit(1);
598			}
599			else{
600			DeltaT = atol(argv[i]);
601			continue;
602			}
603			}
604
605
606			}
607
608			S4_QL(argv[1], outDir, format, DeltaT);
609
610			}
611