HarpoAnalyseEtrack.cxx

Go to the documentation of this file.

// 
// File HarpoAnalyseEtrack.cxx
//
#include "HarpoAnalyseEtrack.h"
#include "HarpoConfig.h"
#include "HarpoDetSet.h"
#include "HarpoDebug.h"
#include "HarpoDccEvent.h"
#include "Pmm2Event.h"
#include "HarpoEvent.h"
#include "HarpoRecoEvent.h"
#include "HarpoSimEvent.h"
#include "TpcSimMCEvent.h"
#include "MakeNiceHisto.h"

#include "TFile.h"
#include "TStyle.h"
#include "TCanvas.h"
#include "TRootEmbeddedCanvas.h"
#include "TLatex.h"
#include "TGraph.h"
#include "TF1.h"
#include "TMath.h"
#include "TGLabel.h"
#include "TSystem.h"

#include <cstdlib>
#include <cstring>
#include <cassert>
#include <fstream>
#include <iostream>

ClassImp(HarpoAnalyseEtrack);

void   HarpoAnalyseEtrack::print()
 {
   unsigned int nd; // number of detectors
   HarpoEventHeader *hdr = fEvt->GetHeader();

   assert(hdr != NULL);
   hdr->print();
   
   for (nd = 0; nd < gkNDetectors; nd++) {
     //     if (fEvt->isdataExist(nd)) {
       HarpoDccMap *plane = fEvt->GetDccMap(nd);
       if (plane != NULL )plane->print();
   }
 }

void   HarpoAnalyseEtrack::process()
{
  //  Bool_t drawEvent = kFALSE;
  nEvents++;
  Info("PROCESS","START %ld",nEvents);
  //???????????????????????????????
  if(fEvt->GetRecoEvent()){  // Reconstructed data (clusters, tracks, matched tracks)
    HarpoRecoEvent* reco = fEvt->GetRecoEvent();
    // Clusters in the event (both X and Y directions)
    TClonesArray* trArray = reco->GetKalmanTracksArray();
    Int_t nTr = trArray->GetEntries();
    // 2D Vertexes in the event (both X and Y directions)
    TClonesArray* vArray = reco->GetVertexArray();
    //Int_t nV = vArray->GetEntries();
    TClonesArray* v3DArray = reco->GetVertex3DArray();
    Int_t nV3D = v3DArray->GetEntries();
    //std::cout << nV3D << " " << nV << " " << nTr << std::endl;
    for(Int_t iV3D = 0; iV3D<nV3D; iV3D++){
      HarpoRecoVertex3D* vertex3D = (HarpoRecoVertex3D*)v3DArray->At(iV3D);
      Int_t iVx = vertex3D->GetVertexX();
      Int_t iVy = vertex3D->GetVertexX();
      HarpoRecoVertex* vertexX = (HarpoRecoVertex*)vArray->At(iVx);
      HarpoRecoVertex* vertexY = (HarpoRecoVertex*)vArray->At(iVy);
      Int_t idTrack[2][2];
      idTrack[XPLANE][0] = vertexX->GetTr1Vertex1();
      idTrack[XPLANE][1] = vertexX->GetTr2Vertex2();
      idTrack[YPLANE][0] = vertexY->GetTr1Vertex1();
      idTrack[YPLANE][1] = vertexY->GetTr2Vertex2();
      Double_t chi2[8][2][2];
      Int_t ncl[2][2];
      for(Int_t plane = 0; plane<2; plane++){
        for(Int_t i = 0; i<2; i++){
          for(Int_t j = 0; j<8; j++)
            chi2[j][plane][i] = -1;
          ncl[plane][i] = -1;
        }
      }
      // for(Int_t plane = 0; plane<2; plane++){
      //        for(Int_t i = 0; i<2; i++){
      //          Info("process","%d  %d    %d",plane,i,idTrack[plane][i]);
      //        }
      // }
      for(Int_t itr = 0; itr<nTr; itr++){
        HarpoRecoKalmanTracks* track = (HarpoRecoKalmanTracks*)trArray->At(itr);
        for(Int_t plane = 0; plane<2; plane++){
          for(Int_t i = 0; i<2; i++){
            if(track->GetPlane()==plane && track->GetNtrack()==idTrack[plane][i]){
              ncl[plane][i] = track->GetNclusterTrack();
              chi2[0][plane][i] = track->GetXstart();
              chi2[1][plane][i] = track->GetZstart();
              chi2[2][plane][i] = track->GetPxStart();
              chi2[3][plane][i] = track->GetPzStart();
              chi2[4][plane][i] = track->GetXend();
              chi2[5][plane][i] = track->GetZend();
              chi2[6][plane][i] = track->GetPxEnd();
              chi2[7][plane][i] = track->GetPzEnd();

              //Info("process","%d  %d    %g",plane,idTrack[plane][i],chi2[plane][i]);
            }
          }
        }
      }

      hChi2_XY->Fill(chi2[0][0][0],chi2[0][1][0]);
      hChi2_XY->Fill(chi2[0][0][1],chi2[0][1][1]);
      hChi2_3D->Fill(0.5*(chi2[0][0][1]+chi2[0][1][1]));
      hChi2_3D->Fill(0.5*(chi2[0][0][0]+chi2[0][1][0]));
      for(Int_t i = 0; i<2; i++){
        Double_t data[18];
        for(Int_t j = 0; j<18; j++) data[j] = -2;
        data[0] = ncl[0][i];
        data[1] = ncl[1][i];
        for(Int_t j = 0; j<8; j++){
          //if(ncl[0][i]<j+4 || ncl[1][i] < j+4) break;
          data[2*j+2] = chi2[j][0][i];
          data[2*j+3] = chi2[j][1][i];
        }
        std::cerr << nEvents << "   ";
        for(Int_t j = 0; j<18; j++)
          std::cerr << data[j] << " ";
        std::cerr << std::endl;
        fNtuple->Fill(data);

      }
    }
  }




   // if (gHDetSet->isExist(SIMDET)){
   //   HarpoSimEvent *simEvt = (HarpoSimEvent *)(fEvt->GetDetEvent(SIMDET));
   //   Int_t nIonTr = simEvt->GetNtracks();
   //   for(Int_t i = 0; i<nIonTr; i++){
   //     TpcSimIonisationTrack* tr = simEvt->GetIonisationTrack(i);
   //     Int_t nPoints = tr->GetNpoints();
   //     for(Int_t j = 0; j<nPoints; j++){
   //    TpcSimIonisationPoint* point = tr->GetPoint(j);
   //    hNeSim->Fill(point->GetNe());
   //     }
   //   }
   //   TpcSimMCEvent* mcEvent = simEvt->GetMCEvent();
   //   Int_t nMCtr = mcEvent->GetNtracks();
   //   data[4] = nMCtr;
   //   for(Int_t i = 0; i<nMCtr; i++){
   //     TpcSimMCTrack* track = mcEvent->GetTrack(i);
   //     hPsim->Fill(track->GetPx(),track->GetPy(),track->GetPz());
   //   }
   // }


}

void   HarpoAnalyseEtrack::Init()
{
  
  fChooseNbins = 0;

  // Initialise histograms here

  hChi2_XY = new TH2F("hChi2_XY","",500,0,TMath::Pi(),500,0,TMath::Pi());
  hChi2_3D = new TH1F("hChi2_3D","",500,0,TMath::Pi());

  //  fNtuple = new TNtupleD("fNtuple","","chi2X:chi2Y:nClX:nClY");
  fNtuple = new TNtupleD("fNtuple","","nClX:nClY:chi2X:chi2Y:chi2X2:chi2Y2:chi2X3:chi2Y3:chi2X4:chi2Y4:chi2X5:chi2Y5:chi2X6:chi2Y6:chi2X7:chi2Y7:chi2X8:chi2Y8");

 }

 void   HarpoAnalyseEtrack::Save(char * /* mode */)
 {

   OpenHistFile("template");

   // Save histograms here
   hChi2_XY->Write();
   hChi2_3D->Write();
   fNtuple->Write();
 }



void HarpoAnalyseEtrack::DisplayAnalysis(TRootEmbeddedCanvas* ecTab, TGListBox* /* infobox */)
{
  // in hrecomonitorgui: fills the canvas on the "Analysis" tab
  TCanvas* c = ecTab->GetCanvas();
  c->GetPad(1)->Delete();
  c->GetPad(2)->Delete();
  c->Divide(2);

  for(Int_t plane = 0; plane<2;plane++){
    TVirtualPad* cMap = c->GetPad(plane+1);
    HarpoDccMap* m = fEvt->GetDccMap(plane);
    TH2F* h = new TH2F("h","",512,0,512,288,0,288);
    for(Int_t i=1;i<NADC;i++){ //Time bins
      for(Int_t j=0;j<NCHAN;j++){ // Channels
        Double_t q = m->GetData(j,i);
        h->SetBinContent(i+1,j+1,q);
      }
    } 
    h->SetMinimum(0);
    MakeNiceHisto(h,cMap,"colz");
    h->Delete();
  }

  c->Update();
}


void HarpoAnalyseEtrack::ConfigFrame(TGMainFrame* fMain, Int_t id)
{
  // in hrecomonitorgui: Creates a popup window for analysis configuration
  // You must define SetConfig() properly

  UInt_t xsize = 200;
  UInt_t ysize2 = 20;
  UInt_t ysize = 10*ysize2;
  TGTransientFrame* main = new TGTransientFrame(gClient->GetRoot(), fMain, xsize, ysize);
  main->Connect("CloseWindow()", "HarpoRecoMonitorGui", main, "CloseWindow()");
  main->DontCallClose(); // to avoid double deletions.
  
  // use hierarchical cleaning
  main->SetCleanup(kDeepCleanup);

  TGVerticalFrame* frame = new TGVerticalFrame(main,xsize,ysize,kVerticalFrame);

  // Object layout options
  TGLayoutHints* fLayout1 = new TGLayoutHints(kLHintsLeft ,5,5,5,5);
  TGLayoutHints* fLayout2 = new TGLayoutHints(kLHintsRight ,5,5,5,5);
  TGLayoutHints* fLayout3 = new TGLayoutHints(kLHintsTop | kLHintsExpandX ,5,5,5,5);
  TGLayoutHints* fLayout4 = new TGLayoutHints(kLHintsExpandX | kLHintsExpandY,5,5,5,5);

  //________ DO NOT MODIFY ABOVE _____________________________




  // Title of the analysis
  TGLabel* fAnalysisLabel = new TGLabel(frame, "Template analysis");

  // Create a line for choosing value of parameter
  TGHorizontalFrame* nBinsFrame = new TGHorizontalFrame(frame,xsize,ysize2,kHorizontalFrame);
  TGLabel* nBinsLabel = new TGLabel(nBinsFrame, "fNbins");
  fChooseNbins    = new TGNumberEntry(nBinsFrame);
  //  fChooseNbins->SetNumber(fNbins);

  main->AddFrame(frame,fLayout4);
  frame->AddFrame(fAnalysisLabel,fLayout3);
  frame->AddFrame(nBinsFrame,fLayout3);
  nBinsFrame->AddFrame(nBinsLabel,fLayout1);
  nBinsFrame->AddFrame(fChooseNbins,fLayout2);



  //________ DO NOT MODIFY BELOW _____________________________
  // Button to validate configuration
  TGTextButton* setConf = new TGTextButton(frame,"Save Config",id);
  setConf->Associate(fMain);

  frame->AddFrame(setConf,fLayout3);

  main->MapSubwindows();
  main->MapWindow();
  main->Resize();
  return;
}



void   HarpoAnalyseEtrack::SetConfig()
{
  // Update the configuration according to the values in the popup window

  if(!fChooseNbins) return;

  //  fNbins = fChooseNbins->GetNumber();
}