ROOT logo
/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */

#include "AlidNdEtaCorrection.h"

#include <AliLog.h>
#include <TCanvas.h>
#include <TH3F.h>
#include <TH2F.h>
#include <TH1D.h>
#include <TDirectory.h>
#include <AliCorrection.h>
#include <AliCorrectionMatrix2D.h>
#include <AliCorrectionMatrix3D.h>

//____________________________________________________________________
ClassImp(AlidNdEtaCorrection)

//____________________________________________________________________
AlidNdEtaCorrection::AlidNdEtaCorrection()
  : TNamed(),
  fTrack2ParticleCorrection(0),
  fVertexRecoCorrection(0),
  fTriggerBiasCorrectionMBToINEL(0),
  fTriggerBiasCorrectionMBToNSD(0),
  fTriggerBiasCorrectionMBToND(0),
  fTriggerBiasCorrectionMBToOnePart(0)
{
  // default constructor
}

//____________________________________________________________________
AlidNdEtaCorrection::AlidNdEtaCorrection(const Char_t* name, const Char_t* title, AliPWG0Helper::AnalysisMode analysis)
  : TNamed(name, title),
  fTrack2ParticleCorrection(0),
  fVertexRecoCorrection(0),
  fTriggerBiasCorrectionMBToINEL(0),
  fTriggerBiasCorrectionMBToNSD(0),
  fTriggerBiasCorrectionMBToND(0),
  fTriggerBiasCorrectionMBToOnePart(0)
{
  //
  // constructor
  //

  fTrack2ParticleCorrection = new AliCorrection("Track2Particle", "Track2Particle", analysis);
  fVertexRecoCorrection     = new AliCorrection("VertexReconstruction", "VertexReconstruction", analysis);

  fTriggerBiasCorrectionMBToINEL = new AliCorrection("TriggerBias_MBToINEL", "TriggerBias_MBToINEL", analysis);
  fTriggerBiasCorrectionMBToNSD  = new AliCorrection("TriggerBias_MBToNSD", "TriggerBias_MBToNSD", analysis);
  fTriggerBiasCorrectionMBToND   = new AliCorrection("TriggerBias_MBToND", "TriggerBias_MBToND", analysis);
  fTriggerBiasCorrectionMBToOnePart   = new AliCorrection("TriggerBias_MBToOnePart", "TriggerBias_MBToOnePart", analysis);
}

//____________________________________________________________________
AlidNdEtaCorrection::~AlidNdEtaCorrection()
{
  // destructor

  if (fTrack2ParticleCorrection) {
    delete fTrack2ParticleCorrection;
    fTrack2ParticleCorrection = 0;
  }

  if (fVertexRecoCorrection) {
    delete fVertexRecoCorrection;
    fVertexRecoCorrection = 0;
  }

  if (fTriggerBiasCorrectionMBToINEL) {
    delete fTriggerBiasCorrectionMBToINEL;
    fTriggerBiasCorrectionMBToINEL = 0;
  }

  if (fTriggerBiasCorrectionMBToNSD) {
    delete fTriggerBiasCorrectionMBToNSD;
    fTriggerBiasCorrectionMBToNSD = 0;
  }

  if (fTriggerBiasCorrectionMBToND) {
    delete fTriggerBiasCorrectionMBToND;
    fTriggerBiasCorrectionMBToND = 0;
  }

  if (fTriggerBiasCorrectionMBToOnePart) {
    delete fTriggerBiasCorrectionMBToOnePart;
    fTriggerBiasCorrectionMBToOnePart = 0;
  }
}

//____________________________________________________________________
void
AlidNdEtaCorrection::Finish() {
  //
  // finish method
  //
  // divide the histograms in the AliCorrectionMatrix2D objects to get the corrections

  fTrack2ParticleCorrection->Divide();
  fVertexRecoCorrection->Divide();
  fTriggerBiasCorrectionMBToINEL->Divide();
  fTriggerBiasCorrectionMBToNSD->Divide();
  fTriggerBiasCorrectionMBToND->Divide();
  fTriggerBiasCorrectionMBToOnePart->Divide();
}

//____________________________________________________________________
Long64_t AlidNdEtaCorrection::Merge(TCollection* list)
{
  // Merge a list of dNdEtaCorrection objects with this (needed for
  // PROOF).
  // Returns the number of merged objects (including this).

  if (!list)
    return 0;

  if (list->IsEmpty())
    return 1;

  TIterator* iter = list->MakeIterator();
  TObject* obj;

  // collections of measured and generated histograms
  TList* collectionNtrackToNparticle    = new TList;
  TList* collectionVertexReco           = new TList;
  TList* collectionTriggerBiasMBToINEL  = new TList;
  TList* collectionTriggerBiasMBToNSD   = new TList;
  TList* collectionTriggerBiasMBToND    = new TList;
  TList* collectionTriggerBiasMBToOnePart    = new TList;

  Int_t count = 0;
  while ((obj = iter->Next())) {

    AlidNdEtaCorrection* entry = dynamic_cast<AlidNdEtaCorrection*> (obj);
    if (entry == 0)
      continue;

    collectionNtrackToNparticle  ->Add(entry->fTrack2ParticleCorrection);
    collectionVertexReco         ->Add(entry->fVertexRecoCorrection);
    collectionTriggerBiasMBToINEL->Add(entry->fTriggerBiasCorrectionMBToINEL);
    collectionTriggerBiasMBToNSD ->Add(entry->fTriggerBiasCorrectionMBToNSD);
    collectionTriggerBiasMBToND  ->Add(entry->fTriggerBiasCorrectionMBToND);
    collectionTriggerBiasMBToOnePart  ->Add(entry->fTriggerBiasCorrectionMBToOnePart);

    count++;
  }
  fTrack2ParticleCorrection      ->Merge(collectionNtrackToNparticle);
  fVertexRecoCorrection          ->Merge(collectionVertexReco);
  fTriggerBiasCorrectionMBToINEL ->Merge(collectionTriggerBiasMBToINEL);
  fTriggerBiasCorrectionMBToNSD  ->Merge(collectionTriggerBiasMBToNSD);
  fTriggerBiasCorrectionMBToND   ->Merge(collectionTriggerBiasMBToND);
  fTriggerBiasCorrectionMBToOnePart->Merge(collectionTriggerBiasMBToOnePart);

  delete collectionNtrackToNparticle;
  delete collectionVertexReco;
  delete collectionTriggerBiasMBToINEL;
  delete collectionTriggerBiasMBToNSD;
  delete collectionTriggerBiasMBToND;
  delete collectionTriggerBiasMBToOnePart;

  return count+1;
}

//____________________________________________________________________
void AlidNdEtaCorrection::Add(AlidNdEtaCorrection* aCorrectionsToAdd, Float_t c) {
  //
  // adds the measured and generated of aCorrectionsToAdd to measured and generated
  // of all corrections in this

  fTrack2ParticleCorrection      ->Add(aCorrectionsToAdd->GetTrack2ParticleCorrection() ,c);
  fVertexRecoCorrection          ->Add(aCorrectionsToAdd->GetVertexRecoCorrection()     ,c);
  fTriggerBiasCorrectionMBToINEL ->Add(aCorrectionsToAdd->GetTriggerBiasCorrectionINEL(),c);
  fTriggerBiasCorrectionMBToNSD  ->Add(aCorrectionsToAdd->GetTriggerBiasCorrectionNSD() ,c);
  fTriggerBiasCorrectionMBToND   ->Add(aCorrectionsToAdd->GetTriggerBiasCorrectionND()  ,c);
  fTriggerBiasCorrectionMBToOnePart   ->Add(aCorrectionsToAdd->GetTriggerBiasCorrectionOnePart()  ,c);
}

//____________________________________________________________________
void AlidNdEtaCorrection::Scale(Float_t c) 
{
  //
  // scales all contained corrections
  // 

  fTrack2ParticleCorrection      ->Scale(c);
  fVertexRecoCorrection          ->Scale(c);
  fTriggerBiasCorrectionMBToINEL ->Scale(c);
  fTriggerBiasCorrectionMBToNSD  ->Scale(c);
  fTriggerBiasCorrectionMBToND   ->Scale(c);
  fTriggerBiasCorrectionMBToOnePart   ->Scale(c);
}

//____________________________________________________________________
void AlidNdEtaCorrection::Reset(void) {
  //
  // reset all corrections
  // 

  fTrack2ParticleCorrection      ->Reset();
  fVertexRecoCorrection          ->Reset();
  fTriggerBiasCorrectionMBToINEL ->Reset();
  fTriggerBiasCorrectionMBToNSD  ->Reset();
  fTriggerBiasCorrectionMBToND   ->Reset();
  fTriggerBiasCorrectionMBToOnePart   ->Reset();
}



//____________________________________________________________________
Bool_t AlidNdEtaCorrection::LoadHistograms(const Char_t* dir)
{
  //
  // loads the histograms
  // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
  //

  if (!dir)
    dir = GetName();

  if (!gDirectory->cd(dir))
    return kFALSE;

  fTrack2ParticleCorrection      ->LoadHistograms();
  fVertexRecoCorrection          ->LoadHistograms();
  fTriggerBiasCorrectionMBToINEL ->LoadHistograms();
  fTriggerBiasCorrectionMBToNSD  ->LoadHistograms();
  fTriggerBiasCorrectionMBToND   ->LoadHistograms();
  fTriggerBiasCorrectionMBToOnePart   ->LoadHistograms();

  gDirectory->cd("..");

  return kTRUE;
}

//____________________________________________________________________
void AlidNdEtaCorrection::SaveHistograms()
{
  //
  // save the histograms
  //

  gDirectory->mkdir(fName.Data());
  gDirectory->cd(fName.Data());

  fTrack2ParticleCorrection     ->SaveHistograms();
  fVertexRecoCorrection         ->SaveHistograms();
  fTriggerBiasCorrectionMBToINEL->SaveHistograms();
  fTriggerBiasCorrectionMBToNSD ->SaveHistograms();
  fTriggerBiasCorrectionMBToND  ->SaveHistograms();
  fTriggerBiasCorrectionMBToOnePart->SaveHistograms();

  gDirectory->cd("..");
}

//____________________________________________________________________
void AlidNdEtaCorrection::DrawHistograms()
{
  //
  // call the draw histogram method of the corrections
  //

  fTrack2ParticleCorrection     ->DrawHistograms();
  fVertexRecoCorrection         ->DrawHistograms();
  fTriggerBiasCorrectionMBToINEL->DrawHistograms();
  fTriggerBiasCorrectionMBToNSD ->DrawHistograms();
  fTriggerBiasCorrectionMBToND  ->DrawHistograms();
  fTriggerBiasCorrectionMBToOnePart  ->DrawHistograms();
}

//____________________________________________________________________
void AlidNdEtaCorrection::DrawOverview(const char* canvasName)
{
  //
  // call the DrawOverview histogram method of the corrections
  //

  fTrack2ParticleCorrection     ->DrawOverview(canvasName);
  fVertexRecoCorrection         ->DrawOverview(canvasName);
  fTriggerBiasCorrectionMBToINEL->DrawOverview(canvasName);
  fTriggerBiasCorrectionMBToNSD ->DrawOverview(canvasName);
  fTriggerBiasCorrectionMBToND  ->DrawOverview(canvasName);
  fTriggerBiasCorrectionMBToOnePart  ->DrawOverview(canvasName);
}

//____________________________________________________________________
void AlidNdEtaCorrection::FillMCParticle(Float_t vtx, Float_t eta, Float_t pt, Bool_t trigger, Bool_t vertex, Int_t processType)
{
  // fills a particle in the corrections
  // it is filled in generated or measured depending of the flags

  fTriggerBiasCorrectionMBToINEL->GetTrackCorrection()->FillGene(vtx, eta, pt);

  if ((processType & AliPWG0Helper::kSD) == 0)
    fTriggerBiasCorrectionMBToNSD->GetTrackCorrection()->FillGene(vtx, eta, pt);

  if (processType & AliPWG0Helper::kND )
    fTriggerBiasCorrectionMBToND->GetTrackCorrection()->FillGene(vtx, eta, pt);

  if (processType & AliPWG0Helper::kOnePart)
    fTriggerBiasCorrectionMBToOnePart->GetTrackCorrection()->FillGene(vtx, eta, pt);

  if (!trigger)
    return;

  fTriggerBiasCorrectionMBToINEL->GetTrackCorrection()->FillMeas(vtx, eta, pt);
  fTriggerBiasCorrectionMBToNSD->GetTrackCorrection()->FillMeas(vtx, eta, pt);
  fTriggerBiasCorrectionMBToND->GetTrackCorrection()->FillMeas(vtx, eta, pt);
  fTriggerBiasCorrectionMBToOnePart->GetTrackCorrection()->FillMeas(vtx, eta, pt);
  fVertexRecoCorrection->GetTrackCorrection()->FillGene(vtx, eta, pt);

  if (!vertex)
    return;

  fVertexRecoCorrection->GetTrackCorrection()->FillMeas(vtx, eta, pt);
  fTrack2ParticleCorrection->GetTrackCorrection()->FillGene(vtx, eta, pt);
}

//____________________________________________________________________
void AlidNdEtaCorrection::FillTrackedParticle(Float_t vtx, Float_t eta, Float_t pt, Double_t weight)
{
  // fills a tracked particle in the corrections

  fTrack2ParticleCorrection->GetTrackCorrection()->FillMeas(vtx, eta, pt, weight);
}

//____________________________________________________________________
void AlidNdEtaCorrection::FillEvent(Float_t vtx, Float_t n, Bool_t trigger, Bool_t vertex, Int_t processType)
{
  // fills an event int he correction
  // it is filled in generated or measured depending of the flags

  fTriggerBiasCorrectionMBToINEL->GetEventCorrection()->FillGene(vtx, n);

  if ((processType & AliPWG0Helper::kSD) == 0)
    fTriggerBiasCorrectionMBToNSD->GetEventCorrection()->FillGene(vtx, n);

  if (processType & AliPWG0Helper::kND )
    fTriggerBiasCorrectionMBToND->GetEventCorrection()->FillGene(vtx, n);

  if (processType & AliPWG0Helper::kOnePart)
    fTriggerBiasCorrectionMBToOnePart->GetEventCorrection()->FillGene(vtx, n);

  if (!trigger)
    return;

  fTriggerBiasCorrectionMBToINEL->GetEventCorrection()->FillMeas(vtx, n);
  fTriggerBiasCorrectionMBToNSD->GetEventCorrection()->FillMeas(vtx, n);
  fTriggerBiasCorrectionMBToND->GetEventCorrection()->FillMeas(vtx, n);
  fTriggerBiasCorrectionMBToOnePart->GetEventCorrection()->FillMeas(vtx, n);
  fVertexRecoCorrection->GetEventCorrection()->FillGene(vtx, n);

  if (!vertex)
    return;

  fVertexRecoCorrection->GetEventCorrection()->FillMeas(vtx, n);
}

//____________________________________________________________________
Float_t AlidNdEtaCorrection::GetMeasuredFraction(CorrectionType correctionType, Float_t ptCutOff, Float_t eta, Int_t vertexBegin, Int_t vertexEnd, Bool_t debug)
{
  // calculates the fraction of particles measured (some are missed due to the pt cut off)
  //
  // uses the generated particle histogram from the correction passed, e.g. pass GetTrack2ParticleCorrection()

  if (!GetCorrection(correctionType))
    return -1;

  const TH3* generated = GetCorrection(correctionType)->GetTrackCorrection()->GetGeneratedHistogram();

  // find eta borders, if eta is negative assume -0.8 ... 0.8
  Int_t etaBegin = 0;
  Int_t etaEnd = 0;
  const TAxis * xax = generated->GetXaxis();
  const TAxis * yax = generated->GetYaxis();
  if (eta < -99)
  {
    etaBegin = yax->FindFixBin(-0.8);
    etaEnd = yax->FindFixBin(0.8);
  }
  else
  {
    etaBegin = yax->FindFixBin(eta);
    etaEnd = etaBegin;
  }

  if (vertexBegin == -1)
    vertexBegin = xax->FindFixBin(-9.99);

  if (vertexEnd == -1)
    vertexEnd = xax->FindFixBin(9.99);

  TH1D* ptProj = dynamic_cast<TH1D*> (generated->ProjectionZ(Form("%s_pt", generated->GetName()), vertexBegin, vertexEnd, etaBegin, etaEnd));
  //printf("GetMeasuredFraction: bin range %d %d %d %d\n", vertexBegin, vertexEnd, etaBegin, etaEnd);
  ptProj->GetXaxis()->SetTitle(generated->GetZaxis()->GetTitle());

  Int_t ptBin = ptProj->FindBin(ptCutOff);
  //printf("GetMeasuredFraction: bin range %d %d\n", ptBin, ptProj->GetNbinsX());
  Float_t abovePtCut = ptProj->Integral(ptBin, ptProj->GetNbinsX()+1);
  Float_t all = ptProj->Integral(1, ptProj->GetNbinsX()+1);

  if (all == 0)
    return -1;

  Float_t fraction = abovePtCut / all;

  //printf("GetMeasuredFraction: all %f above %f fraction %f\n", all, abovePtCut, fraction);

  if (debug)
  {
    new TCanvas;
    ptProj->Draw();
  }
  else
    delete ptProj;

  if (debug)
    printf("AlidNdEtaCorrection::GetMeasuredFraction: pt cut off = %f, eta = %f, => fraction = %f\n", ptCutOff, eta, fraction);

  return fraction;
}

//____________________________________________________________________
TH1* AlidNdEtaCorrection::GetMeasuredEventFraction(CorrectionType correctionType, Int_t multCut)
{
  // calculates the fraction of events above multCut (but including it)
  //
  // uses the generated event histogram from the correction passed, e.g. pass GetTrack2ParticleCorrection()

  if (!GetCorrection(correctionType))
    return 0;

  const TH2* generated = GetCorrection(correctionType)->GetEventCorrection()->GetGeneratedHistogram();

  TH1* allEvents = generated->ProjectionX(Form("%s_all", generated->GetName()), 1, generated->GetNbinsY());
  TH1* aboveEvents = generated->ProjectionX(Form("%s_above", generated->GetName()), generated->GetYaxis()->FindFixBin(multCut), generated->GetNbinsY());
  
  aboveEvents->Divide(aboveEvents, allEvents, 1, 1, "B");

  return aboveEvents;  
}

//____________________________________________________________________
void AlidNdEtaCorrection::ReduceInformation()
{
  // this function deletes the measured and generated histograms from the corrections to reduce the amount of data
  // in memory

  // these are needed for GetMeasuredFraction(): fTrack2ParticleCorrection->ReduceInformation();
  fVertexRecoCorrection          ->ReduceInformation();
  fTriggerBiasCorrectionMBToINEL ->ReduceInformation();
  fTriggerBiasCorrectionMBToNSD  ->ReduceInformation();
  fTriggerBiasCorrectionMBToND   ->ReduceInformation();
  fTriggerBiasCorrectionMBToOnePart   ->ReduceInformation();
}

//____________________________________________________________________
AliCorrection* AlidNdEtaCorrection::GetCorrection(CorrectionType correctionType)
{
  // returns the given correction

  switch (correctionType)
  {
    case kNone : return 0;
    case kTrack2Particle : return fTrack2ParticleCorrection;
    case kVertexReco :     return fVertexRecoCorrection;
    case kINEL :           return fTriggerBiasCorrectionMBToINEL;
    case kNSD :            return fTriggerBiasCorrectionMBToNSD;
    case kND :             return fTriggerBiasCorrectionMBToND;
    case kOnePart :             return fTriggerBiasCorrectionMBToOnePart;
  }

  return 0;
}
 AlidNdEtaCorrection.cxx:1
 AlidNdEtaCorrection.cxx:2
 AlidNdEtaCorrection.cxx:3
 AlidNdEtaCorrection.cxx:4
 AlidNdEtaCorrection.cxx:5
 AlidNdEtaCorrection.cxx:6
 AlidNdEtaCorrection.cxx:7
 AlidNdEtaCorrection.cxx:8
 AlidNdEtaCorrection.cxx:9
 AlidNdEtaCorrection.cxx:10
 AlidNdEtaCorrection.cxx:11
 AlidNdEtaCorrection.cxx:12
 AlidNdEtaCorrection.cxx:13
 AlidNdEtaCorrection.cxx:14
 AlidNdEtaCorrection.cxx:15
 AlidNdEtaCorrection.cxx:16
 AlidNdEtaCorrection.cxx:17
 AlidNdEtaCorrection.cxx:18
 AlidNdEtaCorrection.cxx:19
 AlidNdEtaCorrection.cxx:20
 AlidNdEtaCorrection.cxx:21
 AlidNdEtaCorrection.cxx:22
 AlidNdEtaCorrection.cxx:23
 AlidNdEtaCorrection.cxx:24
 AlidNdEtaCorrection.cxx:25
 AlidNdEtaCorrection.cxx:26
 AlidNdEtaCorrection.cxx:27
 AlidNdEtaCorrection.cxx:28
 AlidNdEtaCorrection.cxx:29
 AlidNdEtaCorrection.cxx:30
 AlidNdEtaCorrection.cxx:31
 AlidNdEtaCorrection.cxx:32
 AlidNdEtaCorrection.cxx:33
 AlidNdEtaCorrection.cxx:34
 AlidNdEtaCorrection.cxx:35
 AlidNdEtaCorrection.cxx:36
 AlidNdEtaCorrection.cxx:37
 AlidNdEtaCorrection.cxx:38
 AlidNdEtaCorrection.cxx:39
 AlidNdEtaCorrection.cxx:40
 AlidNdEtaCorrection.cxx:41
 AlidNdEtaCorrection.cxx:42
 AlidNdEtaCorrection.cxx:43
 AlidNdEtaCorrection.cxx:44
 AlidNdEtaCorrection.cxx:45
 AlidNdEtaCorrection.cxx:46
 AlidNdEtaCorrection.cxx:47
 AlidNdEtaCorrection.cxx:48
 AlidNdEtaCorrection.cxx:49
 AlidNdEtaCorrection.cxx:50
 AlidNdEtaCorrection.cxx:51
 AlidNdEtaCorrection.cxx:52
 AlidNdEtaCorrection.cxx:53
 AlidNdEtaCorrection.cxx:54
 AlidNdEtaCorrection.cxx:55
 AlidNdEtaCorrection.cxx:56
 AlidNdEtaCorrection.cxx:57
 AlidNdEtaCorrection.cxx:58
 AlidNdEtaCorrection.cxx:59
 AlidNdEtaCorrection.cxx:60
 AlidNdEtaCorrection.cxx:61
 AlidNdEtaCorrection.cxx:62
 AlidNdEtaCorrection.cxx:63
 AlidNdEtaCorrection.cxx:64
 AlidNdEtaCorrection.cxx:65
 AlidNdEtaCorrection.cxx:66
 AlidNdEtaCorrection.cxx:67
 AlidNdEtaCorrection.cxx:68
 AlidNdEtaCorrection.cxx:69
 AlidNdEtaCorrection.cxx:70
 AlidNdEtaCorrection.cxx:71
 AlidNdEtaCorrection.cxx:72
 AlidNdEtaCorrection.cxx:73
 AlidNdEtaCorrection.cxx:74
 AlidNdEtaCorrection.cxx:75
 AlidNdEtaCorrection.cxx:76
 AlidNdEtaCorrection.cxx:77
 AlidNdEtaCorrection.cxx:78
 AlidNdEtaCorrection.cxx:79
 AlidNdEtaCorrection.cxx:80
 AlidNdEtaCorrection.cxx:81
 AlidNdEtaCorrection.cxx:82
 AlidNdEtaCorrection.cxx:83
 AlidNdEtaCorrection.cxx:84
 AlidNdEtaCorrection.cxx:85
 AlidNdEtaCorrection.cxx:86
 AlidNdEtaCorrection.cxx:87
 AlidNdEtaCorrection.cxx:88
 AlidNdEtaCorrection.cxx:89
 AlidNdEtaCorrection.cxx:90
 AlidNdEtaCorrection.cxx:91
 AlidNdEtaCorrection.cxx:92
 AlidNdEtaCorrection.cxx:93
 AlidNdEtaCorrection.cxx:94
 AlidNdEtaCorrection.cxx:95
 AlidNdEtaCorrection.cxx:96
 AlidNdEtaCorrection.cxx:97
 AlidNdEtaCorrection.cxx:98
 AlidNdEtaCorrection.cxx:99
 AlidNdEtaCorrection.cxx:100
 AlidNdEtaCorrection.cxx:101
 AlidNdEtaCorrection.cxx:102
 AlidNdEtaCorrection.cxx:103
 AlidNdEtaCorrection.cxx:104
 AlidNdEtaCorrection.cxx:105
 AlidNdEtaCorrection.cxx:106
 AlidNdEtaCorrection.cxx:107
 AlidNdEtaCorrection.cxx:108
 AlidNdEtaCorrection.cxx:109
 AlidNdEtaCorrection.cxx:110
 AlidNdEtaCorrection.cxx:111
 AlidNdEtaCorrection.cxx:112
 AlidNdEtaCorrection.cxx:113
 AlidNdEtaCorrection.cxx:114
 AlidNdEtaCorrection.cxx:115
 AlidNdEtaCorrection.cxx:116
 AlidNdEtaCorrection.cxx:117
 AlidNdEtaCorrection.cxx:118
 AlidNdEtaCorrection.cxx:119
 AlidNdEtaCorrection.cxx:120
 AlidNdEtaCorrection.cxx:121
 AlidNdEtaCorrection.cxx:122
 AlidNdEtaCorrection.cxx:123
 AlidNdEtaCorrection.cxx:124
 AlidNdEtaCorrection.cxx:125
 AlidNdEtaCorrection.cxx:126
 AlidNdEtaCorrection.cxx:127
 AlidNdEtaCorrection.cxx:128
 AlidNdEtaCorrection.cxx:129
 AlidNdEtaCorrection.cxx:130
 AlidNdEtaCorrection.cxx:131
 AlidNdEtaCorrection.cxx:132
 AlidNdEtaCorrection.cxx:133
 AlidNdEtaCorrection.cxx:134
 AlidNdEtaCorrection.cxx:135
 AlidNdEtaCorrection.cxx:136
 AlidNdEtaCorrection.cxx:137
 AlidNdEtaCorrection.cxx:138
 AlidNdEtaCorrection.cxx:139
 AlidNdEtaCorrection.cxx:140
 AlidNdEtaCorrection.cxx:141
 AlidNdEtaCorrection.cxx:142
 AlidNdEtaCorrection.cxx:143
 AlidNdEtaCorrection.cxx:144
 AlidNdEtaCorrection.cxx:145
 AlidNdEtaCorrection.cxx:146
 AlidNdEtaCorrection.cxx:147
 AlidNdEtaCorrection.cxx:148
 AlidNdEtaCorrection.cxx:149
 AlidNdEtaCorrection.cxx:150
 AlidNdEtaCorrection.cxx:151
 AlidNdEtaCorrection.cxx:152
 AlidNdEtaCorrection.cxx:153
 AlidNdEtaCorrection.cxx:154
 AlidNdEtaCorrection.cxx:155
 AlidNdEtaCorrection.cxx:156
 AlidNdEtaCorrection.cxx:157
 AlidNdEtaCorrection.cxx:158
 AlidNdEtaCorrection.cxx:159
 AlidNdEtaCorrection.cxx:160
 AlidNdEtaCorrection.cxx:161
 AlidNdEtaCorrection.cxx:162
 AlidNdEtaCorrection.cxx:163
 AlidNdEtaCorrection.cxx:164
 AlidNdEtaCorrection.cxx:165
 AlidNdEtaCorrection.cxx:166
 AlidNdEtaCorrection.cxx:167
 AlidNdEtaCorrection.cxx:168
 AlidNdEtaCorrection.cxx:169
 AlidNdEtaCorrection.cxx:170
 AlidNdEtaCorrection.cxx:171
 AlidNdEtaCorrection.cxx:172
 AlidNdEtaCorrection.cxx:173
 AlidNdEtaCorrection.cxx:174
 AlidNdEtaCorrection.cxx:175
 AlidNdEtaCorrection.cxx:176
 AlidNdEtaCorrection.cxx:177
 AlidNdEtaCorrection.cxx:178
 AlidNdEtaCorrection.cxx:179
 AlidNdEtaCorrection.cxx:180
 AlidNdEtaCorrection.cxx:181
 AlidNdEtaCorrection.cxx:182
 AlidNdEtaCorrection.cxx:183
 AlidNdEtaCorrection.cxx:184
 AlidNdEtaCorrection.cxx:185
 AlidNdEtaCorrection.cxx:186
 AlidNdEtaCorrection.cxx:187
 AlidNdEtaCorrection.cxx:188
 AlidNdEtaCorrection.cxx:189
 AlidNdEtaCorrection.cxx:190
 AlidNdEtaCorrection.cxx:191
 AlidNdEtaCorrection.cxx:192
 AlidNdEtaCorrection.cxx:193
 AlidNdEtaCorrection.cxx:194
 AlidNdEtaCorrection.cxx:195
 AlidNdEtaCorrection.cxx:196
 AlidNdEtaCorrection.cxx:197
 AlidNdEtaCorrection.cxx:198
 AlidNdEtaCorrection.cxx:199
 AlidNdEtaCorrection.cxx:200
 AlidNdEtaCorrection.cxx:201
 AlidNdEtaCorrection.cxx:202
 AlidNdEtaCorrection.cxx:203
 AlidNdEtaCorrection.cxx:204
 AlidNdEtaCorrection.cxx:205
 AlidNdEtaCorrection.cxx:206
 AlidNdEtaCorrection.cxx:207
 AlidNdEtaCorrection.cxx:208
 AlidNdEtaCorrection.cxx:209
 AlidNdEtaCorrection.cxx:210
 AlidNdEtaCorrection.cxx:211
 AlidNdEtaCorrection.cxx:212
 AlidNdEtaCorrection.cxx:213
 AlidNdEtaCorrection.cxx:214
 AlidNdEtaCorrection.cxx:215
 AlidNdEtaCorrection.cxx:216
 AlidNdEtaCorrection.cxx:217
 AlidNdEtaCorrection.cxx:218
 AlidNdEtaCorrection.cxx:219
 AlidNdEtaCorrection.cxx:220
 AlidNdEtaCorrection.cxx:221
 AlidNdEtaCorrection.cxx:222
 AlidNdEtaCorrection.cxx:223
 AlidNdEtaCorrection.cxx:224
 AlidNdEtaCorrection.cxx:225
 AlidNdEtaCorrection.cxx:226
 AlidNdEtaCorrection.cxx:227
 AlidNdEtaCorrection.cxx:228
 AlidNdEtaCorrection.cxx:229
 AlidNdEtaCorrection.cxx:230
 AlidNdEtaCorrection.cxx:231
 AlidNdEtaCorrection.cxx:232
 AlidNdEtaCorrection.cxx:233
 AlidNdEtaCorrection.cxx:234
 AlidNdEtaCorrection.cxx:235
 AlidNdEtaCorrection.cxx:236
 AlidNdEtaCorrection.cxx:237
 AlidNdEtaCorrection.cxx:238
 AlidNdEtaCorrection.cxx:239
 AlidNdEtaCorrection.cxx:240
 AlidNdEtaCorrection.cxx:241
 AlidNdEtaCorrection.cxx:242
 AlidNdEtaCorrection.cxx:243
 AlidNdEtaCorrection.cxx:244
 AlidNdEtaCorrection.cxx:245
 AlidNdEtaCorrection.cxx:246
 AlidNdEtaCorrection.cxx:247
 AlidNdEtaCorrection.cxx:248
 AlidNdEtaCorrection.cxx:249
 AlidNdEtaCorrection.cxx:250
 AlidNdEtaCorrection.cxx:251
 AlidNdEtaCorrection.cxx:252
 AlidNdEtaCorrection.cxx:253
 AlidNdEtaCorrection.cxx:254
 AlidNdEtaCorrection.cxx:255
 AlidNdEtaCorrection.cxx:256
 AlidNdEtaCorrection.cxx:257
 AlidNdEtaCorrection.cxx:258
 AlidNdEtaCorrection.cxx:259
 AlidNdEtaCorrection.cxx:260
 AlidNdEtaCorrection.cxx:261
 AlidNdEtaCorrection.cxx:262
 AlidNdEtaCorrection.cxx:263
 AlidNdEtaCorrection.cxx:264
 AlidNdEtaCorrection.cxx:265
 AlidNdEtaCorrection.cxx:266
 AlidNdEtaCorrection.cxx:267
 AlidNdEtaCorrection.cxx:268
 AlidNdEtaCorrection.cxx:269
 AlidNdEtaCorrection.cxx:270
 AlidNdEtaCorrection.cxx:271
 AlidNdEtaCorrection.cxx:272
 AlidNdEtaCorrection.cxx:273
 AlidNdEtaCorrection.cxx:274
 AlidNdEtaCorrection.cxx:275
 AlidNdEtaCorrection.cxx:276
 AlidNdEtaCorrection.cxx:277
 AlidNdEtaCorrection.cxx:278
 AlidNdEtaCorrection.cxx:279
 AlidNdEtaCorrection.cxx:280
 AlidNdEtaCorrection.cxx:281
 AlidNdEtaCorrection.cxx:282
 AlidNdEtaCorrection.cxx:283
 AlidNdEtaCorrection.cxx:284
 AlidNdEtaCorrection.cxx:285
 AlidNdEtaCorrection.cxx:286
 AlidNdEtaCorrection.cxx:287
 AlidNdEtaCorrection.cxx:288
 AlidNdEtaCorrection.cxx:289
 AlidNdEtaCorrection.cxx:290
 AlidNdEtaCorrection.cxx:291
 AlidNdEtaCorrection.cxx:292
 AlidNdEtaCorrection.cxx:293
 AlidNdEtaCorrection.cxx:294
 AlidNdEtaCorrection.cxx:295
 AlidNdEtaCorrection.cxx:296
 AlidNdEtaCorrection.cxx:297
 AlidNdEtaCorrection.cxx:298
 AlidNdEtaCorrection.cxx:299
 AlidNdEtaCorrection.cxx:300
 AlidNdEtaCorrection.cxx:301
 AlidNdEtaCorrection.cxx:302
 AlidNdEtaCorrection.cxx:303
 AlidNdEtaCorrection.cxx:304
 AlidNdEtaCorrection.cxx:305
 AlidNdEtaCorrection.cxx:306
 AlidNdEtaCorrection.cxx:307
 AlidNdEtaCorrection.cxx:308
 AlidNdEtaCorrection.cxx:309
 AlidNdEtaCorrection.cxx:310
 AlidNdEtaCorrection.cxx:311
 AlidNdEtaCorrection.cxx:312
 AlidNdEtaCorrection.cxx:313
 AlidNdEtaCorrection.cxx:314
 AlidNdEtaCorrection.cxx:315
 AlidNdEtaCorrection.cxx:316
 AlidNdEtaCorrection.cxx:317
 AlidNdEtaCorrection.cxx:318
 AlidNdEtaCorrection.cxx:319
 AlidNdEtaCorrection.cxx:320
 AlidNdEtaCorrection.cxx:321
 AlidNdEtaCorrection.cxx:322
 AlidNdEtaCorrection.cxx:323
 AlidNdEtaCorrection.cxx:324
 AlidNdEtaCorrection.cxx:325
 AlidNdEtaCorrection.cxx:326
 AlidNdEtaCorrection.cxx:327
 AlidNdEtaCorrection.cxx:328
 AlidNdEtaCorrection.cxx:329
 AlidNdEtaCorrection.cxx:330
 AlidNdEtaCorrection.cxx:331
 AlidNdEtaCorrection.cxx:332
 AlidNdEtaCorrection.cxx:333
 AlidNdEtaCorrection.cxx:334
 AlidNdEtaCorrection.cxx:335
 AlidNdEtaCorrection.cxx:336
 AlidNdEtaCorrection.cxx:337
 AlidNdEtaCorrection.cxx:338
 AlidNdEtaCorrection.cxx:339
 AlidNdEtaCorrection.cxx:340
 AlidNdEtaCorrection.cxx:341
 AlidNdEtaCorrection.cxx:342
 AlidNdEtaCorrection.cxx:343
 AlidNdEtaCorrection.cxx:344
 AlidNdEtaCorrection.cxx:345
 AlidNdEtaCorrection.cxx:346
 AlidNdEtaCorrection.cxx:347
 AlidNdEtaCorrection.cxx:348
 AlidNdEtaCorrection.cxx:349
 AlidNdEtaCorrection.cxx:350
 AlidNdEtaCorrection.cxx:351
 AlidNdEtaCorrection.cxx:352
 AlidNdEtaCorrection.cxx:353
 AlidNdEtaCorrection.cxx:354
 AlidNdEtaCorrection.cxx:355
 AlidNdEtaCorrection.cxx:356
 AlidNdEtaCorrection.cxx:357
 AlidNdEtaCorrection.cxx:358
 AlidNdEtaCorrection.cxx:359
 AlidNdEtaCorrection.cxx:360
 AlidNdEtaCorrection.cxx:361
 AlidNdEtaCorrection.cxx:362
 AlidNdEtaCorrection.cxx:363
 AlidNdEtaCorrection.cxx:364
 AlidNdEtaCorrection.cxx:365
 AlidNdEtaCorrection.cxx:366
 AlidNdEtaCorrection.cxx:367
 AlidNdEtaCorrection.cxx:368
 AlidNdEtaCorrection.cxx:369
 AlidNdEtaCorrection.cxx:370
 AlidNdEtaCorrection.cxx:371
 AlidNdEtaCorrection.cxx:372
 AlidNdEtaCorrection.cxx:373
 AlidNdEtaCorrection.cxx:374
 AlidNdEtaCorrection.cxx:375
 AlidNdEtaCorrection.cxx:376
 AlidNdEtaCorrection.cxx:377
 AlidNdEtaCorrection.cxx:378
 AlidNdEtaCorrection.cxx:379
 AlidNdEtaCorrection.cxx:380
 AlidNdEtaCorrection.cxx:381
 AlidNdEtaCorrection.cxx:382
 AlidNdEtaCorrection.cxx:383
 AlidNdEtaCorrection.cxx:384
 AlidNdEtaCorrection.cxx:385
 AlidNdEtaCorrection.cxx:386
 AlidNdEtaCorrection.cxx:387
 AlidNdEtaCorrection.cxx:388
 AlidNdEtaCorrection.cxx:389
 AlidNdEtaCorrection.cxx:390
 AlidNdEtaCorrection.cxx:391
 AlidNdEtaCorrection.cxx:392
 AlidNdEtaCorrection.cxx:393
 AlidNdEtaCorrection.cxx:394
 AlidNdEtaCorrection.cxx:395
 AlidNdEtaCorrection.cxx:396
 AlidNdEtaCorrection.cxx:397
 AlidNdEtaCorrection.cxx:398
 AlidNdEtaCorrection.cxx:399
 AlidNdEtaCorrection.cxx:400
 AlidNdEtaCorrection.cxx:401
 AlidNdEtaCorrection.cxx:402
 AlidNdEtaCorrection.cxx:403
 AlidNdEtaCorrection.cxx:404
 AlidNdEtaCorrection.cxx:405
 AlidNdEtaCorrection.cxx:406
 AlidNdEtaCorrection.cxx:407
 AlidNdEtaCorrection.cxx:408
 AlidNdEtaCorrection.cxx:409
 AlidNdEtaCorrection.cxx:410
 AlidNdEtaCorrection.cxx:411
 AlidNdEtaCorrection.cxx:412
 AlidNdEtaCorrection.cxx:413
 AlidNdEtaCorrection.cxx:414
 AlidNdEtaCorrection.cxx:415
 AlidNdEtaCorrection.cxx:416
 AlidNdEtaCorrection.cxx:417
 AlidNdEtaCorrection.cxx:418
 AlidNdEtaCorrection.cxx:419
 AlidNdEtaCorrection.cxx:420
 AlidNdEtaCorrection.cxx:421
 AlidNdEtaCorrection.cxx:422
 AlidNdEtaCorrection.cxx:423
 AlidNdEtaCorrection.cxx:424
 AlidNdEtaCorrection.cxx:425
 AlidNdEtaCorrection.cxx:426
 AlidNdEtaCorrection.cxx:427
 AlidNdEtaCorrection.cxx:428
 AlidNdEtaCorrection.cxx:429
 AlidNdEtaCorrection.cxx:430
 AlidNdEtaCorrection.cxx:431
 AlidNdEtaCorrection.cxx:432
 AlidNdEtaCorrection.cxx:433
 AlidNdEtaCorrection.cxx:434
 AlidNdEtaCorrection.cxx:435
 AlidNdEtaCorrection.cxx:436
 AlidNdEtaCorrection.cxx:437
 AlidNdEtaCorrection.cxx:438
 AlidNdEtaCorrection.cxx:439
 AlidNdEtaCorrection.cxx:440
 AlidNdEtaCorrection.cxx:441
 AlidNdEtaCorrection.cxx:442
 AlidNdEtaCorrection.cxx:443
 AlidNdEtaCorrection.cxx:444
 AlidNdEtaCorrection.cxx:445
 AlidNdEtaCorrection.cxx:446
 AlidNdEtaCorrection.cxx:447
 AlidNdEtaCorrection.cxx:448
 AlidNdEtaCorrection.cxx:449
 AlidNdEtaCorrection.cxx:450
 AlidNdEtaCorrection.cxx:451
 AlidNdEtaCorrection.cxx:452
 AlidNdEtaCorrection.cxx:453
 AlidNdEtaCorrection.cxx:454
 AlidNdEtaCorrection.cxx:455
 AlidNdEtaCorrection.cxx:456
 AlidNdEtaCorrection.cxx:457
 AlidNdEtaCorrection.cxx:458
 AlidNdEtaCorrection.cxx:459
 AlidNdEtaCorrection.cxx:460
 AlidNdEtaCorrection.cxx:461
 AlidNdEtaCorrection.cxx:462
 AlidNdEtaCorrection.cxx:463
 AlidNdEtaCorrection.cxx:464
 AlidNdEtaCorrection.cxx:465
 AlidNdEtaCorrection.cxx:466
 AlidNdEtaCorrection.cxx:467
 AlidNdEtaCorrection.cxx:468
 AlidNdEtaCorrection.cxx:469
 AlidNdEtaCorrection.cxx:470
 AlidNdEtaCorrection.cxx:471
 AlidNdEtaCorrection.cxx:472
 AlidNdEtaCorrection.cxx:473
 AlidNdEtaCorrection.cxx:474
 AlidNdEtaCorrection.cxx:475
 AlidNdEtaCorrection.cxx:476
 AlidNdEtaCorrection.cxx:477
 AlidNdEtaCorrection.cxx:478
 AlidNdEtaCorrection.cxx:479
 AlidNdEtaCorrection.cxx:480
 AlidNdEtaCorrection.cxx:481
 AlidNdEtaCorrection.cxx:482
 AlidNdEtaCorrection.cxx:483
 AlidNdEtaCorrection.cxx:484
 AlidNdEtaCorrection.cxx:485
 AlidNdEtaCorrection.cxx:486
 AlidNdEtaCorrection.cxx:487