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$ */

/////////////////////////////////////////////////////////////////////////////
//
//  Contain parametrizations to generate atmospheric muons, and also
//  to generate single muons and muon bundles at surface level.
//
//Begin_Html
/*
<img src="picts/AliGenACORDEClass.gif">
</pre>
<br clear=left>
<font size=+2 color=red>
<p>The responsible person for this module is
<a href="mailto:Enrique.Gamez.Flores@cern.ch">Enrique Gamez</a>.
</font>
<pre>
*/
//End_Html
//
/////////////////////////////////////////////////////////////////////////////

#include "AliGenACORDE.h"

#include <TMCProcess.h>
#include <TPDGCode.h>
#include <TClonesArray.h>
#include <TF1.h>
#include <TH1F.h>

#include "AliRun.h"
#include "AliConst.h"

ClassImp(AliGenACORDE)

//_____________________________________________________________________________
AliGenACORDE::AliGenACORDE()
  : AliGenerator(),
    fIpart(0),
    fCRMode(kSingleMuons),
    fCRModeName(0),
    fXwidth(0),
    fNx(1),
    fZwidth(0),
    fNz(1),
    fMuonGrid(kFALSE),
    fZenithMin(0),
    fZenithMax(0),
    fAzimuthMin(0),
    fAzimuthMax(0),
    fPRange(0),
    fPResolution(1),
    fAp(0),
    fMomentumDist(0),
    fUnfoldedMomentumDist(0),
    fZenithDist(0),
    fPDist(0),
    fNParticles(0)
{
  //
  // Default ctor.
  //
}

//_____________________________________________________________________________
AliGenACORDE::AliGenACORDE(Int_t npart) 
  : AliGenerator(npart),
    fIpart(kMuonMinus),
    fCRMode(kSingleMuons),
    fCRModeName(0),
    fXwidth(0),
    fNx(1),
    fZwidth(0),
    fNz(1),
    fMuonGrid(kFALSE),
    fZenithMin(0),
    fZenithMax(0),
    fAzimuthMin(0),
    fAzimuthMax(0),
    fPRange(0),
    fPResolution(1),
    fAp(0),
    fMomentumDist(0),
    fUnfoldedMomentumDist(0),
    fZenithDist(0),
    fPDist(0),
    fNParticles(0)
{
  //
  // Standard ctor.
  //
  fName = "ACORDE";
  fTitle = "Cosmic Muons generator";

  // Set the origin above the vertex, on the surface.
  fOrigin[0] = 0.;
  fOrigin[1] = AliACORDEConstants::Instance()->Depth(); // At the surface by default.
  fOrigin[2] = 0.;
}

//_____________________________________________________________________________
AliGenACORDE::~AliGenACORDE()
{
  //
  // Default dtor.
  //
  if ( fPDist ) {fPDist->Delete(); delete fPDist; fPDist = 0;}
  if ( fUnfoldedMomentumDist ) delete fUnfoldedMomentumDist;
  if ( fMomentumDist ) delete fMomentumDist;
  if ( fAp )           delete fAp;
  if ( fCRModeName )   delete fCRModeName;
}

//_____________________________________________________________________________
void AliGenACORDE::Generate()
{
  //
  // Generate on one trigger
  // Call the respective method inside the loop for the number
  // of tracks per trigger.

  for (Int_t i = 0; i < fNParticles; i++ ) {

    if ( fCRMode == kMuonBundle ) {
      this->GenerateOneMuonBundle();

    } else if ( fCRMode == kSingleMuons ) {
      this->GenerateOneSingleMuon(kTRUE);

    } else {
      // Generate only single muons following the parametrizations
      // for atmospheric muons.
      this->GenerateOneSingleMuon();

    }

  }
}

//_____________________________________________________________________________
void AliGenACORDE::Init()
{
  //
  // Initialize some internal methods.
  //


   printf("**************************************************************\n");
   printf("<<< *** Starting the AliGenACORDE cosmic generator ******** >>>\n");
   printf("<<< *** No. of muons generated at the surface of P2: %d,  * >>>\n",fNParticles);
   printf("**************************************************************\n");

  // Determine some specific data members.
  fPRange = TMath::Abs(fPMax-fPMin);

  if ( fCRMode == kSingleMuons ) {
    fCRModeName = new TString("Single Muons");
    // Initialisation, check consistency of selected ranges
    if(TestBit(kPtRange)&&TestBit(kMomentumRange)) 
      Fatal("Init","You should not set the momentum range and the pt range!");
    
    if((!TestBit(kPtRange))&&(!TestBit(kMomentumRange))) 
      Fatal("Init","You should set either the momentum or the pt range!");
    
  } else if ( fCRMode == kMuonBundle ) {
    fCRModeName = new TString("Muon Bundles");

  } else if ( fCRMode == kMuonFlux ) {
    fCRModeName = new TString("Muon Fluxes");
    // Initialize the ditribution functions.
    this->InitMomentumGeneration();
    this->InitZenithalAngleGeneration();
    
  } else {
    Fatal("Generate", "Generation Mode unknown!\n");

  }

}

//____________________________________________________________________________
void AliGenACORDE::GenerateOneSingleMuon(Bool_t withFlatMomentum)
{
  //
  // Generate One Single Muon
  // This method will generate only one muon.
  // The momentum will be randomly flat distributed if
  // the paremeter "withFlatMomentum" is set to kTRUE,
  // otherwise the momentum will generate acordingly the parametrization
  // given by 
  // and adpted from Bruno Alessandro's implementation with the
  // CERNLIB to AliRoot.
  // The "withFlatMomentum" parameter also will be used to generate
  // the muons with a flat Zenithal angle distribution.
  // Do the smearing here, so that means per track.

  Float_t polar[3]= {0,0,0}; // Polarization parameters
  Float_t origin[3];
  Int_t nt;
  Float_t p[3];
  Float_t pmom, pt;
  Float_t random[6];

  // Take the azimuth random.
  Rndm(random, 2);
  Float_t azimuth = fAzimuthMin + (fAzimuthMax-fAzimuthMin)*random[0];
  Float_t zenith = fZenithMin + (fZenithMax - fZenithMin)*random[1];

  if ( withFlatMomentum ) {
    Rndm(random,3);
    if(TestBit(kMomentumRange)) {
      pmom = -( fPMin + random[0]*(fPMax - fPMin) ); // always downwards.
      pt = pmom*TMath::Sin(zenith*kDegrad);
    } else {
      pt = -( fPtMin + random[1]*(fPtMax - fPtMin)); // always downwards.
      pmom = pt/TMath::Sin(zenith*kDegrad);
    }

  } else {
    if ( fMomentumDist ) {
      pmom = -this->GetMomentum(); // Always downwards.
    } else {
      pmom = -fPMin;
    }
    zenith = this->GetZenithAngle(pmom);  // In degrees
    pt = pmom*TMath::Sin(zenith*kDegrad);
  }

  p[0] = pt*TMath::Sin(azimuth*kDegrad);
  p[1] = pmom*TMath::Cos(zenith*kDegrad);
  p[2] = pt*TMath::Cos(azimuth*kDegrad);

  // Finaly the origin, with the smearing
  Rndm(random,6);
  origin[0] = AliACORDEConstants::Instance()->Depth()*TMath::Tan(zenith*kDegrad)*
    TMath::Sin(azimuth*kDegrad)
    + fOsigma[0]* TMath::Cos(2*random[0]*TMath::Pi())*TMath::Sqrt(-2*TMath::Log(random[1]));

  origin[1] = AliACORDEConstants::Instance()->Depth();

  origin[2] = AliACORDEConstants::Instance()->Depth()*TMath::Tan(zenith*kDegrad)*
    TMath::Cos(azimuth*kDegrad)
    + fOsigma[2]* TMath::Cos(2*random[2]*TMath::Pi())*TMath::Sqrt(-2*TMath::Log(random[3]));

  // Put the track on the stack.
  PushTrack(fTrackIt,-1,fIpart,p,origin,polar,0,kPPrimary,nt);

}

//____________________________________________________________________________
void AliGenACORDE::GenerateOneMuonBundle()
{
  //
  // Generate One Muon Bundle method
  // This method will generate a bunch of muons following the
  // procedure of the AliGenScan class.
  // These muons will be generated with flat momentum.

  Float_t polar[3]= {0,0,0}; // Polarization parameters
  Float_t origin[3];
  Float_t p[3];
  Int_t nt;
  Float_t pmom;
  Float_t random[6];

  Rndm(random, 3);
  Float_t zenith = fZenithMin + (fZenithMax - fZenithMin)*random[1];
  Float_t azimuth = fAzimuthMin + (fAzimuthMax-fAzimuthMin)*random[2];
  //Float_t zenith = 10;
  //Float_t azimuth = 30;

  // Generate the kinematics a la AliGenScan (Andreas Morchs)
  Float_t dx, dz;
  if ( fNx > 0 ) {
    dx = fXwidth/fNx;
  } else {
    dx = 1e10;
    //dx = 100.;
  }

  if ( fNz > 0 ) {
    dz = fZwidth/fNz;
  } else {
    dz = 1e10;
    //dz = 100.;
  }

  origin[0] = AliACORDEConstants::Instance()->Depth()*TMath::Tan(zenith*kDegrad)*
              TMath::Sin(azimuth*kDegrad);
  //origin[0] = 0.;
  origin[1] = AliACORDEConstants::Instance()->Depth();
  //origin[1] = 900;
  origin[2] = AliACORDEConstants::Instance()->Depth()*TMath::Tan(zenith*kDegrad)*
              TMath::Cos(azimuth*kDegrad);
    //origin[2] = 0.;

  for (Int_t ix = 0; ix < fNx; ix++ ) {
    for (Int_t iz = 0; iz < fNz; iz++ ) {
      Rndm(random,6);
      origin[0]+=ix*dx+2*(random[1]-0.5)*fOsigma[0];
      origin[2]+=iz*dz+2*(random[2]-0.5)*fOsigma[2];
      if ( random[4] < 0.5 ) {
        origin[0] = -1*origin[0];
      }
      if ( random[5] < 0.5 ) {
        origin[2] = -1*origin[2];
      }

      pmom = -(fPMin + random[3] *(fPMax - fPMax) ); // Always downwards
      p[0] = TMath::Sin(zenith*kDegrad)*TMath::Sin(azimuth*kDegrad)*pmom;
      p[1] = TMath::Cos(zenith*kDegrad)*pmom;
      p[2] = TMath::Sin(zenith*kDegrad)*TMath::Cos(azimuth*kDegrad)*pmom;

      PushTrack(fTrackIt, -1, fIpart, p, origin, polar, 0, kPPrimary, nt);
    }

  }

}

//____________________________________________________________________________
void AliGenACORDE::SetGridRange(Int_t nx,Float_t xwidth, Int_t nz, Float_t zwidth)
{
  //
  // Define the grid
  // This data shuold be used for Muon bundles generation.
  //
  fXwidth=xwidth;
  fNx=nx;
  fZwidth=zwidth;
  fNz=nz;

  // Print a message  about the use, if the Mode has not been set, or
  // it has to a different Mode.
  if ( fCRMode != kMuonBundle ) {
    Warning("SetRange","You have been specified a grid to generate muon bundles, but seems that you haven't choose this generation mode, or you have already select a different one");
    fMuonGrid = kTRUE;
  }
}

//____________________________________________________________________________
void AliGenACORDE::InitApWeightFactors()
{
  //
  // This factors will be  to correct the zenithal angle distribution
  // acording the momentum

  //
  // Fill the array for the flux zenith angle dependence.
  // at the index 0 of fAp[] will be the "factor" if we have a muon
  // of 0 GeV.
  Float_t a[6] = {-1.61, -1.50, -1.28, -0.94, -0.61, -0.22};
  Float_t p[6] = { 0., 10., 30., 100., 300., 1000.};

  // Get the information from the momentum
  Int_t pEnd  = TMath::Abs(TMath::Nint(fPMax/fPResolution)) + 1;
  // Initialize the Array of floats to hold the a(p) factors.
  fAp = new TArrayF(pEnd);
  
  Int_t index = 0;

  for (Int_t i = 0; i < pEnd; i++ ) {
    Float_t currentP = ((Float_t)i)*fPResolution;
    if ( currentP < p[1] )                          index = 0;
    else if ( currentP >= p[1] && currentP < p[2] ) index = 1;
    else if ( currentP >= p[2] && currentP < p[3] ) index = 2;
    else if ( currentP >= p[3] && currentP < p[4] ) index = 3;
    else if ( currentP >= p[4] )                    index = 4;

    Float_t ap = (currentP -p[index])*(a[index+1] - a[index])/
                 (p[index+1] - p[index]) + a[index];
    fAp->AddAt(ap, i);
  }

}

//___________________________________________________________________________
void AliGenACORDE::InitMomentumGeneration()
{
  //
  // Initialize a funtion to generate the momentum randomly
  // acording this function.
  //

  // Check if we nned to initialize the function
  if ( fPMin != fPMax ) {

    // Check also if the function have been defined yet.
    if ( !fMomentumDist ) {

      // If not, use this function
      const char* y      = "log10(x)";
      
      const char* h1Coef = "[0]*( %s*%s*%s/2 - (5*%s*%s/2) + 3*%s )";
      const char* h2Coef = "[1]*( (-2*%s*%s*%s/3) + (3*%s*%s) - 10*%s/3 + 1 )";
      const char* h3Coef = "[2]*( %s*%s*%s/6 - %s*%s/2 + %s/3 )";
      const char* s2Coef = "[3]*( %s*%s*%s/3 - 2*%s*%s + 11*%s/3 - 2 )";
      
      const char* h = "%s + %s + %s + %s";
      const char* flux = "pow(10., %s)";
      const char* normalizedFlux = "0.86*x*x*x*pow(10., %s)";
      const char* paramNames[4] = {"H1", "H2", "H3", "S1"};
      
      char buffer1[1024];
      char buffer2[1024];
      char buffer3[1024];
      char buffer4[1024];
      char buffer5[1024];
      char buffer6[1024];
      char buffer7[1024];

      snprintf(buffer1, 1023, h1Coef, y, y, y, y, y, y);
      snprintf(buffer2, 1023, h2Coef, y, y, y, y, y, y);
      snprintf(buffer3, 1023, h3Coef, y, y, y, y, y, y);
      snprintf(buffer4, 1023, s2Coef, y, y, y, y, y, y);
      
      snprintf(buffer5, 1023, h, buffer1, buffer2, buffer3, buffer4);
      
      snprintf(buffer6, 1023, flux, buffer5);
      
      fMomentumDist = new TF1("fMomentumDist", buffer6, fPMin, fPMax);
      snprintf(buffer7, 1023, normalizedFlux, buffer5);
      fUnfoldedMomentumDist = new TF1("fUnfoldedMomentumDist", buffer7, fPMin, fPMax);
      for (Int_t i = 0; i < 4; i++ ) {
	fMomentumDist->SetParName(i, paramNames[i]);
	fUnfoldedMomentumDist->SetParName(i, paramNames[i]);
      }
      
      fMomentumDist->SetParameter(0, 0.133);
      fMomentumDist->SetParameter(1, -2.521);
      fMomentumDist->SetParameter(2, -5.78);
      fMomentumDist->SetParameter(3, -2.11);

      fUnfoldedMomentumDist->SetParameter(0, 0.133);
      fUnfoldedMomentumDist->SetParameter(1, -2.521);
      fUnfoldedMomentumDist->SetParameter(2, -5.78);
      fUnfoldedMomentumDist->SetParameter(3, -2.11);
      
    }

  }

}

//____________________________________________________________________________
void AliGenACORDE::InitZenithalAngleGeneration()
{
  //
  // Initalize a distribution function for the zenith angle.
  // This angle will be obtained randomly acording this function.
  // The generated angles  will been in degrees.

  // Check if we need to create the function.
  if ( fZenithMin != fZenithMax ) {

    // Check also if another function have been defined.
    if ( !fZenithDist ) {
      
      // initialize the momentum dependent coefficients, a(p) 
      this->InitApWeightFactors();

      Int_t pEnd  = TMath::Abs(TMath::Nint(fPMax/fPResolution)) + 1;
      char name[26];
      char title[52];
      fPDist = new TClonesArray("TH1F", pEnd);
      TClonesArray &mom = *fPDist;
      TH1F* zenith = 0;
      Float_t weight = 0;
      for ( Int_t i = 0; i < pEnd; i++ ) {
	// Fill the distribution
	snprintf(name, 25, "zenith%d", i+1);
	snprintf(title, 51, "Zenith distribution, p=%f", fPMin+(Float_t)i);
	zenith = new(mom[i]) TH1F(name, title, TMath::Abs(TMath::Nint(fZenithMax-fZenithMin)), TMath::Cos(fZenithMax*TMath::Pi()/180), TMath::Cos(fZenithMin*TMath::Pi()/180));

	// Make a loop for the angle and fill the histogram for the weight
	Int_t steps = 1000;
	Float_t value = 0;
	for (Int_t j = 0; j < steps; j++ ) {
	  value = TMath::Cos(fZenithMin*TMath::Pi()/180) + (Float_t)j * ( TMath::Cos(fZenithMax*TMath::Pi()/180) - TMath::Cos(fZenithMin*TMath::Pi()/180))/1000;
	  weight = 1 + fAp->At(i)*(1 - value);
	  zenith->Fill(value, weight);
	}

      }

    } 

  }

}

//____________________________________________________________________________
Float_t AliGenACORDE::GetZenithAngle(Float_t mom) const
{

  Float_t zenith = 0.;
  // Check if you need to generate a constant zenith angle.
  if ( !fZenithDist ) {
    // Check if you have defined an array of momentum functions
    if ( fPDist ) {
      Int_t pIndex = TMath::Abs(TMath::Nint(mom));
      TH1F* cosZenithAngle = (TH1F*)fPDist->UncheckedAt(pIndex);
      Float_t tmpzenith = TMath::ACos(cosZenithAngle->GetRandom());
      // Correct the value
      zenith = kRaddeg*tmpzenith;
      return zenith;
    } else {

      if ( fCRMode != kMuonFlux ) {
	// If you aren't generating muons obeying any ditribution
	// only generate a flat zenith angle, acording the input settings
	Float_t random[2];
	Rndm(random, 2);
	zenith = fZenithMin + (fZenithMax - fZenithMin)*random[0];

      } else {
	// Even if you are generating muons acording some distribution,
	// but you don't want to ...
	zenith = fZenithMin;

      }

    }
  } else {
    zenith = fZenithDist->GetRandom();
  }

  return zenith;
}

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