Universita' di GenovaINFN Sezione di Genova  
AIRWATCH / EUSO Genova

EAtmosphereHistoPainter - source file

// $Id: EAtmosphereHistoPainter.cc,v 1.18 2005/04/18 14:37:01 moreggia Exp $
// Author: Anne Stutz   2004/10/19

/*****************************************************************************
 * ESAF: Euso Simulation and Analysis Framework                              *
 *                                                                           *
 *  Id: EAtmosphereHistoPainter                                                           *
 *  Package: <packagename>                                                   *
 *  Coordinator: <coordinator>                                               *
 *                                                                           *
 *****************************************************************************/

//_____________________________________________________________________________
//
// EAtmosphereHistoPainters
//
// <extensive class description>
//
//   Config file parameters
//   ======================
//
//   <parameter name>: <parameter description>
//   -Valid options: <available options>
//

#include "EAtmosphereHistoPainter.hh"
#include "EAtmosphere.hh"
#include "EBunchPhotons.hh"
#include "ESystemOfUnits.hh"

#include <TH1F.h>
#include <TF1.h>
#include <TTimer.h>
#include <TPad.h>
#include <TPaveText.h>
#include <TGaxis.h>
#include <TVector3.h>
#include <TMath.h>

ClassImp(EAtmosphereHistoPainter)


// TOOL function
extern Double_t Zv(const TVector3&);
//    return sqrt( pow(6.370949e6*m + v.Z(),2) + v.Perp2() ) - 6.370949e6*m;
//}


//_____________________________________________________________________________
EAtmosphereHistoPainter::EAtmosphereHistoPainter(EAtmosphere* atm) {
    //
    // Constructor
    //

    fAtmosphere = atm;
    fHistos = new TList();

    fTimeBins = 0;
    fAltBins = 0;

    Build();
}

//_____________________________________________________________________________
EAtmosphereHistoPainter::~EAtmosphereHistoPainter() {
    //
    // Destructor
    //

    if ( fTimeBins ) delete[] fTimeBins;
    if ( fAltBins ) delete[] fAltBins;
    if ( fHistos ) fHistos->Delete();
}

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

    if ( !fAtmosphere ) {
        Info("Build()","EAtmosphere object is NULL. Painter made zombie.");
        MakeZombie();
        return;
    }
        
    fNphMax = 0; 
    Int_t n = fAtmosphere->GetNumBunches();
    Int_t ns = fAtmosphere->GetNumSingles();

    if ( n <= 0 ) {
        Info("Build()","No Bunches in Atmosphere ( NumBunch = 0 ). Painter made zombie.");
        MakeZombie();
        return;
    }

    // determine bins for time development
    Double_t temp[n+1];
    fNBins = 0;
    for (Int_t i=0; i<n; i++) {
      if (fAtmosphere->GetBunch(i)->GetType() == 0 ) {
        temp[fNBins] = fAtmosphere->GetBunch(i)->GetDatei()/ms;  
        temp[fNBins+1] = fAtmosphere->GetBunch(i)->GetDatef()/ms;  

        if( fNphMax < fAtmosphere->GetBunch(i)->GetWeight() )
            fNphMax = fAtmosphere->GetBunch(i)->GetWeight();
        fNBins++;
      }
    }

    // determine bins for development w.r.t altitude
    fAltBins = new Double_t[fNBins+1];
    Int_t incrm = 0;
    for (Int_t i=0; i<n; i++) {
      if (fAtmosphere->GetBunch(i)->GetType() == 0 ) {
        fAltBins[fNBins - incrm] = Zv(fAtmosphere->GetBunch(i)->GetShowerPosi())/km;  
        fAltBins[fNBins - 1 - incrm] = Zv(fAtmosphere->GetBunch(i)->GetShowerPosf())/km; 
	incrm++; 
      }
    }

    // determine bins for time on pupil
    fTimeMin = 3000.*microsecond;
    fTimeMax = 0;
    fTimeMaxF = 0;
    for (Int_t i=0; i<ns; i++) {
      //      if ( fAtmosphere->GetSingle(i)->GetType() == 0 && !(fAtmosphere->GetSingle(i)->IsAbsorbed()) ) {
     if ( !(fAtmosphere->GetSingle(i)->IsAbsorbed()) ) {
         if ( fTimeMin > fAtmosphere->GetSingle(i)->GetTof()+fAtmosphere->GetSingle(i)->GetDate() )
            fTimeMin = fAtmosphere->GetSingle(i)->GetTof()+fAtmosphere->GetSingle(i)->GetDate();
         if ( fTimeMax < fAtmosphere->GetSingle(i)->GetTof()+fAtmosphere->GetSingle(i)->GetDate() )
            fTimeMax = fAtmosphere->GetSingle(i)->GetTof()+fAtmosphere->GetSingle(i)->GetDate();
         if ( fAtmosphere->GetSingle(i)->GetType()==0 ) {
             if ( fTimeMaxF < fAtmosphere->GetSingle(i)->GetTof()+fAtmosphere->GetSingle(i)->GetDate() )
                 fTimeMaxF = fAtmosphere->GetSingle(i)->GetTof()+fAtmosphere->GetSingle(i)->GetDate();
         }
      }
    }

    Printf("From EAtmosphere %d bins in time to draw",fNBins);

    fTimeBins = new Double_t[fNBins+1];
    for (Int_t i=0; i<fNBins+1; i++) {
        fTimeBins[i] = temp[i];
    }
    // build histograms from bunches vs time developpement 
    TH1F *th0 = GetTimeHisto("Bunch_Nph", "Nph vs time");
    TH1F *th1 = GetTimeHisto("Bunch_Nph_F", "Nph_Fluo vs time");
    TH1F *th2 = GetTimeHisto("Bunch_Nph_C", "Nph_Cerenkov vs time");
    TH1F *th3 = GetTimeHisto("Bunch_Yield_F","Fluo Yield vs time");
    TH1F *th4 = GetTimeHisto("Bunch_Yield_C","Cerenkov Yield vs time");
    TH1F *th5 = GetTimeHisto("Bunch_Nph_F_SA", "Nph_Fluo_in_Solide_Angle vs time");
    th0->SetYTitle("Nph");
    th1->SetYTitle("Nph Fluo");
    th2->SetYTitle("Nph Cerenkov");
    th3->SetYTitle("Nph/m Fluo");
    th4->SetYTitle("Nph/m Cerenkov");
    th5->SetYTitle("Nph Fluo");

    
    // build histograms from singles vs time developpement
    TH1F *th33 = GetTimeHisto("Single_Nph_F", "Nph_Fluo_in_Solide_Angle vs time");
    TH1F *th34 = GetTimeHisto("Single_Nph_F_P", "Nph_Fluo_on_Pupil vs time");
    TH1F *th35 = GetTimeHisto("Single_Nph_CB", "Nph_Cerenkov_backscttered_in_Solide_Angle vs time");
    TH1F *th35_2 = GetTimeHisto("Single_Nph_ckov_CloudBackscat", "Nph_Cerenkov_CloudBackscat_in_Solide_Angle vs time");
    TH1F *th35_tot = GetTimeHisto("Single_Nph_CB_tot", "Nph_Cerenkov_Total_backscattered_in_Solide_Angle vs time");  // not used so far
    TH1F *th36 = GetTimeHisto("Single_Nph_CB_P", "Nph_Cerenkov_backscattered_on_Pupil vs time");
    TH1F *th36_2 = GetTimeHisto("Single_Nph_ckov_CloudBackscat_P", "Nph_Cerenkov_CloudBackscat_on_Pupil vs time");
    TH1F *th36_tot = GetTimeHisto("Single_Nph_CB_tot_P", "Nph_Cerenkov_Total_backscattered_on_Pupil vs time");   // not used so far
    TH1F *th37 = GetTimeHisto("Single_Nph_CR", "Nph_Cerenkov_reflected_in_Solide_Angle vs time");
    TH1F *th38 = GetTimeHisto("Single_Nph_CR_P", "Nph_Cerenkov_reflected_on_Pupil vs time");
    th33->SetYTitle("Nph Fluo");
    th34->SetYTitle("Nph Fluo");
    th35->SetYTitle("Nph Cerenkov backscattered");
    th35_2->SetYTitle("Nph Cerenkov backscattered by clouds");
    th35_tot->SetYTitle("Nph Cerenkov total backscattered");
    th36->SetYTitle("Nph Cerenkov backscattered ");
    th36_2->SetYTitle("Nph Cerenkov backscattered by clouds");
    th36_tot->SetYTitle("Nph Cerenkov total backscattered");
    th37->SetYTitle("Nph Cerenkov reflecteded");
    th38->SetYTitle("Nph Cerenkov refelected");
    
    // build histo vs altitude
    TH1F *alt0 = GetAltHisto("alt_Bunch_Nph", "Nph vs altitude");
    TH1F *alt1 = GetAltHisto("alt_Bunch_Nph_F", "Nph_Fluo vs altitude");
    TH1F *alt2 = GetAltHisto("alt_Bunch_Nph_C", "Nph_Cerenkov vs altitude");
    TH1F *alt3 = GetAltHisto("alt_Bunch_Yield_F","Fluo Yield vs altitude");
    TH1F *alt4 = GetAltHisto("alt_Bunch_Yield_C","Cerenkov Yield vs altitude");
    alt0->SetYTitle("Nph");
    alt1->SetYTitle("Nph Fluo");
    alt2->SetYTitle("Nph Cerenkov");
    alt3->SetYTitle("Nph/m Fluo");
    alt4->SetYTitle("Nph/m Cerenkov");

    //build histogram vs wavelenght
    TH1F *th11 = GetWlHisto("Bunch_Wl_F", "Nph_F vs Lambda");
    TH1F *th12 = GetWlHisto("Bunch_Wl_C", "Nph_C vs Lambda");

    TH1F *th13 = GetWlHisto("Single_Wl_F", "Nph_Fluo_in_Solid_Angle vs Lambda");
    TH1F *th14 = GetWlHisto("Single_Wl_F_P", "Nph_Fluo_on_pupil vs Lambda");
    TH1F *th15 = GetWlHisto("Single_Wl_CB", "Nph_Cerenkov_back_in_Solid_Angle vs Lambda");
    TH1F *th15_2 = GetWlHisto("Single_Wl_ckov_CloudBackscat", "Nph_Cerenkov_Cloudback_in_Solid_Angle vs Lambda");  // not used so far
    TH1F *th15_tot = GetWlHisto("Single_Wl_CB_tot", "Nph_Cerenkov_Total_back_in_Solid_Angle vs Lambda");
    TH1F *th16 = GetWlHisto("Single_Wl_CB_P", "Nph_Cerenkov_back_on_pupil vs Lambda");
    TH1F *th16_2 = GetWlHisto("Single_Wl_ckov_CloudBackscat_P", "Nph_Cerenkov_Cloudback_on_pupil vs Lambda");
    TH1F *th16_tot = GetWlHisto("Single_Wl_CB_tot_P", "Nph_Cerenkov_Total_back_on_pupil vs Lambda");  // not used so far
    TH1F *th17 = GetWlHisto("Single_Wl_CR","Nph_Cerenkov refl in Solid Angle vs lambda");
    TH1F *th18 = GetWlHisto("Single_Wl_CR_P","Nph_Cerenkov refl on pupil vs lambda");
    th11->SetYTitle("Nph Fluo");
    th12->SetYTitle("Nph Cerenkov");
    th13->SetYTitle("Nph Fluo");
    th14->SetYTitle("Nph Fluo");
    th15->SetYTitle("Nph Cerenkov");
    th15_2->SetYTitle("Nph clouds Cerenkov");
    th15_tot->SetYTitle("Nph total Cerenkov");
    th16->SetYTitle("Nph Cerenkov");
    th16_2->SetYTitle("Nph clouds Cerenkov");
    th16_tot->SetYTitle("Nph total Cerenkov");
    th17->SetYTitle("Nph Cerenkov");
    th18->SetYTitle("Nph Cerenkov");

    // build histogram vs time on pupil
    TH1F *th20 = GetTofHisto("Single_Nph_P_t","Nph_Total_on_Pupil vs time");
    TH1F *th24 = GetTofHisto("Single_Nph_F_P_t","Nph_Fluo_on_Pupil vs time");
    TH1F *th26 = GetTofHisto("Single_Nph_CB_P_t","Nph_Cerenkov back on pupil vs time");
    TH1F *th26_2 = GetTofHisto("Single_Nph_ckov_CloudBackscat_P_t","Nph_Cerenkov cloudback on pupil vs time");
    TH1F *th26_tot = GetTofHisto("Single_Nph_ckov_Total_Backscat_P_t","Nph_Cerenkov total back on pupil vs time");
    TH1F *th28 = GetTofHisto("Single_Nph_CR_P_t","Nph_Cerenkov refl on pupil vs time");
    th20->SetYTitle("Nph ");
    th24->SetYTitle("Nph Fluo");
    th26->SetYTitle("Nph Backscattered Cerenkov");
    th26_2->SetYTitle("Nph CloudBackscattered Cerenkov");
    th26_tot->SetYTitle("Nph total Backscattered Cerenkov");
    th28->SetYTitle("Nph Refelected Cerenkov");

    //fill histograms 
    TVector3 euso(0,0,430*km);
    Double_t radius = 1250 * mm;
    Info("EAtmosphereHistoPainter::Build()","USE pupil radius = %f mm, euso altitude = %f km",radius/mm,euso.Z()/km);

    // fill histograms from bunches
    for (Int_t i=0; i<n; i++) { 
        Double_t time = ( fAtmosphere->GetBunch(i)->GetDatei()  
                        + fAtmosphere->GetBunch(i)->GetDatef() )/2/ms;
	Double_t alt =( Zv(fAtmosphere->GetBunch(i)->GetShowerPosi())  
                        + Zv(fAtmosphere->GetBunch(i)->GetShowerPosf()) )/2/km; 
        Double_t Nph =      fAtmosphere->GetBunch(i)->GetWeight(); 
        Double_t Yield =    fAtmosphere->GetBunch(i)->GetYield()*m;
        Double_t NumWl =    fAtmosphere->GetBunch(i)->GetNumWavelengths();
        const Float_t* lambda = 0;
        const Float_t* wlweight = 0;
        wlweight = fAtmosphere->GetBunch(i)->GetTable("weight");
        lambda = fAtmosphere->GetBunch(i)->GetTable("lambda");

	th0->Fill(time,Nph); 
	alt0->Fill(alt,Nph);

        if (fAtmosphere->GetBunch(i)->GetType()==0) {
	   // fluorescence bunches
           th1->Fill(time,Nph);  
           th3->Fill(time,Yield);
	   alt1->Fill(alt,Nph);
	   alt3->Fill(alt,Yield);
           for ( Int_t j=0; j<NumWl; j++ ) {
	       th11->Fill( lambda[j]/nm , wlweight[j]*Nph );  
	   }

           TVector3 dist = euso - (fAtmosphere->GetBunch(i)->GetShowerPosi()+
                               fAtmosphere->GetBunch(i)->GetShowerPosf())*0.5;
           Double_t omega = TMath::Pi() *radius*radius *cos( dist.Theta() ) / dist.Mag2();
           th5->Fill(time,Nph*omega/(4*TMath::Pi()));
        } 

        if (fAtmosphere->GetBunch(i)->GetType()==1) { 
	   // Cerenkov bunches
           th2->Fill(time,Nph); 
           th4->Fill(time,Yield);
	   alt2->Fill(alt,Nph);
	   alt4->Fill(alt,Yield); 
           for ( Int_t j=0; j<NumWl; j++ ) {
	       th12->Fill( lambda[j]/nm , wlweight[j]*Nph );  
	   }
        }  
    } 

    // fill histograms from singles
    for (Int_t i=0; i<ns; i++) { 
 
           Double_t time =  fAtmosphere->GetSingle(i)->GetDate()/ms;
           Double_t tof = (fAtmosphere->GetSingle(i)->GetTof() + 
                           fAtmosphere->GetSingle(i)->GetDate())/microsecond; 
           Double_t wl = fAtmosphere->GetSingle(i)->GetWl()/nm;

           if ( !(fAtmosphere->GetSingle(i)->IsAbsorbed()) ) th20->Fill(tof,1);

           if ( fAtmosphere->GetSingle(i)->GetType() == 0 ) { 
              // for fluorescence photons
	      th33->Fill(time,1); 
              th13->Fill(wl,1);
              if ( !(fAtmosphere->GetSingle(i)->IsAbsorbed()) ) {
                  th34->Fill(time,1);
                  th14->Fill(wl,1);
                  th24->Fill(tof,1);
              }
           }
           if ( fAtmosphere->GetSingle(i)->GetType() == 1 ) { 
	    // for Backscattered Cerenkov photons 
              if ( fAtmosphere->GetSingle(i)->GetHistory()==2) {
		  th35->Fill(time,1);
		  th35_tot->Fill(time,1);
		  th15->Fill(wl,1);
		  th15_tot->Fill(wl,1);
                  if ( !(fAtmosphere->GetSingle(i)->IsAbsorbed()) ) { 
		      th36->Fill(time,1);
		      th36_tot->Fill(time,1);
		      th16->Fill(wl,1);
		      th16_tot->Fill(wl,1);
                      th26->Fill(tof,1);
                      th26_tot->Fill(tof,1);
                  }
              }
              if ( fAtmosphere->GetSingle(i)->GetHistory()==3) {
		  th35_2->Fill(time,1);
		  th35_tot->Fill(time,1);
		  th15_2->Fill(wl,1);
		  th15_tot->Fill(wl,1);
                  if ( !(fAtmosphere->GetSingle(i)->IsAbsorbed()) ) { 
		      th36_2->Fill(time,1);
		      th36_tot->Fill(time,1);
		      th16_2->Fill(wl,1);
		      th16_tot->Fill(wl,1);
                      th26_2->Fill(tof,1);
                      th26_tot->Fill(tof,1);
                  }
              }
              // for reflected Cerenkov photons
              if ( fAtmosphere->GetSingle(i)->GetHistory()==1 ) {
		  th37->Fill(time,1);
		  th17->Fill(wl,1);
                  if ( !(fAtmosphere->GetSingle(i)->IsAbsorbed()) ) { 
		      th38->Fill(time,1);
		      th18->Fill(wl,1);
                      th28->Fill(tof,1);
                  }
              }
           }

	}

} 

//_____________________________________________________________________________
 TH1F *EAtmosphereHistoPainter::GetTimeHisto(const char *name, const char *title) {
    //
    // Get histograms versus time
    //

    if ( title == 0 ) title = name;
    TH1F* th = (TH1F*)fHistos->FindObject(name);

    if (!th) {
        th = new TH1F( name, title, fNBins, fTimeBins);
        fHistos->Add(th);
        th->SetStats(0);
        th->SetXTitle("Development Time, #{m}s");
        th->SetFillColor(9);
    }

    return th;
}

//_____________________________________________________________________________
 TH1F *EAtmosphereHistoPainter::GetAltHisto(const char *name, const char *title) {
    //
    // Get histograms versus altitude
    //

    if ( title == 0 ) title = name;
    TH1F* th = (TH1F*)fHistos->FindObject(name);

    if (!th) {
        th = new TH1F( name, title, fNBins, fAltBins);
        fHistos->Add(th);
        th->SetStats(0);
        th->SetXTitle("Altitude (km)");
        th->SetFillColor(9);
    }

    return th;
}

//_____________________________________________________________________________
 TH1F *EAtmosphereHistoPainter::GetWlHisto(const char *name, const char *title) {
    //
    // Get wavelength histograms 
    //

    if ( title == 0 ) title = name;
    TH1F* th = (TH1F*)fHistos->FindObject(name);

    if (!th) {
      th = new TH1F( name, title, 150, 300., 450.);
        fHistos->Add(th);
        th->SetStats(0);
        th->SetXTitle("Wavelength, #nm");
        th->SetFillColor(9);
    }

    return th;
}

//_____________________________________________________________________________
 TH1F *EAtmosphereHistoPainter::GetTofHisto(const char *name, const char *title) {
    //
    // Get wavelength histograms 
    //

    if ( title == 0 ) title = name;
    TH1F* th = (TH1F*)fHistos->FindObject(name);

    Float_t min = (Int_t)(fTimeMin/microsecond)-1;
    Int_t nbins = (Int_t)((fTimeMax-fTimeMin)/(2.5*microsecond)) + 20;    

    if (!th) {
      th = new TH1F( name, title, nbins, min, min+nbins*2.5);
        fHistos->Add(th);
        th->SetStats(0);
        th->SetXTitle("Time, #mus");
        th->SetFillColor(9);
    }

    return th;
}

//______________________________________________________________________________
 void EAtmosphereHistoPainter::Draw(Option_t *opt) {
    //
    //
    //
    
    if ( IsZombie() ) return;

    TString option(opt);
    TH1F* th(0);
    TH1F* th1(0);
    TH1F* th2(0);
    TH1F* th3(0);
    TH1F* th4(0);

    if ( option == "" ) option = "Nph";

    if ( option == "Nph" ) {
        th2 = GetTimeHisto("Bunch_Nph_C");
        if ( th2 ) {
            th2->SetFillColor(39);
            th2->Draw();
        }       
        th1 = GetTimeHisto("Bunch_Nph_F");        
        if ( th1 ) {
           th1->SetFillColor(29);
           th1->Draw("SAME"); 
        }  
    } 

    if ( option == "Nph_F" ) {
        th = GetTimeHisto("Bunch_Nph_F");        
        if ( th ) {
           th->SetFillColor(29);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           sprintf(txt,"Total:%e",th->Integral());
           t->DrawText(0.005,th->GetMaximum()*0.94,txt);
        }  
    }
    
    // fluo photons created vs/ altitude
    if ( option == "alt_Nph_F" ) {
        th = GetAltHisto("alt_Bunch_Nph_F");        
        if ( th ) {
           th->SetFillColor(29);
           TF1* mygaus = new TF1("mygaus","gaus",0,10);
           th->Fit("mygaus","R");
           Double_t par[3];
           mygaus->GetParameters(par);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           sprintf(txt,"Total:%e",th->Integral());
           t->DrawText(5,th->GetMaximum()*0.94,txt);
           sprintf(txt,"Alt max : %e",par[1]);
           t->DrawText(5,th->GetMaximum()*0.84,txt);
           sprintf(txt,"Nph max : %e",par[0]);
           t->DrawText(5,th->GetMaximum()*0.74,txt);
        }  
    }
    // ckov created, or fluo yield, or ckov yield
    if ( option == "alt_Nph_C" ) {
        th = GetTimeHisto("alt_Bunch_Nph_C");        
        if ( th ) {
           th->SetFillColor(39);
           TF1* mygaus = new TF1("mygaus","gaus",0,10);
	   th->Fit("mygaus","R");
           Double_t par[3];
           mygaus->GetParameters(par);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           sprintf(txt,"Total:%e",th->Integral());
           t->DrawText(5,th->GetMaximum()*0.94,txt);
           sprintf(txt,"Alt max : %e",par[1]);
           t->DrawText(5,th->GetMaximum()*0.84,txt);
           sprintf(txt,"Nph max : %e",par[0]);
           t->DrawText(5,th->GetMaximum()*0.74,txt);
        } 
    }
    else if ( option == "alt_Yield_F" ) {
        th = GetTimeHisto("alt_Bunch_Yield_F");
        if ( th ) {
           th->SetFillColor(29);
	   th->SetMinimum(0.9*th->GetMinimum());
           th->Draw();
	}
    }            
    else if ( option == "alt_Yield_C" ) {
        th = GetTimeHisto("alt_Bunch_Yield_C");  
        if ( th ) {
           th->SetFillColor(39);
           th->Draw();
	}
    } 

    if ( option == "Nph_F_SA" ) {
        th = GetTimeHisto("Bunch_Nph_F_SA");        
        if ( th ) {
           th->SetFillColor(29);
            
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           sprintf(txt,"Total:%e",th->Integral());
           t->DrawText(0.005,th->GetMaximum()*0.94,txt);
        }  
    }


    if ( option == "Nph_C" ) {
        th = GetTimeHisto("Bunch_Nph_C");        
        if ( th ) {
           th->SetFillColor(39);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           sprintf(txt,"Total:%e",th->Integral());
           t->DrawText(0.005,th->GetMaximum()*0.94,txt);
        } 
    }
    else if ( option == "Yield_F" ) {
        th = GetTimeHisto("Bunch_Yield_F");
        if ( th ) {
           th->SetFillColor(29);
	   th->SetMinimum(0.9*th->GetMinimum());
           th->Draw();
	}
    }            
    else if ( option == "Yield_C" ) {
        th = GetTimeHisto("Bunch_Yield_C");  
        if ( th ) {
           th->SetFillColor(39);
           th->Draw();
	}
    } 



    if ( option == "Sp_F" ) {
        th = GetWlHisto("Bunch_Wl_F");        
        if ( th ) {
           th->SetFillColor(29);
           th->Draw(""); 
        }  
    }

    if ( option == "Sp_C" ) {
        th = GetWlHisto("Bunch_Wl_C");        
        if ( th ) {
           th->SetFillColor(39);
           th->Draw(""); 
        }  
    }

    if ( option == "Nphs_F" ) {
        th = GetTimeHisto("Single_Nph_F"); 
        th1 = GetTimeHisto("Single_Nph_F_P");       
        if ( th ) {
           th->SetFillColor(29);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           t->SetTextColor(29);
           sprintf(txt,"Total:%d",(int)th->Integral());
           t->DrawText(0.005,th->GetMaximum()*0.94,txt);
        
           if ( th1 ) {
              th1->SetFillColor(2);
              th1->Draw("SAME"); 
              t->SetTextColor(2);
              sprintf(txt,"Total:%d",(int)th1->Integral());
              t->DrawText(0.005,th1->GetMaximum()*0.94,txt);
          }  
        }  
    }

    // Fluorescence photons vs time on pupil
    if ( option == "Nphs_F_t" ) { 
        th = GetTofHisto("Single_Nph_F_P_t");
        Double_t tx = fTimeMin/microsecond+10;       
        if ( th ) {
            th->SetFillColor(2);
	    TF1* mygaus = new TF1("mygaus","gaus",fTimeMin/microsecond,fTimeMaxF/microsecond);
	    th->Fit("mygaus","R");
	    Double_t par[3];
	    mygaus->GetParameters(par);
            th->Draw(""); 
            TText *t = new TText();
            char txt[200];
            t->SetTextFont(62);  
            t->SetTextSize(0.06);
            t->SetTextColor(2);
            sprintf(txt,"Total:%d",(int)th->Integral());
            t->DrawText(tx,th->GetMaximum()*0.94,txt);
            t->SetTextColor(1);
	    sprintf(txt,"DTmax : %d",(int)(par[1]-fTimeMin/microsecond));
            t->DrawText(tx,th->GetMaximum()*0.84,txt);
	    sprintf(txt,"Nph max : %d",(int)par[0]);
            t->DrawText(tx,th->GetMaximum()*0.74,txt);
        }    
    }

    // fluorescence spectrum
    if ( option == "Sps_F" ) {
        th = GetWlHisto("Single_Wl_F"); 
        th1 = GetWlHisto("Single_Wl_F_P");       
        if ( th ) {
           th->SetFillColor(29);
           th->Draw();   
           if( th1 ) { 
              th1->SetFillColor(2); 
              th1->Draw("SAME"); 
           }
        }  
    }

    // Cerenkov photons
    if ( option == "Nphs_CR" ) {
        th = GetTimeHisto("Single_Nph_CR"); 
        th1 = GetTimeHisto("Single_Nph_CR_P");       
        if ( th ) {
           th->SetFillColor(39);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           t->SetTextColor(39);
           sprintf(txt,"Total:%d",(int)th->Integral());
           t->DrawText(0.005,th->GetMaximum()*0.94,txt);
        
           if ( th1 ) {
              th1->SetFillColor(50);
              th1->Draw("SAME"); 
              t->SetTextColor(50);
              sprintf(txt,"Total:%d",(int)th1->Integral());
              t->DrawText(0.005,th1->GetMaximum()*0.94,txt);
          }  
        }  
    }
    if ( option == "Nphs_CB" ) {
        //th = GetTimeHisto("Single_Nph_CB"); 
	th = GetTimeHisto("Single_Nph_CB_tot");  // total ckov backscat before last transmission
        th1 = GetTimeHisto("Single_Nph_CB_P");   // air backscat on pupil
        th2 = GetTimeHisto("Single_Nph_ckov_CloudBackscat_P");  // clouds backscat on pupil
        if ( th ) {
           th->SetFillColor(39);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           t->SetTextColor(39);
           sprintf(txt,"Total:%d",(int)th->Integral());
           t->DrawText(0.005,th->GetMaximum()*0.94,txt);
           if ( th1 ) {
              th1->SetFillColor(50);
              th1->Draw("SAME"); 
              t->SetTextColor(50);
              sprintf(txt,"Total:%d",(int)th1->Integral());
              t->DrawText(0.005,th1->GetMaximum()*0.94,txt);
          }
           if ( th2 ) {
              th2->SetFillColor(5);
              th2->Draw("SAME"); 
              t->SetTextColor(5);
              sprintf(txt,"Total:%d",(int)th2->Integral());
              t->DrawText(0.005,th2->GetMaximum()*0.94,txt);
          }  
  
        } 
    }

    if ( option == "Nphs_CR_t" ) {
        th1 = GetTofHisto("Single_Nph_CR_P_t"); 

        Double_t tx = fTimeMin/microsecond+10;    
        if ( th1 ) {
           th1->SetFillColor(50);
           th1->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           t->SetTextColor(50);
           sprintf(txt,"Total:%d",(int)th1->Integral());
           t->DrawText(tx,th1->GetMaximum()*0.94,txt);       
        }  
    }
    if ( option == "Nphs_CB_t" ) {
        th = GetTofHisto("Single_Nph_ckov_Total_Backscat_P_t"); 
	th2 = GetTofHisto("Single_Nph_CB_P_t"); 
        th1 = GetTofHisto("Single_Nph_ckov_CloudBackscat_P_t"); 

        Double_t tx = fTimeMin/microsecond+10;    
        if ( th ) {
           th->SetFillColor(39);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           t->SetTextColor(39);
           sprintf(txt,"Total:%d",(int)th->Integral());
           t->DrawText(0.005,th->GetMaximum()*0.94,txt);
           if ( th1 ) {
              th1->SetFillColor(5);
              th1->Draw("SAME"); 
              t->SetTextColor(5);
              sprintf(txt,"Cloud scattered:%d",(int)th2->Integral());
              t->DrawText(0.005,th2->GetMaximum()*0.94,txt);
          }  
           if ( th2 ) {
              th2->SetFillColor(50);
              th2->Draw("SAME"); 
              t->SetTextColor(50);
              sprintf(txt,"Air scattered:%d",(int)th2->Integral());
              t->DrawText(0.005,th2->GetMaximum()*0.94,txt);
          }  
        }  
    }

    // Cerenkov spectrum
    if ( option == "Sps_CR" ) {
        th = GetWlHisto("Single_Wl_CR"); 
        th1 = GetWlHisto("Single_Wl_CR_P"); 
        if ( th ) {
           th->SetFillColor(39);
           th->Draw();   
           if( th1 ) { 
              th1->SetFillColor(50); 
              th1->Draw("SAME"); 
           }
        }  
    }
    if ( option == "Sps_CB" ) {
        //th = GetWlHisto("Single_Wl_CB"); 
        th = GetWlHisto("Single_Wl_CB_tot"); 
        th1 = GetWlHisto("Single_Wl_CB_P"); 
        th2 = GetWlHisto("Single_Wl_ckov_CloudBackscat_P"); 
        if ( th ) {
           th->SetFillColor(39);
           th->Draw();   
           if( th1 ) { 
              th1->SetFillColor(50); 
              th1->Draw("SAME"); 
           }
           if( th2 ) { 
              th2->SetFillColor(5); 
              th2->Draw("SAME"); 
           }
        }  
    }

    // Photons on Pupil
    if ( option == "Nphs_t" ) {
        th = GetTofHisto("Single_Nph_P_t");
        th1 = GetTofHisto("Single_Nph_F_P_t"); 
        th2 = GetTofHisto("Single_Nph_CR_P_t"); 
        th3 = GetTofHisto("Single_Nph_CB_P_t"); 
        th4 = GetTofHisto("Single_Nph_ckov_CloudBackscat_P_t"); 
        Double_t tx = fTimeMin/microsecond+10;    
        if ( th ) {
           th->SetFillColor(9);
           th->Draw(""); 
           TText *t = new TText();
           char txt[200];
           t->SetTextFont(62);  
           t->SetTextSize(0.06);
           t->SetTextColor(9);
           sprintf(txt,"Total:%d",(int)th->Integral());
           t->DrawText(tx,th->GetMaximum()*0.94,txt);
           if ( th1 ) {
              th1->SetFillColor(2);
              th1->Draw("SAME"); 
              t->SetTextColor(2);
              sprintf(txt,"Fluo:%d",(int)th1->Integral());
              t->DrawText(tx,th1->GetMaximum()*0.94,txt);
           } 
           if ( th3 ) {
              th3->SetFillColor(3);
              th3->Draw("SAME"); 
              t->SetTextColor(3);
              sprintf(txt,"Air scattered:%d",(int)th3->Integral());
              t->DrawText(tx,th3->GetMaximum()*0.94,txt);
           }       
           if ( th2 ) {
              th2->SetFillColor(6);
              th2->Draw("SAME"); 
              t->SetTextColor(6);
              sprintf(txt,"Reflected:%d",(int)th2->Integral());
              t->DrawText(tx,th2->GetMaximum()*0.5,txt);
           }       
           if ( th4 ) {
              th4->SetFillColor(5);
              th4->Draw("SAME"); 
              t->SetTextColor(5);
              sprintf(txt,"Cloud scattered:%d",(int)th4->Integral());
              t->DrawText(tx,th4->GetMaximum()*0.5,txt);
           }       
        }  
    }
    
}

//_____________________________________________________________________________
 void EAtmosphereHistoPainter::NextFrame() {
    //
    // Switch to the next frame and update histograms
    // 

    fFrame++;
    UpdateHistos();
}

//______________________________________________________________________________
 void EAtmosphereHistoPainter::UpdateHistos() {
    //
    // Update hitograms at the current frame.
    //

    if ( IsZombie() ) return;

    Int_t n;
    if ( fFrame == -1 )
        n = fAtmosphere->GetNumBunches();
    else
        n = (Int_t)(((Float_t)fAtmosphere->GetNumBunches()/fNumFrames)*fFrame);


    TH1F* thNph    = GetTimeHisto("Bunch_Nph");
    // TH1F* thNph_F   = GetTimeHisto("Bunch_Nph_F");
    // TH1F* thNph_C = GetTimeHisto("Bunch_Nph_C");

    for (Int_t i=0; i<fAtmosphere->GetNumBunches() ; i++) {
        Double_t Nph = fAtmosphere->GetBunch(i)->GetWeight();

        thNph->SetBinContent(i+1, i < n ? Nph : 0 );
	//	thNph_F->SetBinContent(i+1, i < n ? Nph_F : 0 );
        //thNph_C->SetBinContent(i+1, i < n ? Nph_C : 0 );
    }

    thNph->SetMaximum( fNphMax*1.1 );
    //thNph_F->SetMaximum( fNphFMax*1.1 );
    //thNph_C->SetMaximum( fNphCMax*1.1 );

}
About Us | EUSO Official Website | Web pages created by Roberto Pesce and Alessandro Thea - Last Update 14-May-2005 21:31