Universita' di GenovaINFN Sezione di Genova  
AIRWATCH / EUSO Genova

EnergyViewer - source file

// $Id: EnergyViewer.cc,v 1.24 2005/11/14 17:02:49 naumov Exp $
// Author: Dmitry.Naumov   2005/09/25

/*****************************************************************************
 * ESAF: Euso Simulation and Analysis Framework                              *
 *                                                                           *
 *  Id: EnergyViewer                                                         *
 *  Package:  Energy                                                         *
 *  Coordinator: Dmitry.Naumov and Artem Chukanov                            *
 *                                                                           *
 *****************************************************************************/

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

#include "EnergyViewer.hh"
#include "EnergyModule.hh"
#include "RecoPixelData.hh"
#include "RecoEvent.hh"
#include "Etypes.hh"
#include "EPmtGeo.hh"

#include <TCanvas.h>
#include <TPad.h>
#include <TGraph.h>
#include <TPolyLine.h>
#include <TPolyLine3D.h>
#include <TPolyMarker.h>
#include <TStyle.h>
#include <TColor.h>
#include <TText.h>
#include <TF1.h>
#include <TView.h>
#include <TProfile.h>
#include <TSystem.h>

using namespace TMath;

ClassImp(EnergyViewer)
  
//_____________________________________________________________________________
 EnergyViewer::EnergyViewer(EnergyModule* fEM) : EsafMsgSource() {
  //
  // Constructor
  //
    fCanvas           = NULL;
    fEnergyModule     = fEM;
    fColNum           = 200;
    fColorDefined     = kFALSE;
    fFSViewBuild      = kFALSE; 
    fDefinePixelsView = kFALSE;

    DefineColors();
    gStyle->SetCanvasColor(kWhite);
}

//_____________________________________________________________________________
 EnergyViewer::~EnergyViewer() {
    //
    // Destructor
    //
    SafeDelete(fColor);
    SafeDelete(fCanvas);
}


//_____________________________________________________________________________
 void EnergyViewer::InitFilePS() {
  fCanvas = new TCanvas("EnergyModule");
  fPsFile = Form("output/debug.run%d.event%d.ps",
		 fEnergyModule->fEvent->GetHeader().GetRun(),
		 fEnergyModule->fEvent->GetHeader().GetNum());
  fCanvas->Print(fPsFile+"[");
}

//_____________________________________________________________________________
 void EnergyViewer::CloseFilePS() {

  if (!fCanvas) return;
  
  fCanvas->Print(fPsFile+"]");
  gSystem->Exec(Form("rm -rf %s.gz ; gzip -9 %s", fPsFile.Data(), fPsFile.Data()));
  Msg(EsafMsg::Info) << fPsFile << " is gzipped" << MsgDispatch;
  fCanvas->Clear();
}

//_____________________________________________________________________________
 void EnergyViewer::DisplayFSTrack() {

    if (!fCanvas) return;

    fEnergyModule->fGraphTheta->GetFunction("pol2")->SetLineColor(kMagenta);
    fEnergyModule->fGraphPhi->GetFunction("pol2")->SetLineColor(kMagenta);
    fEnergyModule->fGraphX->GetFunction("pol1")->SetLineColor(kMagenta);
    fEnergyModule->fGraphY->GetFunction("pol1")->SetLineColor(kMagenta);
    fEnergyModule->fGraphXY->GetFunction("pol1")->SetLineColor(kMagenta);

    fCanvas->Clear(); 
    fCanvas->Divide(1,2); fCanvas->cd(1);
    SetHistoStyle(fEnergyModule->fHistoTheta);
    fEnergyModule->fHistoTheta->Draw();
    fEnergyModule->fGraphTheta->Draw("X");
    fEnergyModule->fHistoTheta->GetYaxis()->SetTitle("#Theta, [deg]");
    fEnergyModule->fHistoTheta->GetXaxis()->SetTitle("Gtu Hit");
    
    fCanvas->cd(2);
    SetHistoStyle(fEnergyModule->fHistoPhi);
    fEnergyModule->fHistoPhi->Draw();
    fEnergyModule->fGraphPhi->Draw("X");
    fEnergyModule->fHistoPhi->GetYaxis()->SetTitle("#Phi, [deg]");
    fEnergyModule->fHistoPhi->GetXaxis()->SetTitle("Gtu Hit");
    
    fCanvas->Print(fPsFile);
    
    fCanvas->Clear();
    fCanvas->Divide(1,2); fCanvas->cd(1);
    SetHistoStyle(fEnergyModule->fHistoX);
    fEnergyModule->fHistoX->Draw();
    fEnergyModule->fGraphX->Draw("X");
    fEnergyModule->fHistoX->GetYaxis()->SetTitle("X, [mm]");
    fEnergyModule->fHistoX->GetXaxis()->SetTitle("Gtu Hit");
    
    fCanvas->cd(2);
    SetHistoStyle(fEnergyModule->fHistoY);
    fEnergyModule->fHistoY->Draw();
    fEnergyModule->fGraphY->Draw("X");
    fEnergyModule->fHistoY->GetYaxis()->SetTitle("Y, [mm]");
    fEnergyModule->fHistoY->GetXaxis()->SetTitle("Gtu Hit");

    fCanvas->Print(fPsFile);

    fCanvas->Clear();
    SetHistoStyle(fEnergyModule->fHistoXY);
    fEnergyModule->fHistoXY->SetMarkerColor(kGreen);
    fEnergyModule->fHistoXY->SetMarkerStyle(kMultiply);
    fEnergyModule->fHistoXY->SetMarkerSize(1.e-3*Min(gPad->GetWh(),gPad->GetWw()));
    fEnergyModule->fHistoXY->Draw();
    fEnergyModule->fGraphXY->Draw("X");
    fCanvas->Print(fPsFile);
    fCanvas->Clear();
}


//_____________________________________________________________________________
 void EnergyViewer::DisplayRecoveryDebug(TH1* Psf, Float_t IntegratedEfficacy, Float_t time) {

  if (!fCanvas) return;
  if (ceil(time) - time) return;

  DefinePixelsView();

  Pixels_t OneGtuPixels;
  map <Int_t, TPolyLine>::iterator itPmt;
  Pixels_t::iterator itPix;

  TPaveText *cl = NULL;
  TText *text = new TText();

  fCanvas->Clear();

  if (gPad) {
    Float_t OffSet = 5; 
    gPad->Range(fXmin-OffSet,fYmin-OffSet,fXmax+OffSet,fYmax+OffSet);
  }

  for (itPmt = fMyPixelLines.begin(); itPmt != fMyPixelLines.end(); itPmt++)
    itPmt->second.Draw();

  fEnergyModule->fHistoXY->SetStats(kFALSE);
  fEnergyModule->fHistoXY->Draw("same");
  Psf->Draw("samecolz");
  if (fEnergyModule->fGraphXY) fEnergyModule->fGraphXY->Draw("X");

  OneGtuPixels.clear();
  OneGtuPixels = fEnergyModule->fGtuMyPixels[Int_t(time)];

  for (itPix = OneGtuPixels.begin(); itPix != OneGtuPixels.end(); itPix++) {
    Int_t Count = (*itPix)->GetCounts();
    Int_t PixUid = (*itPix)->GetPixelId();
    if (Count > 0) {
      fMyPixelLines[PixUid].SetFillColor(fPalette[Count]);
      fMyPixelLines[PixUid].Draw("f");
      fMyPixelLines[PixUid].Draw();
      TVector3 PixelCenter = fEnergyModule->fRunpars->PixelCenter(PixUid);
      Float_t x = PixelCenter.X()-0.2*fEnergyModule->fRunpars->GetPmtData().GetPadSide();
      Float_t y = PixelCenter.Y()-0.2*fEnergyModule->fRunpars->GetPmtData().GetPadSide();
      Float_t padMin = Min(gPad->GetWh(),gPad->GetWw());
      text->SetTextSize(5.e-5*padMin);
      text->SetTextColor(5);
      text->DrawText(x,y,Form("%d",Count));
    }
  }

  TString pavetext = Form("Efficacies: total = %f mm^{2}, integrated %f mm^{2} at %3.2f ms",Psf->Integral(), IntegratedEfficacy, time);

  DrawLabel(cl,pavetext);

  fCanvas->Print(fPsFile);

  delete cl;
  delete text;
}

//_____________________________________________________________________________
 void EnergyViewer::DisplayRecovery() {
    //
    //
    //
    if (!fCanvas) return;
    fCanvas->Clear();
    fCanvas->Divide(2,2); fCanvas->cd(1);

    fEnergyModule->fHistoOverlap->Draw("histo");
    fEnergyModule->fHistoOverlap->GetYaxis()->SetTitle("Psf Integrated over Gtu Overlap with pixels, [mm^{2}]");
    fEnergyModule->fHistoOverlap->GetXaxis()->SetTitle("Gtu Hit");
    
    fCanvas->cd(2);
    fEnergyModule->fHistoIdealOverlap->Draw("histo");
    fEnergyModule->fHistoIdealOverlap->GetYaxis()->SetTitle("Psf Ideal Overlap with pixels, [mm^{2}]");
    fEnergyModule->fHistoIdealOverlap->GetXaxis()->SetTitle("Gtu Hit");
    
//     map <Int_t, TPolyLine>::iterator itPix;
//     fCanvas->cd(3);
//     for (itPix = fMyPixelLines.begin(); itPix != fMyPixelLines.end(); itPix++) {
//          itPix->second.Draw("f");
//          itPix->second.Draw();
//     }
//     if (fEnergyModule->fGraphXY) fEnergyModule->fGraphXY->Draw("X");
   
    fCanvas->Print(fPsFile);
}

//_____________________________________________________________________________
 void EnergyViewer::DefineColors() {

  if (fColorDefined) return;

  Int_t Red, Green, Blue;

  fPalette[0] = 1;

  for (Int_t iColor = 1; iColor < 1000; iColor++) {
    if (iColor <= fColNum) {

      if (iColor < 0.3*fColNum) {Blue = 1; Green = Red = 0;}
      else if (iColor > 0.7*fColNum) {Red = 1; Green = Blue = 0;}
      else {Green = 1; Red = Blue = 0;}
      
      if(!gROOT->GetColor(300+iColor))
	fColor = new TColor(300+iColor,
			    Red*(0.6+0.4*iColor/fColNum),
			    Green*(0.6+0.5*iColor/fColNum),
			    Blue*(0.6+1.*iColor/fColNum),"");
      else {
	fColor = gROOT->GetColor(300+iColor);
	fColor->SetRGB(Red*(0.6+0.4*iColor/fColNum),
		       Green*(0.6+0.5*iColor/fColNum),
		       Blue*(0.6+1.*iColor/fColNum));
      }
      fPalette[iColor] = 300 + iColor;
    }
    else 
      fPalette[iColor] = 300 + fColNum;
  }

  fColorDefined = kTRUE;
}

//_____________________________________________________________________________
 void EnergyViewer::DefinePixelsView() {

    if (!fCanvas || fDefinePixelsView) return;

    PMTs_t::iterator itPmt;

    fMyPixelLines.clear();
    fPixelLines.clear(); fPixelLines3D.clear(); 
    fPmtLines.clear();   fPmtLines3D.clear();

    fXmin =  fEnergyModule->fXmin;
    fXmax =  fEnergyModule->fXmax;
    fYmin =  fEnergyModule->fYmin;
    fYmax =  fEnergyModule->fYmax;
    fZmin =  kHuge;
    fZmax = -kHuge;
  
    Float_t x[5], y[5], p[15];

    // make list of lines to draw around each pixels with signal

    for (itPmt = fEnergyModule->fMyPixels.begin(); itPmt != fEnergyModule->fMyPixels.end(); itPmt++) {
      for(Int_t iCorner = 0; iCorner < 4; iCorner++) {
	TVector3 vPix = fEnergyModule->fRunpars->PixelCorner(*itPmt,iCorner);
	x[iCorner] = vPix.X(); y[iCorner] = vPix.Y();
      }
      x[4] = x[0]; y[4] = y[0];
      fMyPixelLines[*itPmt].SetPolyLine(5,x,y);
    }
    
    // make list of lines to draw around each PMT
  
    for (itPmt = fEnergyModule->fPmts.begin(); itPmt != fEnergyModule->fPmts.end(); itPmt++) {
      for(Int_t iCorner(0); iCorner < 4; iCorner++) {
	TVector3 vPmt = fEnergyModule->fRunpars->PmtCorner(*itPmt,iCorner);

	p[3*iCorner]   =  x[iCorner] =  vPmt.X();
	p[3*iCorner+1] =  y[iCorner] =  vPmt.Y();
	p[3*iCorner+2] = -vPmt.Z();

	if (fXmin >  vPmt.X()) fXmin =  vPmt.X();
	if (fYmin >  vPmt.Y()) fYmin =  vPmt.Y();
	if (fZmin > -vPmt.Z()) fZmin = -vPmt.Z();
	if (fXmax <  vPmt.X()) fXmax =  vPmt.X();
	if (fYmax <  vPmt.Y()) fYmax =  vPmt.Y();
	if (fZmax < -vPmt.Z()) fZmax = -vPmt.Z();
      }

      p[12] = x[4] = x[0]; p[13] = y[4] = y[0]; p[14] = p[2];

      fPmtLines[*itPmt].SetPolyLine(5,x,y);
      fPmtLines[*itPmt].SetLineColor(20);

      fPmtLines3D[*itPmt].SetPolyLine(5,p);
      fPmtLines3D[*itPmt].SetLineColor(kBlack);
    }

    // make list of lines to draw around each pixels

    for (itPmt = fEnergyModule->fPixelsBgr.begin(); itPmt != fEnergyModule->fPixelsBgr.end(); itPmt++) {
      for(Int_t iCorner(0); iCorner < 4; iCorner++) {
	TVector3 vPmt = fEnergyModule->fRunpars->PixelCorner(*itPmt,iCorner);
	p[3*iCorner]   =  x[iCorner] =  vPmt.X();
	p[3*iCorner+1] =  y[iCorner] =  vPmt.Y();
	p[3*iCorner+2] = -vPmt.Z();
      }
      
      p[12] = x[4] = x[0]; p[13] = y[4] = y[0]; p[14] = p[2];
      
      fPixelLines[*itPmt].SetPolyLine(5,x,y);
      fPixelLines[*itPmt].SetLineColor(10);
      fPixelLines[*itPmt].SetFillColor(fPalette[0]);
    
      fPixelLines3D[*itPmt].SetPolyLine(5,p);
      fPixelLines3D[*itPmt].SetLineColor(kWhite);
    }
    fDefinePixelsView = kTRUE;
}

//_____________________________________________________________________________
 void EnergyViewer::DefineFSView() {

  if (!fCanvas || fFSViewBuild) return;

  Float_t x[5], y[5], p[15];
  Float_t PhiMin, PhiMax, PhiMin1, PhiMax1;

  PMTs_t::iterator itPmt;

  fPmtLinesFS3D.clear();

  fGXmin = fGYmin = fGZmin = PhiMin = PhiMin1 =  kHuge;
  fGXmax = fGYmax = fGZmax = PhiMax = PhiMax1 = -kHuge;

  for (itPmt = fEnergyModule->fPmts.begin(); itPmt != fEnergyModule->fPmts.end(); itPmt++) {
    TVector3 vPmt = fEnergyModule->fRunpars->GetPmtGeo(*itPmt)->GetCenter();
    Float_t Phi = vPmt.Phi();
    if (Phi < PhiMin) PhiMin = Phi;
    if (Phi > PhiMax) PhiMax = Phi;

    if (Phi < 0.) Phi += 2*Pi();

    if (Phi < PhiMin1) PhiMin1 = Phi;
    if (Phi > PhiMax1) PhiMax1 = Phi;
  }

  Int_t Case = 0;

  if (PhiMax > 0.7*Pi() && PhiMin < -0.7*Pi()) {
    Case = 1;
    PhiMin = PhiMin1; PhiMax = PhiMax1;
  }

  for(Int_t i(1); i <= fEnergyModule->fRunpars->GetNumPmts(); i++) {
    TVector3 vPmt = fEnergyModule->fRunpars->GetPmtGeo(i)->GetCenter();
    Float_t phi = vPmt.Phi();
    if (Case && vPmt.Phi() < 0) phi = vPmt.Phi() + 2*Pi();
    if (phi > PhiMin - 0.05*Pi() && phi < PhiMax + 0.05*Pi()) {
      for(Int_t iCorner(0); iCorner < 4; iCorner++) {
	TVector3 vPmt = fEnergyModule->fRunpars->PmtCorner(i,iCorner);
	
	p[3*iCorner]   =  x[iCorner] =  vPmt.X();
	p[3*iCorner+1] =  y[iCorner] =  vPmt.Y();
	p[3*iCorner+2] = -vPmt.Z();
	
	if (fGXmin >  vPmt.X()) fGXmin =  vPmt.X();
	if (fGYmin >  vPmt.Y()) fGYmin =  vPmt.Y();
	if (fGZmin > -vPmt.Z()) fGZmin = -vPmt.Z();
	if (fGXmax <  vPmt.X()) fGXmax =  vPmt.X();
	if (fGYmax <  vPmt.Y()) fGYmax =  vPmt.Y();
	if (fGZmax < -vPmt.Z()) fGZmax = -vPmt.Z();
      }
      
      p[12] = x[4] = x[0]; p[13] = y[4] = y[0]; p[14] = p[2];
      
      fPmtLinesFS3D[i].SetPolyLine(5,p);
      fPmtLinesFS3D[i].SetLineColor(kBlack);
    }
  }

  fFSViewBuild = kTRUE;
}

//_____________________________________________________________________________
 void EnergyViewer::DisplayPixels() {

  if (!fCanvas) return;
  DefinePixelsView();

  Int_t Count, PixUid;
  Pixels_t OneGtuPixels;
  map <Int_t,Int_t> PixelsSignal;
  Pixels_t::iterator itPix;
  GtuPixels_t::iterator
  map <Int_t, TPolyLine>::iterator itPmt;

  PixelsSignal.clear();
  fCanvas->Clear();

  if (gPad) {
    Float_t OffSet = 5; 
    gPad->Range(fXmin-OffSet,fYmin-OffSet,fXmax+OffSet,fYmax+OffSet);
  }

  // Draw PMTs 

  for (itPmt = fPmtLines.begin(); itPmt != fPmtLines.end(); itPmt++) {
    itPmt->second.Draw("f");
    itPmt->second.Draw();
  }

  // Integrate events over all pixels

  for (itPix = fEnergyModule->fActivePixels.begin(); itPix != fEnergyModule->fActivePixels.end(); itPix++){
    PixelsSignal[(*itPix)->GetPixelId()] += (*itPix)->GetCounts();
  }

  // Draw integral signals

  for (itPmt = fPixelLines.begin(); itPmt != fPixelLines.end(); itPmt++) {
    Count = PixelsSignal[itPmt->first];
    itPmt->second.SetFillColor(fPalette[Count]);
    itPmt->second.Draw("f");
    itPmt->second.Draw();
  }

  if (fEnergyModule->fGraphXY) fEnergyModule->fGraphXY->Draw("X");
  fEnergyModule->fHistoXY->SetStats(0);
  fEnergyModule->fHistoXY->Draw("SAME");

  // show the legend on top of the page
  
  TPaveText *cl = NULL;
  TString pavetext = Form("Integral signal over all GTUs in [%d,%d] microsecond interval",fEnergyModule->fGtuMin,fEnergyModule->fGtuMax);

  DrawLabel(cl,pavetext);

  fCanvas->Print(fPsFile);

  SafeDelete(cl);

  // Restore default color for inactive pixels

  for (itPmt = fPixelLines.begin(); itPmt != fPixelLines.end(); itPmt++) 
    itPmt->second.SetFillColor(fPalette[0]);

  if (fEnergyModule->fDebugPs >= 2) {

    // Loop for all Gtu and fill the map of all active pixels
  
    for (itGtu = fEnergyModule->fGtuMyPixels.begin(); itGtu != fEnergyModule->fGtuMyPixels.end(); itGtu++) {

      fCanvas->Clear();

      // Draw PMTs 

      for (itPmt = fPmtLines.begin(); itPmt != fPmtLines.end(); itPmt++) {
	itPmt->second.Draw("f");
	itPmt->second.Draw();
      }

      // Define colors and signals in each pixels per one Gtu

      OneGtuPixels.clear(); PixelsSignal.clear();
      OneGtuPixels = itGtu->second;

      for (itPix = OneGtuPixels.begin(); itPix != OneGtuPixels.end(); itPix++) {
	Count = (*itPix)->GetCounts();
	PixUid = (*itPix)->GetPixelId();
	fPixelLines[PixUid].SetFillColor(fPalette[Count]);
        fPixelLines[PixUid].SetLineColor(kRed);
	PixelsSignal[PixUid] = Count;
      }

      // Draw Pixels 

      TText *text = new TText();

      for (itPmt = fPixelLines.begin(); itPmt != fPixelLines.end(); itPmt++) {
	itPmt->second.Draw("f");
	itPmt->second.Draw();
	Count = PixelsSignal[itPmt->first];
	if (Count > 0) {
	  TVector3 PixelCenter = fEnergyModule->fRunpars->PixelCenter(itPmt->first);
	  Float_t x = PixelCenter.X()-0.2*fEnergyModule->fRunpars->GetPmtData().GetPadSide();
	  Float_t y = PixelCenter.Y()-0.2*fEnergyModule->fRunpars->GetPmtData().GetPadSide();
	  Float_t padMin = Min(gPad->GetWh(),gPad->GetWw());
	  text->SetTextSize(5.e-5*padMin);
	  text->SetTextColor(5);
	  text->DrawText(x,y,Form("%d",Count));
	}
      }

      if (fEnergyModule->fGraphXY) fEnergyModule->fGraphXY->Draw("X");

      TString pavetext = Form("Signal at GTU = %d",itGtu->first);
      DrawLabel(cl,pavetext);

      fCanvas->Print(fPsFile);

      SafeDelete(text);
      SafeDelete(cl);

      // Restore default color for inactive pixels

      for (itPmt = fPixelLines.begin(); itPmt != fPixelLines.end(); itPmt++) {
	itPmt->second.SetFillColor(fPalette[0]);
	itPmt->second.SetLineColor(kWhite);
      }
    }
  }
}

//_____________________________________________________________________________
 void EnergyViewer::DisplayPixels3D(TString DrawSurf) {

    if (!fCanvas) return;
    Int_t Count;
    Float_t Xmin, Xmax, Ymin, Ymax, Zmin, Zmax;
    map <Int_t,Int_t> PixelsSignal;
    map<Int_t,TPolyLine3D> Pmts3D;
    Pixels_t::iterator itPix;
    map <Int_t, TPolyLine3D>::iterator itPmt;

    fCanvas->Clear();
    PixelsSignal.clear();
    Pmts3D.clear();
    
    if (DrawSurf == "FocalSurf") {
      DefineFSView();
      Xmin = fGXmin; Xmax = fGXmax; Ymin = fGYmin; Ymax = fGYmax;
      Zmin = fGZmin; Zmax = fGZmax;
      Pmts3D = fPmtLinesFS3D;
    }
    else if (DrawSurf == "ClusterSurf") {
      DefinePixelsView();
      Xmin = fXmin; Xmax = fXmax; Ymin = fYmin; Ymax = fYmax;
      Zmin = fZmin; Zmax = fZmax;
      Pmts3D = fPmtLines3D;
    }

    TView *view = new TView(1);
    view->ShowAxis();
    view->SetRange(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);

    Float_t xm = 0.5*(Xmax+Xmin);
    Float_t ym = 0.5*(Ymax+Ymin);

    if (xm > 0 && ym > 0) view->RotateView(-135,60);
    else if (xm < 0 && ym < 0) view->RotateView(45,60);
    else if (xm < 0 && ym > 0) view->RotateView(-45,60);
    else if (xm > 0 && ym < 0) view->RotateView(135,60);

    // Draw PMTs

    for (itPmt = Pmts3D.begin(); itPmt != Pmts3D.end(); itPmt++)
      itPmt->second.Draw();

    // Calculate number of counts for each pixel

    for (itPix = fEnergyModule->fActivePixels.begin(); 
 	 itPix != fEnergyModule->fActivePixels.end(); itPix++)
      PixelsSignal[(*itPix)->GetPixelId()] += (*itPix)->GetCounts();

    // Draw integral signals

    for (itPmt = fPixelLines3D.begin(); itPmt != fPixelLines3D.end(); itPmt++) {
      Count = PixelsSignal[itPmt->first];
      if (Count) {
 	itPmt->second.SetLineColor(fPalette[Count]);
	itPmt->second.Draw();
      }
    }

    // show the legend on top of the page

    TPaveText *cl = NULL;
    TString pavetext = Form("Integral signal over all GTUs in [%d,%d] microsecond interval",fEnergyModule->fGtuMin,fEnergyModule->fGtuMax);
    
    DrawLabel(cl,pavetext);

    fCanvas->Print(fPsFile);

    // Restore default color for inactive pixels

    for (itPmt = fPixelLines3D.begin(); itPmt != fPixelLines3D.end(); itPmt++)
      itPmt->second.SetLineColor(kWhite);

    SafeDelete(cl);
    SafeDelete(view);
}

//_____________________________________________________________________________
 void EnergyViewer::DisplayCherenkov(TH1F *htime, Int_t nPoints, Float_t x[], Float_t y[]) {

  if (!fCanvas) return;
  fCanvas->Clear();

  TList *functions = htime->GetListOfFunctions();                        
  TPolyMarker *pm = (TPolyMarker*)functions->FindObject("TPolyMarker");

  htime->Draw("histo");

  TGraph *gr1 = new TGraph(nPoints,x,y);
  gr1->SetMarkerColor(kBlue);
  gr1->SetMarkerStyle(21);
  gr1->Draw("P");
  if (pm) pm->Draw();
      
  // fCanvas->Print(fPsFileTime+"(");
  fCanvas->Print(fPsFile);
  SafeDelete(gr1);
}


//_____________________________________________________________________________
 void EnergyViewer::DisplayTimeFit(TH1F *htime, TF1 *func, TF1 *bg)
{
    if (!fCanvas) return;
    fCanvas->Clear();

    htime->SetStats(0);
    htime->Draw("PE");
    func->Draw("same");
    bg->Draw("same");

    TPaveText *cl = new TPaveText(0.7,0.85,1.2,1.,"brNDC");
    cl->SetTextFont(72); cl->SetTextSize(0.02);
    cl->SetFillColor(0); cl->SetTextAlign(12);

    // Event Information

    cl->AddText("Truth Information");
    cl->AddText(Form("Energy = %.2E MeV",
		     fEnergyModule->fEvent->GetHeader().GetTrueEnergy()));
    cl->AddText(Form("#Theta = %.2f deg",
		     fEnergyModule->fEvent->GetHeader().GetTrueTheta()*RadToDeg()));
    cl->AddText(Form("#phi   = %.2f deg",
		     fEnergyModule->fEvent->GetHeader().GetTruePhi()*RadToDeg()));
    cl->AddText("Reco Information");
    if (!IsNaN(fEnergyModule->fEnergy)) 
      cl->AddText(Form("Energy = %.2E #pm %2.E MeV",
		       fEnergyModule->fEnergy,fEnergyModule->fEnergyError));
    if (!IsNaN(fEnergyModule->fTheta)) 
      cl->AddText(Form("#Theta = %.2f deg",fEnergyModule->fTheta*RadToDeg()));
    if (!IsNaN(fEnergyModule->fPhi))
      cl->AddText(Form("#phi   = %.2f deg",fEnergyModule->fPhi*RadToDeg()));
    cl->Draw();

    fCanvas->Print(fPsFile);

    SafeDelete(cl);

}


//_____________________________________________________________________________
 void EnergyViewer::DrawHisto(TH1F *histo)
{    
  if (!fCanvas) return;
  fCanvas->Clear();

  histo->Draw();
  fCanvas->Print(fPsFile);
}


//_____________________________________________________________________________
 void EnergyViewer::Clean() {

    fMyPixelLines.clear();
    fPixelLines.clear(); fPixelLines3D.clear();
    fPmtLines.clear();   fPmtLines3D.clear();
    fPmtLinesFS3D.clear();

    fFSViewBuild      = kFALSE; 
    fDefinePixelsView = kFALSE;
}
 
 void EnergyViewer::SetHistoStyle(TH2F *histo)
{
  histo->GetXaxis()->SetLabelFont(63);
  histo->GetXaxis()->SetLabelSize(15);
  histo->GetYaxis()->SetLabelFont(63);
  histo->GetYaxis()->SetLabelSize(15);
  
  histo->GetXaxis()->SetTitleFont(73);
  histo->GetXaxis()->SetTitleSize(15);
  histo->GetYaxis()->SetTitleFont(73);
  histo->GetYaxis()->SetTitleSize(15);
  
  histo->SetTitleOffset(1.3,"X");
  histo->SetTitleOffset(1.3,"Y");

  histo->GetXaxis()->SetNdivisions(507,kTRUE);
  histo->GetYaxis()->SetNdivisions(505,kTRUE);

  histo->SetMarkerStyle(kFullDotMedium);
  histo->SetMarkerColor(kBlue);
}


 void EnergyViewer::DrawLabel(TPaveText *cl, TString pavetext)
{
  cl = new TPaveText(0.1,0.96,0.9,1.,"brNDC");
  cl->SetTextFont(72);
  cl->SetTextSize(0.03);
  cl->SetFillColor(kGreen);
  cl->SetTextAlign(22);
  cl->SetTextColor(kBlue);
  cl->AddText(pavetext);
  cl->Draw();
}
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