Universita' di GenovaINFN Sezione di Genova  
AIRWATCH / EUSO Genova

EDetectorHistoPainter - source file

// ESAF : Euso Simulation and Analysis Framework
// $Id: EDetectorHistoPainter.cc,v 1.7 2005/11/04 00:17:00 thea Exp $
//  created May, 27 2004

#include "EDetectorHistoPainter.hh"
#include "EDetector.hh"
#include "EDetPhoton.hh"
#include "EElectronics.hh"
#include "EFee.hh"
#include "EAFee.hh"
#include "EMacroCellHit.hh"

#include "Etypes.hh"
#include "ERunParameters.hh"

#include "Riostream.h"
#include "TCanvas.h"
#include "TLegend.h"
#include "TList.h"
#include "TH1F.h"
#include "TH2F.h"
#include "THStack.h"
#include "TStyle.h"
#include "TVector3.h"
#include "TPolyLine.h"


ClassImp(EDetectorHistoPainter)
   
//_____________________________________________________________________________
 EDetectorHistoPainter::EDetectorHistoPainter() {
    // 
    // Constructor
    // 
    
    fDetector = 0;
    fRunPars = 0;

    Build();
}

//_____________________________________________________________________________
 EDetectorHistoPainter::EDetectorHistoPainter( EDetector *det, ERunParameters* rp ) {
    // 
    // Constructor
    // 
    
    fDetector = 0;
    fRunPars = 0;

    SetDetector(det);
    SetRunPars(rp);
    Build();
}


//_____________________________________________________________________________
 EDetectorHistoPainter::~EDetectorHistoPainter() {
    // 
    // Destuctor
    // 

    fHistos->Delete();
}

//_____________________________________________________________________________
 void EDetectorHistoPainter::Build() {
    // 
    // Build 
    // 

    if ( !fDetector || !fRunPars ) {
        Info("Build()","Either EDetector or EShower object are NULL. Painter made zombie.");
        MakeZombie();
        return;
    }

    fHistos = new TList;

    gStyle->SetPalette(1);

    if ( fDetector ) {
        fGtuLength = fDetector->GetElectronics()->GetGtuLength();
        fFirstGtu = TMath::Floor((fDetector->GetFirstTime()
                    -fDetector->GetElectronics()->GetTimeGtu0() ) / fGtuLength);

        fLastGtu = TMath::Ceil((fDetector->GetLastTime()
                    -fDetector->GetElectronics()->GetTimeGtu0() ) / fGtuLength);

        fNumGtus = (Int_t) (fLastGtu -fFirstGtu );
        fGtuLength = fDetector->GetElectronics()->GetGtuLength();
    }


}
//_____________________________________________________________________________
 void EDetectorHistoPainter::Clear( Option_t *opt ) {
    TString option(opt);

    if (option.Contains("C"))
        fHistos->Delete();
    else 
        fHistos->Clear();
    
}

//_____________________________________________________________________________
 TH1F *EDetectorHistoPainter::GetGtuHisto(const char *name, const char *title){
    TH1F* th = (TH1F*)fHistos->FindObject(name);

    if ( !title ) title = name;

    if (!th) {
        
        th = new TH1F( name, title, fNumGtus, fFirstGtu, fLastGtu);
        th->SetDirectory(0);
        fHistos->Add(th);
    }
 
    return th;
}

//_____________________________________________________________________________
 void EDetectorHistoPainter::DoGtuPhotonsHisto(TVirtualPad *c) {
    // Time distribution of photons at various stage of Euso fDetector
    
    if ( !fDetector ) {
        cout << "No Event defined" << endl;
        return;
    }

    if ( !c && TPad::Pad() ) c = TPad::Pad()->GetCanvas();
    if ( !c ) c = new TCanvas();  
    
//     c->Clear();
//     c->cd();

    TPad *pad1 = new TPad("Pad1","Pad1",0.01,0.51,0.49,0.99);
    TPad *pad2 = new TPad("Pad2","Pad2",0.51,0.51,0.99,0.99);
    TPad *pad3 = new TPad("Pad3","Pad3",0.01,0.01,0.49,0.49);
    TPad *pad4 = new TPad("Pad4","Pad4",0.51,0.01,0.99,0.49);

    pad1->Draw();
    pad2->Draw();
    pad3->Draw();
    pad4->Draw();
    
    // build histograms
    TH1F *th0 = GetGtuHisto("EusoHF_GtuPho0", "Photons vs GTU");
    TH1F *th1 = GetGtuHisto("EusoHF_GtuPho1", "Photons(IFS) vs GTU");
    TH1F *th2 = GetGtuHisto("EusoHF_GtuPho2", "Photons(Signal) vs GTU");

    th0->Reset();
    th1->Reset();
    th2->Reset();

    for( Int_t i=0; i < fDetector->GetNumPhotons(); i++ ) {
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        Float_t time = ( ph.GetTime() - fDetector->GetElectronics()->GetTimeGtu0()) / fGtuLength;
        th0->Fill( time );
        if ( ph.CrossedIFS() ) th1->Fill( time );
        if ( ph.CrossedIFS() && ph.MadeSignal() ) th2->Fill( time );
    }
    pad1->cd();
    // draw histograms with colors and filling
    th0->SetStats(0);
    th0->SetXTitle("GTU number");
    th0->SetYTitle("Number of Photons");

    th0->SetFillColor(9);
    th1->SetFillColor(2);
    th2->SetFillColor(8);

    th0->Draw();
    th1->Draw("same");
    th2->Draw("same");
    
    // write meaning and statistics
    Axis_t xMax = th0->GetXaxis()->GetXmax();
    Axis_t yMax = th0->GetMaximum();
    
    TLegend* leg;
    leg = new TLegend(0.68, 0.68, 0.95, 0.95);
    leg->AddEntry(th0,Form("Ph. on Pupil:%d",(int)th0->GetEntries()), "f");
    leg->AddEntry(th1,Form("Ph. on FS:%d",(int)th1->GetEntries()), "f");
    leg->AddEntry(th2,Form("Ph. detected:%d",(int)th2->GetEntries()), "f");
    leg->SetBorderSize(2);
    leg->Draw();

//    TText *t = new TText();
//    t->SetTextFont(62);
//    t->SetTextAlign(31);
//    t->SetTextColor(9);
//    t->SetTextSize(0.04);
//    sprintf(txt,"Ph. on Pupil:%d",(int)th0->GetEntries());
//    t->DrawText(0.92*xMax,0.94*yMax,txt); 
//    t->SetTextColor(2);
//    sprintf(txt,"Ph. on FS:%d",(int)th1->GetEntries());
//    t->DrawText(0.92*xMax,0.86*yMax,txt); 
//    t->SetTextColor(8);
//    sprintf(txt,"Ph. detected:%d",(int)th2->GetEntries());
//    t->DrawText(0.92*xMax,0.78*yMax,txt);
    
    c->Update();
   
    // second pad 
    pad2->cd();
    TH1F *th3 = GetGtuHisto( "EusoHF_GtuPho3", "Photons vs GTU");
    TH1F *th4 = GetGtuHisto( "EusoHF_GtuPho4", "Photons vs GTU");
    TH1F *th5 = GetGtuHisto( "EusoHF_GtuPho5", "Photons vs GTU");
    TH1F *th6 = GetGtuHisto( "EusoHF_GtuPho6", "Photons vs GTU");

    th3->Reset();
    th4->Reset();
    th5->Reset();
    th6->Reset();

    for( Int_t i=0; i < fDetector->GetNumPhotons(); i++ ) {
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        Float_t time = ( ph.GetTime() - fDetector->GetElectronics()->GetTimeGtu0()) / fGtuLength;
        if ( ph.CrossedIFS() ) {
            th3->Fill( time );
            if ( ph.MadeSignal() ) th4->Fill( time ); 
            if ( ph.MadeCount()  ) th5->Fill( time );
            if ( ph.MadeFastOR() ) th6->Fill( time );
        }
    }
    
    th3->SetStats(0);
    th3->SetXTitle("GTU number");
    th3->SetYTitle("Number of Photons");

    th3->Draw();
    th4->Draw("same");
    th5->Draw("same");
    th6->Draw("same");

    th3->SetFillColor(9);
    th4->SetFillColor(6);
    th5->SetFillColor(2);
    th6->SetFillColor(8);

    xMax = th3->GetXaxis()->GetXmax();
    yMax = th3->GetMaximum();

    leg = new TLegend(0.68, 0.68, 0.95, 0.95);
    leg->AddEntry(th3,Form("Photons on FS:%d",(int)th3->GetEntries()), "f");
    leg->AddEntry(th4,Form("Pmt signals:%d",(int)th4->GetEntries()), "f");
    leg->AddEntry(th5,Form("Frontend hits:%d",(int)th5->GetEntries()), "f");
    leg->AddEntry(th6,Form("Fast OR counts:%d",(int)th6->GetEntries()), "f");
    leg->SetBorderSize(2);
    leg->Draw();
//
//    
//    TText *t2 = new TText();
//    t2->SetTextFont(62);
//    t2->SetTextAlign(31);
//    t2->SetTextColor(9);
//    t2->SetTextSize(0.04);
//    sprintf(txt,"Photons on FS:%d",(int)th3->GetEntries());
//    t2->DrawText(0.92*xMax,0.94*yMax,txt); 
//    t2->SetTextColor(6);
//    sprintf(txt,"Pmt signals:%d",(int)th4->GetEntries());
//    t2->DrawText(0.92*xMax,0.86*yMax,txt); 
//    t2->SetTextColor(2);
//    sprintf(txt,"Frontend hits:%d",(int)th5->GetEntries());
//    t2->DrawText(0.92*xMax,0.78*yMax,txt);
//    t2->SetTextColor(8);
//    sprintf(txt,"Fast OR counts:%d",(int)th6->GetEntries());
//    t2->DrawText(0.92*xMax,0.70*yMax,txt);
     
    c->Update();

    // third pad: charge
    pad3->cd();

    TH1F *th7 = GetGtuHisto("EusoHF_GtuPho7", "Charge vs GTU");
    TH1F *th8 = GetGtuHisto("EusoHF_GtuPho9", "Photons vs GTU");
    TH1F *th9 = GetGtuHisto("EusoHF_GtuPho8", "Photons(IFS) vs GTU");
    
    Float_t gain;
    gain = (fRunPars->GetPmtData()).GetGain();
    for( Int_t i=0; i < fDetector->GetNumPhotons(); i++ ) {
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        Float_t time = ( ph.GetTime() - fDetector->GetElectronics()->GetTimeGtu0() ) / fGtuLength;
        if ( ph.CrossedIFS() && ph.MadeCount() ) th8->Fill( time );
        if ( ph.CrossedIFS() && ph.MadeFastOR() ) th9->Fill( time );
    }
    Float_t total_afee_charge=0.;
    for(Int_t i=0; i < fDetector->GetNumAFee(); i++ ) {
        EAFee *afee = fDetector->GetAFee( i );
        Float_t gtu = (Float_t)afee->GetGtu();
        Float_t charge = afee->GetDyCharge();
        Int_t numpe = (Int_t)(charge/1.6e-19/gain);   // convert into NOMINAL photoelectrons, gain 3.E5
        if ( gtu > fFirstGtu && gtu < fLastGtu ) {
            Int_t gtubin = th7->FindBin( gtu );
            Int_t old = (Int_t)(th7->GetBinContent( gtubin ));
            th7->SetBinContent( gtubin, old + numpe );
            total_afee_charge += numpe;
        }
    }

    th7->SetStats(0);
    th7->SetXTitle("GTU number");
    th7->SetYTitle("Number of Photons");

    th7->SetFillColor(8);
    th8->SetFillColor(9);
    th9->SetFillColor(6);
    
    th7->Draw();
    th8->Draw("same");
    th9->Draw("same");


    xMax = th7->GetXaxis()->GetXmax();
    yMax = th7->GetMaximum();

    leg = new TLegend(0.68, 0.68, 0.95, 0.95);
    leg->AddEntry(th7,Form("Ph. detected:%d",(int)th7->GetEntries()), "f");
    leg->AddEntry(th8,Form("Ph. in DFEE:%d",(int)th8->GetEntries()), "f");
    leg->AddEntry(th9,Form("Charge (pC,dynode):%d",(int)th9->GetEntries()), "f");
    leg->SetBorderSize(2);
    leg->Draw();

//    TText *t3 = new TText();
//    t3->SetTextFont(62);
//    t3->SetTextAlign(31);
//    t3->SetTextColor(9);
//    t3->SetTextSize(0.04);
//    sprintf(txt,"Ph. detected:%d",(Int_t)th8->GetEntries());
//    t3->DrawText(0.92*xMax,0.94*yMax,txt); 
//    t3->SetTextColor(6);
//    sprintf(txt,"Ph. in DFEE:%d",(Int_t)th9->GetEntries());
//    t3->DrawText(0.92*xMax,0.86*yMax,txt); 
//    t3->SetTextColor(8);
//    sprintf(txt,"Charge (pC,dynode):%.2f",total_afee_charge);
//    t3->DrawText(0.92*xMax,0.78*yMax,txt);
    

    // fourth pad
    pad4->cd();
    TH1F *th10 = GetGtuHisto("EDFrontEndHits","FrontEnd Hits");
    TH1F *th11 = GetGtuHisto("EDFrontEndSignals","FrontEnd Signal");
    TH1F *th12 = GetGtuHisto("EDFrontEndAnalogHits","FrontEnd Signal");

    for( Int_t i(0); i<fDetector->GetNumFee(); i++ ) {
        EFee* fe = fDetector->GetFee( i );
        Int_t nSignals = fe->GetNumSignals();
        Int_t nHits = fe->GetNumHits();

        if ( nSignals > 0 )
            while ( nHits-- > 0 ) 
                th10->Fill((Float_t)fe->GetGtu());

        while ( nSignals-- > 0 ) 
            th11->Fill((Float_t)fe->GetGtu());
    }
    
    for(Int_t i=0; i < fDetector->GetNumAFee(); i++ ) {
        EAFee *afee = fDetector->GetAFee( i );
        Int_t nPe = (Int_t)(afee->GetDyCharge()/1.6e-19/gain);   // convert into NOMINAL photoelectrons, gain 3.E5

        while ( nPe-- > 0 ) 
            th12->Fill((Float_t)afee->GetGtu());
    }

    th10->SetStats(0);
    th10->SetXTitle("GTU number");
    th10->SetYTitle("Hits");
    
    th10->SetFillColor(216);
    th11->SetFillColor(220);
    th12->SetLineColor(kRed);

    th10->Draw();
    th11->Draw("same");
    th12->Draw("same");

    xMax = th10->GetXaxis()->GetXmax();
    yMax = th10->GetMaximum();

    leg = new TLegend(0.68, 0.68, 0.95, 0.95);
    leg->AddEntry(th10,Form("DFee Hits:%d",(int)th10->GetEntries()), "f");
    leg->AddEntry(th11,Form("DFee Signal:%d",(int)th11->GetEntries()), "f");
    leg->AddEntry(th12,Form("AFee Hits:%d",(int)th12->GetEntries()), "f");
    leg->SetBorderSize(2);
    leg->Draw();

//    TText *t4 = new TText();
//    t4->SetTextFont(62);
//    t4->SetTextAlign(31);
//    t4->SetTextColor(216);
//    t4->SetTextSize(0.04);
//    sprintf(txt,"DFee Hits:%d",(Int_t)th10->GetEntries());
//    t4->DrawText(0.92*xMax,0.94*yMax,txt); 
//    t4->SetTextColor(220);
//    sprintf(txt,"DFee Pure Signal Hits:%d",(Int_t)th11->GetEntries());
//    t4->DrawText(0.92*xMax,0.86*yMax,txt); 
//    t4->SetTextColor(kRed);
//    sprintf(txt,"AFee Hits:%d",(Int_t)th12->GetEntries());
//    t4->DrawText(0.92*xMax,0.78*yMax,txt);
//    c->Update();
}

//_____________________________________________________________________________
 void EDetectorHistoPainter::DoGtuPhotonsHisto2(TVirtualPad *c) {
    // Time distribution of photons at various stage of Euso detector
    
    if ( !fDetector ) {
        cout << "No Event defined" << endl;
        return;
    }

    if ( !c && TPad::Pad() ) c = TPad::Pad()->GetCanvas();
    if ( !c ) c = new TCanvas();  
    
//     c->Clear();
//     c->cd();

    TPad *pad1 = new TPad("Pad1","Pad1",0.01,0.51,0.49,0.99);
    TPad *pad2 = new TPad("Pad2","Pad2",0.51,0.51,0.99,0.99);
    TPad *pad3 = new TPad("Pad3","Pad3",0.01,0.01,0.49,0.49);
    TPad *pad4 = new TPad("Pad4","Pad4",0.51,0.01,0.99,0.49);

    pad1->Draw();
    pad2->Draw();
    pad3->Draw();
    pad4->Draw();
    
    pad1->cd();

    // build histograms
    //
    //look for hs in fHistos list.
    
    THStack *hs1 = (THStack*)fHistos->FindObject("EStackGtuPho1");
    if ( !hs1 ) {
       hs1 = new THStack("EStackGtuPho1","Photon's types");
//       fHistos->Add( hs1 );
    }
    
    TH1F *th0 = GetGtuHisto("EusoHF_GtuPho10", "#gamma vs GTU");
    TH1F *th1 = GetGtuHisto("EusoHF_GtuPho11", "Fluorescence #gamma vs GTU");
    TH1F *th2 = GetGtuHisto("EusoHF_GtuPho12", "Direct Cherenkov #gamma vs GTU");
    TH1F *th3 = GetGtuHisto("EusoHF_GtuPho13", "BackScattered Cherenkov #gamma  vs GTU");

    th0->Reset();
    th1->Reset();
    th2->Reset();
    th3->Reset();

    for( Int_t i=0; i < fDetector->GetNumPhotons(); i++ ) {
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        
        Float_t time = ( ph.GetTime() - fDetector->GetFirstTime() ) / fGtuLength;
        th0->Fill(time);
        if ( ph.GetType() == kFluorescence ) th1->Fill(time);
        if ( ph.GetType() == kDirectCherenkov ) th2->Fill(time);
        if ( ph.GetType() == kBackCherenkov ) th3->Fill(time);
    }

    th1->SetFillColor(kBlue);
    th2->SetFillColor(kRed);
    th3->SetFillColor(kGreen);

    //hs->Add(th0);
    hs1->Add(th1);
    hs1->Add(th2);
    hs1->Add(th3);
//    th0->Draw();
//    th1->Draw("same");
//    th2->Draw("same");
//    th3->Draw("same");
    hs1->Draw();

    Float_t xMax, yMax;
    
    xMax = hs1->GetXaxis()->GetXmax();
    yMax = hs1->GetMaximum();
    
    TLegend* leg;
    leg = new TLegend(0.68, 0.78, 0.95, 0.95);
    leg->AddEntry(th1,Form("Fluo : %d",(int)th1->GetEntries()), "f");
    leg->AddEntry(th2,Form("Dir. Cher. : %d",(int)th2->GetEntries()), "f");
    leg->AddEntry(th3,Form("Bck. Cher. : %d",(int)th3->GetEntries()), "f");
    leg->SetBorderSize(2);
    leg->Draw();

    TString txt;
//    TText *t1 = new TText();
//    t1->SetTextFont(62);
//    t1->SetTextAlign(31);
//    t1->SetTextSize(0.04);
//
//    t1->SetTextColor(kBlue);
//    txt = Form("Fluorescence:%d",(Int_t)th1->GetEntries());
//    t1->DrawText(0.92*xMax,0.94*yMax,txt); 
//
//    t1->SetTextColor(kRed);
//    txt = Form("Direct Cherenkov:%d",(Int_t)th2->GetEntries());
//    t1->DrawText(0.92*xMax,0.86*yMax,txt); 
//    
//    t1->SetTextColor(kGreen);
//    txt = Form("BackScattered Cherenkov:%d",(Int_t)th3->GetEntries());
//    t1->DrawText(0.92*xMax,0.78*yMax,txt);

    //
    // **************************
    //

    pad2->cd();

    THStack *hs2 = (THStack*)fHistos->FindObject("EStackGtuPho2");
    if ( !hs2 ) {
       hs2 = new THStack("EStackGtuPho2","Photons in the Detector");
  //     fHistos->Add( hs2 );
    }
    TH1F *th4 = GetGtuHisto("EusoHF_GtuPho14", "#gamma in");
    TH1F *th5 = GetGtuHisto("EusoHF_GtuPho15", "#gamma in");
    TH1F *th6 = GetGtuHisto("EusoHF_GtuPho16", "#gamma in");
    TH1F *th7 = GetGtuHisto("EusoHF_GtuPho17", "#gamma in");
    TH1F *th8 = GetGtuHisto("EusoHF_GtuPho18", "#gamma in");
    TH1F *th9 = GetGtuHisto("EusoHF_GtuPho19", "#gamma in");

    th4->Reset();
    th5->Reset();
    th6->Reset();
    th7->Reset();
    th8->Reset();
    th9->Reset();

    for( Int_t i=0; i < fDetector->GetNumPhotons(); i++ ) {
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        
        Float_t time = ( ph.GetTime() - fDetector->GetFirstTime() ) / fGtuLength;
        if ( ph.GetHistory() == kNotInside ) th4->Fill(time);
        if ( ph.GetHistory() == kBaffle ) th5->Fill(time);
        if ( ph.GetHistory() == kOptics ) th6->Fill(time);
        if ( ph.GetHistory() == kFocalPlane ) th7->Fill(time);
        if ( ph.GetHistory() == kWalls ) th8->Fill(time);
        if ( ph.GetHistory() == kEarth ) th9->Fill(time);
    }
    
    th4->SetFillColor(kBlue);
    th5->SetFillColor(kRed);
    th6->SetFillColor(kGreen);
    th7->SetFillColor(kCyan);
    th8->SetFillColor(kMagenta);
    th9->SetFillColor(kYellow);
    
    hs2->Add(th4);
    hs2->Add(th5);
    hs2->Add(th6);
    hs2->Add(th7);
    hs2->Add(th8);
    hs2->Add(th9);
    hs2->Draw();

    xMax = hs2->GetXaxis()->GetXmax();
    yMax = hs2->GetMaximum();

    leg = new TLegend(0.68, 0.58, 0.95, 0.95);
    leg->AddEntry(th4,Form("Not Inside : %d",(int)th4->GetEntries()), "f");
    leg->AddEntry(th5,Form("Baffle : %d",(int)th5->GetEntries()), "f");
    leg->AddEntry(th6,Form("Optics : %d",(int)th6->GetEntries()), "f");
    leg->AddEntry(th7,Form("FocalPlane : %d",(int)th7->GetEntries()), "f");
    leg->AddEntry(th8,Form("Walls : %d",(int)th8->GetEntries()), "f");
    leg->AddEntry(th9,Form("Reflected : %d",(int)th9->GetEntries()), "f");
    leg->SetBorderSize(2);
    leg->Draw();

//    TText *t2 = new TText();
//    t2->SetTextFont(62);
//    t2->SetTextAlign(31);
//    t2->SetTextSize(0.04);
//
//    t2->SetTextColor(kBlue);
//    txt = Form("Not Inside :%d",(Int_t)th4->GetEntries());
//    t2->DrawText(0.92*xMax,0.94*yMax,txt); 
//
//    t2->SetTextColor(kRed);
//    txt = Form("Baffle:%d",(Int_t)th5->GetEntries());
//    t2->DrawText(0.92*xMax,0.86*yMax,txt); 
//    
//    t2->SetTextColor(kGreen);
//    txt = Form("Optics:%d",(Int_t)th6->GetEntries());
//    t2->DrawText(0.92*xMax,0.78*yMax,txt);
//    
//    t2->SetTextColor(kCyan);
//    txt = Form("FocalPlane:%d",(Int_t)th6->GetEntries());
//    t2->DrawText(0.92*xMax,0.70*yMax,txt);
//
//    t2->SetTextColor(kMagenta);
//    txt = Form("Walls:%d",(Int_t)th6->GetEntries());
//    t2->DrawText(0.92*xMax,0.62*yMax,txt);
//
//    t2->SetTextColor(kYellow);
//    txt = Form("Reflected toward Earth:%d",(Int_t)th6->GetEntries());
//    t2->DrawText(0.92*xMax,0.54*yMax,txt);

    //
    // **************************
    //

    pad3->cd();

    THStack *hs3 = (THStack*)fHistos->FindObject("EStackGtuPho3");
    if ( !hs3 ) {
       hs3 = new THStack("EStackGtuPho2","Photons vs Wavelength");
  //     fHistos->Add( hs2 );
    }
    TH1F *th10 = GetGtuHisto("EusoHF_GtuPho20", "#gamma in");
    TH1F *th11 = GetGtuHisto("EusoHF_GtuPho21", "#gamma in");
    TH1F *th12 = GetGtuHisto("EusoHF_GtuPho22", "#gamma in");

    th10->Reset();
    th11->Reset();
    th12->Reset();

    for( Int_t i=0; i < fDetector->GetNumPhotons(); i++ ) {
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        
        Float_t time = ( ph.GetTime() - fDetector->GetFirstTime() ) / fGtuLength;
        if ( ph.GetLambda() < 340E-6 ) th10->Fill(time);
        else if ( ph.GetLambda() > 370E-6 ) th12->Fill(time);
        else th11->Fill(time);
    }
    
    th10->SetFillColor(223);
    th11->SetFillColor(211);
    th12->SetFillColor(94);
    
    hs3->Add(th10);
    hs3->Add(th11);
    hs3->Add(th12);
    hs3->Draw();

    xMax = hs3->GetXaxis()->GetXmax();
    yMax = hs3->GetMaximum();

    leg = new TLegend(0.54, 0.78, 0.95, 0.95);
    leg->AddEntry(th10,Form("#lambda < 340 nm : %d",(int)th10->GetEntries()), "f");
    leg->AddEntry(th11,Form("340 nm < #lambda < 370 nm : %d",(int)th11->GetEntries()), "f");
    leg->AddEntry(th12,Form("#lambda > 370 nm : %d",(int)th12->GetEntries()), "f");
    leg->SetBorderSize(2);
    leg->Draw();
//
//    TText *t3 = new TText();
//    t3->SetTextFont(62);
//    t3->SetTextAlign(31);
//    t3->SetTextSize(0.04);
//
//    t3->SetTextColor(223);
//    txt = Form("#lambda < 340 nm: %d",(Int_t)th6->GetEntries());
//    t3->DrawText(0.92*xMax,0.94*yMax,txt);
//
//    t3->SetTextColor(211);
//    txt = Form("340 nm < #lambda < 370 nm: %d",(Int_t)th6->GetEntries());
//    t3->DrawText(0.92*xMax,0.86*yMax,txt);
//
//    t3->SetTextColor(94);
//    txt = Form("#lambda > 370 nm: %d",(Int_t)th6->GetEntries());
//    t3->DrawText(0.92*xMax,0.78*yMax,txt);
//

    pad4->cd();
    // Aggiungere lo spettro dei fotoni
    pad4->SetLogy();
    

    Float_t lmin = kHuge;
    Float_t lmax = -kHuge;
    
    for( Int_t i=0; i < fDetector->GetNumPhotons(); i++ ) {
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        
        if ( ph.GetLambda() > lmax ) lmax = ph.GetLambda();
        if ( ph.GetLambda() < lmin ) lmin = ph.GetLambda();
    }

    Bool_t status = TH1::AddDirectoryStatus();
    TH1::AddDirectory(kFALSE);

    THStack *hs4 = (THStack*)fHistos->FindObject("EStackWlPho1");
    if ( !hs4 ) {
       hs4 = new THStack("EStackWlPho1","Photons` spectrum");
    }

    TH1F* th13 = (TH1F*)fHistos->FindObject("EusoHF_Wl1");
    if (!th13) {
        th13 = new TH1F( "EusoHF_Wl1", "Photons' spectrum", 100, lmin, lmax);
        fHistos->Add(th13);
    }
    TH1F* th14 = (TH1F*)fHistos->FindObject("EusoHF_Wl2");
    if (!th14) {
        th14 = new TH1F( "EusoHF_Wl2", "Photons' spectrum", 100, lmin, lmax);
        fHistos->Add(th14);
    }
    TH1F* th15 = (TH1F*)fHistos->FindObject("EusoHF_Wl3");
    if (!th15) {
        th15 = new TH1F( "EusoHF_Wl3", "Photons' spectrum", 100, lmin, lmax);
        fHistos->Add(th15);
    }

    th13->Reset();
    TH1::AddDirectory(status);

    for( Int_t i=0; i < fDetector->GetNumPhotons(); i++ ) {
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        if ( ph.GetType() == kFluorescence ) th13->Fill( ph.GetLambda() );
        if ( ph.GetType() == kDirectCherenkov ) th14->Fill( ph.GetLambda() );
        if ( ph.GetType() == kBackCherenkov ) th15->Fill( ph.GetLambda() );
    }
    
    th13->SetFillColor(kBlue);
    th14->SetFillColor(kRed);
    th15->SetFillColor(kGreen);

    hs4->Add(th14);
    hs4->Add(th15);
    hs4->Add(th13);

    hs4->Draw();
    xMax = hs4->GetXaxis()->GetXmax();
    yMax = hs4->GetMaximum();

    leg = new TLegend(0.68, 0.78, 0.95, 0.95);
    leg->AddEntry(th13,Form("Fluo. : %d",(int)th13->GetEntries()), "f");
    leg->AddEntry(th14,Form("Dir. Cher. : %d",(int)th14->GetEntries()), "f");
    leg->AddEntry(th15,Form("Bak. Cher. : %d",(int)th15->GetEntries()), "f");
    leg->SetBorderSize(2);
    leg->Draw();

//    TText *t4 = new TText();
//    t4->SetTextFont(62);
//    t4->SetTextAlign(31);
//    t4->SetTextSize(0.04);
//
//    t4->SetTextColor(kBlue);
//    if (!pad4->GetLogy()) {
//        txt = Form("Fluorescence:%d",(Int_t)th13->GetEntries());
//        t4->DrawText(0.92*xMax,0.94*yMax,txt); 
//
//        t4->SetTextColor(kRed);
//        txt = Form("Direct Cherenkov:%d",(Int_t)th14->GetEntries());
//        t4->DrawText(0.92*xMax,0.86*yMax,txt); 
//    
//        t4->SetTextColor(kGreen);
//        txt = Form("BackScattered Cherenkov:%d",(Int_t)th15->GetEntries());
//        t4->DrawText(0.92*xMax,0.78*yMax,txt);
//    }
//    else {
//       // redefine a bit the layout
//       hs4->SetMaximum(yMax*5);
//       hs4->SetMinimum(yMax*1.e-4);
//       txt = Form("Fluorescence:%d",(Int_t)th13->GetEntries());
//       t4->DrawText(0.92*xMax,5*yMax,txt); 
//
//       t4->SetTextColor(kRed);
//       txt = Form("Direct Cherenkov:%d",(Int_t)th14->GetEntries());
//       t4->DrawText(0.92*xMax,2.5*yMax,txt); 
//    
//       t4->SetTextColor(kGreen);
//       txt = Form("BackScattered Cherenkov:%d",(Int_t)th15->GetEntries());
//       t4->DrawText(0.92*xMax,1*yMax,txt);
//   }
  
}

//______________________________________________________________________________
 void EDetectorHistoPainter::DoPhiCosThetaHisto(TVirtualPad *c){
    // photons on focal surface

    if ( !fDetector ) {
        cout << "No Event defined" << endl;
        return;
    }

    if ( !c && TPad::Pad() ) c = TPad::Pad()->GetCanvas();
    if ( !c ) c = new TCanvas();  

    c->Clear();
    c->Divide(2,2);

    Float_t th0CosThMin(10),th0CosThMax(-10);
    Float_t th0PhiMin(10),th0PhiMax(-10);

    for ( Int_t i=0; i < fDetector->GetNumPhotons(); i++){
        EDetPhoton& ph=*(fDetector->GetPhoton(i));
        if (TMath::Cos(ph.GetTheta()) > th0CosThMax)
            th0CosThMax = TMath::Cos(ph.GetTheta());

        if (ph.GetPhi() > th0PhiMax)
            th0PhiMax = ph.GetPhi();

        if (TMath::Cos(ph.GetTheta()) < th0CosThMin)
            th0CosThMin = TMath::Cos(ph.GetTheta());

        if (ph.GetPhi() < th0PhiMin)
            th0PhiMin = ph.GetPhi();

    }
    Float_t th0Delta = th0CosThMax-th0CosThMin;
    th0CosThMin-=0.1*th0Delta;
    th0CosThMax+=0.1*th0Delta;

    th0Delta = th0PhiMax-th0PhiMin;
    th0PhiMin-=0.1*th0Delta;
    th0PhiMax+=0.1*th0Delta;

    gDirectory->Delete("EusoHF_PhiCosTh*");

    Bool_t status = TH1::AddDirectoryStatus();
    TH1::AddDirectory(kFALSE);
    TH2F *th0 = new TH2F( "EusoHF_PhiCosTh0", "Cos(Theta) vs Phi (CrossedIFS)", 100, 
            th0CosThMin, th0CosThMax, 100, th0PhiMin, th0PhiMax);
    TH2F *th1 = new TH2F( "EusoHF_PhiCosTh1", "Cos(Theta) vs Phi (MadeSignal)", 100, 
            th0CosThMin, th0CosThMax, 100, th0PhiMin, th0PhiMax);
    TH2F *th2 = new TH2F( "EusoHF_PhiCosTh2", "Cos(Theta) vs Phi (MadeCount)", 100, 
            th0CosThMin, th0CosThMax, 100, th0PhiMin, th0PhiMax);
    TH2F *th3 = new TH2F( "EusoHF_PhiCosTh3", "Cos(Theta) vs Phi (MadeFastOR)", 100, 
            th0CosThMin, th0CosThMax, 100, th0PhiMin, th0PhiMax);
    TH1::AddDirectory(status);

    for ( Int_t i=0; i < fDetector->GetNumPhotons(); i++){
        EDetPhoton& ph = *(fDetector->GetPhoton( i ));
        if ( ph.CrossedIFS() ) th0->Fill(TMath::Cos(ph.GetTheta()),ph.GetPhi());
        if ( ph.MadeSignal() ) th1->Fill(TMath::Cos(ph.GetTheta()),ph.GetPhi());
        if ( ph.MadeCount() ) th2->Fill(TMath::Cos(ph.GetTheta()),ph.GetPhi());
        if ( ph.MadeFastOR() ) th3->Fill(TMath::Cos(ph.GetTheta()),ph.GetPhi());
    }

    th0->SetStats(0);
    th1->SetStats(0);
    th2->SetStats(0);
    th3->SetStats(0);

    c->cd(1);
    th0->Draw("colz");
    c->cd(3);
    th1->Draw("colz");
    c->cd(2);
    th2->Draw("colz");
    c->cd(4);
    th3->Draw("colz");

    c->Update();
}

//_____________________________________________________________________________
 void EDetectorHistoPainter::DoXYHisto(TVirtualPad *c){
    //photons on XY plane

    if ( !fDetector ) {
        cout << "No Event defined" << endl;
        return;
    }

    if ( !c && TPad::Pad() ) c = TPad::Pad()->GetCanvas();
    if ( !c ) c = new TCanvas();  

    c->Clear();
    c->Divide(2,2);


    //Get the number of photons
    Int_t numph=fDetector->GetNumPhotons();
    cout << "Total number of photons: " << numph << endl;

    gDirectory->Delete("EusoHF_XY*");


    Float_t *x = new Float_t[numph];
    Float_t *y = new Float_t[numph];
    Float_t *xsh = new Float_t[numph];
    Float_t *ysh = new Float_t[numph];
    Float_t xshmax=0, xshmin=0, yshmax=0, yshmin=0;

    Int_t *nlist = new Int_t[numph];

    //read coordinates x,y for each photon
    for (Int_t i=0; i<numph; i++) {
        const TVector3& pos = fDetector->GetPhoton(i)->GetPos();
        x[i]=pos[0];
        y[i]=pos[1];
    }

    Int_t j=0;

    //roughly determinate which photons belong to the shower
    for (Int_t i=0; i<numph; i++){
        Float_t deltax=TMath::Abs(x[i]-x[i+1]);

        if (deltax < 1.5) {
            Float_t deltay=TMath::Abs(y[i]-y[i+1]);

            if (deltay < 1.5) {
                xsh[j]=x[i];
                ysh[j]=y[i];
                xsh[j+1]=x[i+1];
                ysh[j+1]=y[i+1];
                j++;
            }
        }
    }

    j++; 

    //calculate x,y of the shower minimum and maximum
    for (Int_t i=0; i<j; i++) {
        if (i==0) {
            xshmin=xsh[0];
            xshmax=xsh[0];
            yshmin=ysh[0];
            yshmax=ysh[0];
        }
        if (i>0) {
            if (xsh[i]<xshmin) {xshmin=xsh[i];}
            if (xsh[i]>xshmax) {xshmax=xsh[i];}
            if (ysh[i]<yshmin) {yshmin=ysh[i];}
            if (ysh[i]>yshmax) {yshmax=ysh[i];}
        }
    }

    Float_t deltashx = xshmax - xshmin;
    Float_t deltashy = yshmax - yshmin;
    xshmin -= deltashx / 10.0;
    xshmax += deltashx / 10.0;
    yshmin -= deltashy / 10.0;
    yshmax += deltashy / 10.0;

    cout << "Creating istograms..." << endl;
    cout << "X-Min: " << xshmin << endl;
    cout << "X-Max: " << xshmax << endl;
    cout << "Y-Min: " << yshmin << endl;
    cout << "Y-Max: " << yshmax << endl;

    //Select photons in the interval (xshmin,xshmax);(yshmin,yshmax)
    Int_t l=0;
    for (Int_t i=0; i<numph; i++) {
        if ((((x[i]>xshmin && x[i]<xshmax) && y[i]>yshmin) && y[i]<yshmax)) {
            nlist[l]=i;
            l++;
        }
    }
    l++;

    //Create istograms and fill them

    Bool_t status = TH1::AddDirectoryStatus();
    TH1::AddDirectory(kFALSE);

    TH2F *th1=new TH2F("EusoHF_XY1","Hits on Focal Surface",100,xshmin,xshmax,100,yshmin,yshmax);
    TH2F *th2=new TH2F("EusoHF_XY2","MadeSignal",100,xshmin,xshmax,100,yshmin,yshmax);
    TH2F *th3=new TH2F("EusoHF_XY3","MadeCount",100,xshmin,xshmax,100,yshmin,yshmax);
    TH2F *th4=new TH2F("EusoHF_XY4","MadeFastOR",100,xshmin,xshmax,100,yshmin,yshmax); 

    TH1::AddDirectory(status);

    for (Int_t i=0; i<l; i++) {
        Int_t k=nlist[i];
        EDetPhoton* ph=fDetector->GetPhoton(k);
        const TVector3& pos=(fDetector->GetPhoton(k)->GetPos());

        if (ph->CrossedIFS()) th1->Fill(pos[0],pos[1]);
        if (ph->MadeSignal()) th2->Fill(pos[0],pos[1]);
        if (ph->MadeCount()) th3->Fill(pos[0],pos[1]);
        if (ph->MadeFastOR()) th4->Fill(pos[0],pos[1]);
    }

    gStyle->SetOptStat(10);

    th1->SetStats(kTRUE);
    th1->SetXTitle("x");
    th1->SetYTitle("y");
    th2->SetStats(kTRUE);
    th2->SetXTitle("x");
    th2->SetYTitle("y");
    th3->SetStats(kTRUE);
    th3->SetXTitle("x");
    th3->SetYTitle("y");
    th4->SetStats(kTRUE);
    th4->SetXTitle("x");
    th4->SetYTitle("y");

    c->cd(1);
    th1->Draw("colz"); 
    c->cd(2);   
    th2->Draw("colz");
    c->cd(3);
    th3->Draw("colz");
    c->cd(4);
    th4->Draw("colz");

    c->Update();

    if ( x ) delete x;
    if ( y ) delete y;
    if ( xsh ) delete xsh;
    if ( ysh ) delete ysh;
    if ( nlist ) delete nlist;

}

//______________________________________________________________________________
 void EDetectorHistoPainter::DoWavelengthXYHisto(TVirtualPad *c){
    //
    // photons on the IFS in XY plane at various wavelenghts
    //

    if ( !fDetector ) {
        cout << "No Event defined" << endl;
        return;
    }

    if ( !c && TPad::Pad() ) c = TPad::Pad()->GetCanvas();
    if ( !c ) c = new TCanvas();  

    c->Clear();
    c->Divide(2,2);


    //Get the number of photons
    Int_t numph=fDetector->GetNumPhotons();
    cout << "Total number of photons: " << numph << endl;

    gDirectory->Delete("EusoHF_WLXY*");


    Float_t *x = new Float_t[numph];
    Float_t *y = new Float_t[numph];
    Float_t *xsh = new Float_t[numph];
    Float_t *ysh = new Float_t[numph];
    Float_t xshmax=0, xshmin=0, yshmax=0, yshmin=0;

    Int_t *nlist = new Int_t[numph];

    //read coordinates x,y for each photon
    for (Int_t i=0; i<numph; i++) {
        const TVector3 &pos=(fDetector->GetPhoton(i)->GetPos());
        x[i]=pos[0];
        y[i]=pos[1];
    }

    Int_t j=0;

    //determinate which photons belong to the shower
    for (Int_t i=0; i<numph; i++){
        Float_t deltax=TMath::Abs(x[i]-x[i+1]);

        if (deltax < 1.5) {
            Float_t deltay=TMath::Abs(y[i]-y[i+1]);

            if (deltay < 1.5) {
                xsh[j]=x[i];
                ysh[j]=y[i];
                xsh[j+1]=x[i+1];
                ysh[j+1]=y[i+1];
                j++;
            }
        }
    }

    j++; 

    //calculate x,y of the shower minimum and maximum
    for (Int_t i=0; i<j; i++) {
        if (i==0) {
            xshmin=xsh[0];
            xshmax=xsh[0];
            yshmin=ysh[0];
            yshmax=ysh[0];
        }
        if (i>0) {
            if (xsh[i]<xshmin) {xshmin=xsh[i];}
            if (xsh[i]>xshmax) {xshmax=xsh[i];}
            if (ysh[i]<yshmin) {yshmin=ysh[i];}
            if (ysh[i]>yshmax) {yshmax=ysh[i];}
        }
    }

    Float_t deltashx = xshmax - xshmin;
    Float_t deltashy = yshmax - yshmin;
    xshmin -= deltashx / 10.0;
    xshmax += deltashx / 10.0;
    yshmin -= deltashy / 10.0;
    yshmax += deltashy / 10.0;

    cout << "Creating histograms..." << endl;
    cout << "X-Min: " << xshmin << endl;
    cout << "X-Max: " << xshmax << endl;
    cout << "Y-Min: " << yshmin << endl;
    cout << "Y-Max: " << yshmax << endl;

    //Select photons in the interval (xshmin,xshmax);(yshmin,yshmax)
    Int_t l=0;
    for (Int_t i=0; i<numph; i++) {
        if ((((x[i]>xshmin && x[i]<xshmax) && y[i]>yshmin) && y[i]<yshmax)) {
            nlist[l]=i;
            l++;
        }
    }
    l++;

    //Create istograms and fill them

    TH2F *th1=new TH2F("EusoHF_WLXY1","Hits on IFS, all wavelengths", 100,xshmin,xshmax,100,yshmin,yshmax);
    TH2F *th2=new TH2F("EusoHF_WLXY2","Hits on IFS, #lambda < 350  nm", 100,xshmin,xshmax,100,yshmin,yshmax);
    TH2F *th3=new TH2F("EusoHF_WLXY3","Hits on IFS, 350 nm #leq #lambda #leq 370 nm", 100,xshmin,xshmax,100,yshmin,yshmax);
    TH2F *th4=new TH2F("EusoHF_WLXY4","Hits on IFS, #lambda > 370 nm", 100,xshmin,xshmax,100,yshmin,yshmax); 

    for (Int_t i=0; i<l; i++) {
        Int_t k=nlist[i];
        EDetPhoton* ph=fDetector->GetPhoton(k);
        const TVector3& pos=(fDetector->GetPhoton(k)->GetPos());

        if (ph->CrossedIFS()) th1->Fill(pos[0],pos[1]);
        if (ph->CrossedIFS() && ph->GetLambda()<0.00035) th2->Fill(pos[0],pos[1]);
        if (ph->CrossedIFS() && ph->GetLambda()>=0.00035 && ph->GetLambda()<=0.00037) th3->Fill(pos[0],pos[1]);
        if (ph->CrossedIFS() && ph->GetLambda()>0.00037) th4->Fill(pos[0],pos[1]);
    }

    gStyle->SetOptStat(10);

    th1->SetStats(kTRUE);
    th1->SetXTitle("x");
    th1->SetYTitle("y");
    th2->SetStats(kTRUE);
    th2->SetXTitle("x");
    th2->SetYTitle("y");
    th3->SetStats(kTRUE);
    th3->SetXTitle("x");
    th3->SetYTitle("y");
    th4->SetStats(kTRUE);
    th4->SetXTitle("x");
    th4->SetYTitle("y");


    c->cd(1);
    th1->Draw("colz"); 
    c->cd(2);   
    th2->Draw("colz");
    c->cd(3);
    th3->Draw("colz");
    c->cd(4);
    th4->Draw("colz");

    c->Update();

    if ( x ) delete x;
    if ( y ) delete y;
    if ( xsh ) delete xsh;
    if ( ysh ) delete ysh;
    if ( nlist ) delete nlist;

    cout << "DoWaveLenght ended" << endl;
}

About Us | EUSO Official Website | Web pages created by Roberto Pesce and Alessandro Thea - Last Update Wed Nov 16 16:57:39 2005 Wed Nov 16 16:29:22 2005