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

// Read background particles from a boundary source
// Very specialized generator to simulate background from beam halo.
// The input file is a text file specially prepared 
// for this purpose.
// Author: andreas.morsch@cern.ch

#include <stdlib.h>

#include <TDatabasePDG.h>
#include <TCanvas.h>
#include <TGraph.h>
#include <TPDGCode.h>
#include <TSystem.h>

#include "AliLog.h"
#include "AliGenHaloProtvino.h"
#include "AliRun.h"

ClassImp(AliGenHaloProtvino)

AliGenHaloProtvino::AliGenHaloProtvino()
    :AliGenerator(-1), 
     fFile(0),
     fFileName(0),
     fSide(1),
     fRunPeriod(kY3D90),
     fTimePerEvent(1.e-4),
     fNskip(0),
     fZ1(0),
     fZ2(0),
     fG1(0),
     fG2(0),
     fGPASize(0)
{
// Constructor
    
    fName  = "HaloProtvino";
    fTitle = "Halo from LHC Tunnel";
//
//  Read all particles
    fNpart = -1;
    SetAnalog(0);
}

AliGenHaloProtvino::AliGenHaloProtvino(Int_t npart)
    :AliGenerator(npart),
     fFile(0),
     fFileName(0),
     fSide(1),
     fRunPeriod(kY3D90),
     fTimePerEvent(1.e-4),
     fNskip(0),
     fZ1(0),
     fZ2(0),
     fG1(0),
     fG2(0),
     fGPASize(0)
{
// Constructor
    fName = "Halo";
    fTitle= "Halo from LHC Tunnel";
//
    fNpart   = npart;
//
    SetAnalog(0);
}

//____________________________________________________________
AliGenHaloProtvino::~AliGenHaloProtvino()
{
// Destructor
}

//____________________________________________________________
void AliGenHaloProtvino::Init() 
{
// Initialisation
    fFile = fopen(fFileName,"r");
    if (fFile) {
	printf("\n File %s opened for reading, %p ! \n ",  fFileName.Data(), (void*)fFile);
    } else {
	printf("\n Opening of file %s failed,  %p ! \n ",  fFileName.Data(), (void*)fFile);
    }
//
//
//
//    Read file with gas pressure values
    char *name = 0;
    if (fRunPeriod < 5) {
	name = gSystem->ExpandPathName("$(ALICE_ROOT)/LHC/gasPressure.dat" );
	fGPASize = 21;
	fG1 = new Float_t[fGPASize];
	fG2 = new Float_t[fGPASize];
	fZ1 = new Float_t[fGPASize];
	fZ2 = new Float_t[fGPASize];
    } else if (fRunPeriod == 5) {
	name = gSystem->ExpandPathName("$(ALICE_ROOT)/LHC/pressure_2003_startup.dat");
	fGPASize = 18853;
	fG1 = new Float_t[fGPASize];
	fZ1 = new Float_t[fGPASize];
    } else if (fRunPeriod ==6) {
	name = gSystem->ExpandPathName("$(ALICE_ROOT)/LHC/pressure_2003_conditioned.dat");
	fGPASize = 12719;
	fG1 = new Float_t[fGPASize];
	fZ1 = new Float_t[fGPASize];
    } else {
	Fatal("Init()", "No gas pressure file for given run period !");
    }

    FILE* file = 0;
    if (name) file = fopen(name, "r");
    if (!file) {
	AliError("No gas pressure file");
	return;
    }

    Float_t z;
    Int_t i;
    Float_t p[5];    

    const Float_t kCrossSection = 0.094e-28;      // m^2
    const Float_t kFlux         = 1.e11 / 25.e-9; // protons/s
    Float_t pFlux[5] = {0.2, 0.2, 0.3, 0.3, 1.0};

    Int_t ncols = 0;
    if (fRunPeriod < 5) {
//
//  Ring 1   
// 

	for (i = 0; i < fGPASize; i++)
	{
	    ncols = fscanf(file, "%f %f %f %f %f %f", &z, &p[0], &p[1], &p[2] , &p[3], &p[4]);
	    if (ncols<0) break;

	    fG1[i] = p[fRunPeriod];
	    
	    if (i > 0) {
		fZ1[i] = fZ1[i-1] + z;
	    } else {
		fZ1[i] = 20.;
	    }
	}
//
// Ring 2
//
	for (i = 0; i < fGPASize; i++)
	{
	    ncols = fscanf(file, "%f %f %f %f %f %f", &z, &p[0], &p[1], &p[2] , &p[3], &p[4]);
	    if (ncols<0) break;

	    fG2[i] = p[fRunPeriod];
	    if (i > 0) {
		fZ2[i] = fZ2[i-1] + z;
	    } else {
		fZ2[i] = 20.;
	    }
	}
//
// Interaction rates
//
	for (i = 0; i < fGPASize; i++)  
	{
	    fG1[i] = fG1[i] * kCrossSection * pFlux[fRunPeriod] * kFlux; // 1/m/s 
	    fG2[i] = fG2[i] * kCrossSection * pFlux[fRunPeriod] * kFlux; // 1/m/s
	}

    } else {
	for (i = 0; i < fGPASize; i++) 
	{
	    ncols = fscanf(file, "%f %e %e %e %e %e", &z, &p[0], &p[1], &p[2], &p[3], &p[4]);
	    if (ncols<0) break;

	    z /= 1000.;
	    fG1[i] = p[4] * kCrossSection * kFlux;             // 1/m/s
	    // 1/3 of nominal intensity at startup
	    if (fRunPeriod ==  kLHCPR674Startup) fG1[i] /= 3.;
	    fZ1[i] = z;
	}
    }
    


    
//
//  Transform into interaction rates
//


    

    Float_t sum1 = 0.;
    Float_t sum2 = 0.;
    
    for (Int_t iz = 0; iz < 300; iz++) {
	Float_t zpos = 20. + iz * 1.;
	zpos *= 100;
	Float_t wgt1 = GasPressureWeight( zpos);
	Float_t wgt2 = GasPressureWeight(-zpos);
	sum1 += wgt1;
	sum2 += wgt2;
    }
    sum1/=250.;
    sum2/=250.;
    printf("\n %f %f \n \n", sum1, sum2);
    delete file;
}

//____________________________________________________________
void AliGenHaloProtvino::Generate()
{
// Generate from input file
 
  Float_t polar[3]= {0,0,0};
  Float_t origin[3];
  Float_t p[3], p0;
  Float_t tz, txy;
  Float_t amass;
  //
  Int_t ncols, nt;
  static Int_t nskip = 0;
  Int_t nread = 0;

  Float_t* zPrimary = new Float_t [fNpart];
  Int_t  * inuc     = new Int_t   [fNpart];
  Int_t  * ipart    = new Int_t   [fNpart];
  Float_t* wgt      = new Float_t [fNpart]; 
  Float_t* ekin     = new Float_t [fNpart];
  Float_t* vx       = new Float_t [fNpart];
  Float_t* vy       = new Float_t [fNpart];
  Float_t* tx       = new Float_t [fNpart];
  Float_t* ty       = new Float_t [fNpart];
  
  Float_t zVertexOld = -1.e10;
  Int_t   nInt       = 0;        // Counts number of interactions
  Float_t wwgt = 0.;
  
  while(1) {
//
// Load event into array
//
      ncols = fscanf(fFile,"%f %d %d %f %f %f %f %f %f",
		     &zPrimary[nread], &inuc[nread], &ipart[nread], &wgt[nread], 
		     &ekin[nread], &vx[nread], &vy[nread],
		     &tx[nread], &ty[nread]);
      
      if (ncols < 0) break;
// Skip fNskip events
      nskip++;
      if (fNpart !=-1 && nskip <= fNskip) continue;
// Count interactions
      if (zPrimary[nread] != zVertexOld) {
	  nInt++;
	  zVertexOld = zPrimary[nread];
      }
// Count tracks      
      nread++;
      if (fNpart !=-1 && nread >= fNpart) break;
  }
//
// Mean time between interactions
//

  Float_t dT = 0.;   // sec 
  if (nInt > 0) 
      dT = fTimePerEvent/nInt;   
  Float_t t  = 0;                    // sec
  
//
// Loop over primaries
//
  zVertexOld   = -1.e10;
  Double_t arg = 0.;
  
  for (Int_t nprim = 0; nprim < fNpart; nprim++) 
  {
      amass = TDatabasePDG::Instance()->GetParticle(ipart[nprim])->Mass();

      //
      // Momentum vector
      //
      p0=sqrt(ekin[nprim]*ekin[nprim] + 2.*amass*ekin[nprim]);
      
      txy=TMath::Sqrt(tx[nprim]*tx[nprim]+ty[nprim]*ty[nprim]);
      if (txy == 1.) {
	  tz=0;
      } else {
	  tz=-TMath::Sqrt(1.-txy);
      }
    
      p[0] = p0*tx[nprim];
      p[1] = p0*ty[nprim];
      p[2] =-p0*tz;
      
      origin[0] = vx[nprim];
      origin[1] = vy[nprim];
      origin[2] = -2196.5;

      //
      //
      // Particle weight

      Float_t originP[3] = {0., 0., 0.};
      originP[2] = zPrimary[nprim];
      
      Float_t pP[3] = {0., 0., 0.};
      Int_t ntP;
      
      if (fSide == -1) {
	  originP[2] = -zPrimary[nprim];
	  origin[2]  = -origin[2];
	  p[2]       = -p[2];
      }

      //
      // Time
      //
      if (zPrimary[nprim] != zVertexOld) {
	  while(arg==0.) arg = gRandom->Rndm();
	  t -= dT*TMath::Log(arg);              // (sec)   
	  zVertexOld = zPrimary[nprim];
      }

//    Get statistical weight according to local gas-pressure
//
      fParentWeight=wgt[nprim]*GasPressureWeight(zPrimary[nprim]);

      if (!fAnalog || gRandom->Rndm() < fParentWeight) {
//    Pass parent particle
//
	  PushTrack(0,-1,kProton,pP,originP,polar,t,kPNoProcess,ntP, fParentWeight);
	  KeepTrack(ntP);
	  PushTrack(fTrackIt,ntP,ipart[nprim],p,origin,polar,t,kPNoProcess,nt,fParentWeight);
      }
      //
      // Both sides are considered
      //

      if (fSide > 1) {
	  fParentWeight=wgt[nprim]*GasPressureWeight(-zPrimary[nprim]);
	  if (!fAnalog || gRandom->Rndm() < fParentWeight) {
	      origin[2]  = -origin[2];
	      originP[2] = -originP[2];
	      p[2]=-p[2];
	      PushTrack(0,-1,kProton,pP,originP,polar,t,kPNoProcess,ntP, fParentWeight);
	      KeepTrack(ntP);
	      PushTrack(fTrackIt,ntP,ipart[nprim],p,origin,polar,t,kPNoProcess,nt,fParentWeight);
	  }
      }
      wwgt += fParentWeight;
      
      SetHighWaterMark(nt);
  }
  delete [] zPrimary;
  delete [] inuc;    
  delete [] ipart;   
  delete [] wgt;     
  delete [] ekin;    
  delete [] vx;      
  delete [] vy;      
  delete [] tx;      
  delete [] ty;     
  printf("Total weight %f\n\n", wwgt);
  
}
 

Float_t AliGenHaloProtvino::GasPressureWeight(Float_t zPrimary)
{
//
// Return z-dependent gasspressure weight = interaction rate [1/m/s].
//
    Float_t weight = 0.;
    zPrimary /= 100.;        // m
    if (fRunPeriod < 5) {
	Float_t zAbs = TMath::Abs(zPrimary);
	if (zPrimary > 0.) 
	{
	    if (zAbs > fZ1[20]) {
		weight = 2.e4;
	    } else {
		for (Int_t i = 1; i < 21; i++) {
		    if (zAbs < fZ1[i]) {
			weight = fG1[i];
			break;
		    }
		}
	    }
	} else {
	    if (zAbs > fZ2[20]) {
		weight = 2.e4;
	    } else {
		for (Int_t i = 1; i < 21; i++) {
		    if (zAbs < fZ2[i]) {
		    weight = fG2[i];
		    break;
		    }
		}
	    }
	}
    } else {
	Int_t index = TMath::BinarySearch(fGPASize, fZ1, zPrimary);
	weight = fG1[index];
    }
    return weight;
}

void AliGenHaloProtvino::Draw(Option_t *)
{
// Draws the gas pressure distribution
    Float_t z[400];
    Float_t p[400];
    
    for (Int_t i = 0; i < 400; i++)
    {
	z[i] = -20000. + Float_t(i) * 100;
	p[i] = GasPressureWeight(z[i]);
    }
    
    TGraph*  gr = new TGraph(400, z, p);   
    TCanvas* c1 = new TCanvas("c1","Canvas 1",400,10,600,700);
    c1->cd();
    gr->Draw("AL");
}


/*
# Title:    README file for the sources of IR8 machine induced background
# Author:   Vadim Talanov <Vadim.Talanov@cern.ch>
# Modified: 12-12-2000 

0. Overview

	There are three files, named ring.one.beta.[01,10,50].m, which
	contain the lists of background particles, induced by proton losses
	upstream of IP8 in the LHC ring one, for the beta* values of 1, 10
	and 50 m, respectively.

1. File contents

	Each line in the files contains the coordinates of particle track
	crossing with the infinite plane, positioned at z=-1m, together with
	the physical properties of corresponding particle, namely:

	S  - S coordinate of the primary interaction vertex, cm;
	N  - type of the gas nuclei at interaction, 1 is H, 2 - C and 3 - O;
	I  - particle ID in PDG particle numbering scheme;
	W  - particle weight;
	E  - particle kinetic energy, GeV;
	X  - x coordinate of the crossing point, cm;
	Y  - y coordinate of the crossing point, cm;
	Dx - x direction cosine;
	Dy - y direction cosine.

2. Normalisation

	Each file is given per unity of linear density of proton inelastic
	interactions with the gas nuclei, [1 inelastic interaction/m].

# ~/vtalanov/public/README.mib: the end.

*/




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