EVGEN/AliGenParam.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id$ */
  17
  18 // Class to generate particles from using paramtrized pT and y distributions.
  19 // Distributions are obtained from pointer to object of type AliGenLib.
  20 // (For example AliGenMUONlib)
  21 // Decays are performed using Pythia.
  22 // andreas.morsch@cern.ch
  23
  24 #include <TCanvas.h>
  25 #include <TClonesArray.h>
  26 #include <TDatabasePDG.h>
  27 #include <TF1.h>
  28 #include <TH1F.h>
  29 #include <TLorentzVector.h>
  30 #include <TMath.h>
  31 #include <TParticle.h>
  32 #include <TParticlePDG.h>
  33 #include <TROOT.h>
  34 #include <TVirtualMC.h>
  35
  36 #include "AliDecayer.h"
  37 #include "AliGenMUONlib.h"
  38 #include "AliGenParam.h"
  39 #include "AliMC.h"
  40 #include "AliRun.h"
  41 #include "AliGenEventHeader.h"
  42
  43 ClassImp(AliGenParam)
  44
  45 //------------------------------------------------------------
  46
  47 //Begin_Html
  48 /*
  49     <img src="picts/AliGenParam.gif">
  50 */
  51 //End_Html
  52
  53 //____________________________________________________________
  54 AliGenParam::AliGenParam()
  55 : fPtParaFunc(0),
  56         fYParaFunc(0),
  57         fIpParaFunc(0),
  58   fV2ParaFunc(0),
  59         fPtPara(0),
  60         fYPara(0),
  61   fV2Para(0),
  62   fdNdPhi(0),
  63         fParam(0),
  64         fdNdy0(0.),
  65         fYWgt(0.),
  66         fPtWgt(0.),
  67         fBias(0.),
  68         fTrials(0),
  69         fDeltaPt(0.01),
  70         fSelectAll(kFALSE),
  71         fDecayer(0)
  72 {
  73   // Default constructor
  74 }
  75 //____________________________________________________________
  76 AliGenParam::AliGenParam(Int_t npart, const AliGenLib * Library,  Int_t param, const char* tname)
  77     :AliGenMC(npart),
  78      fPtParaFunc(Library->GetPt(param, tname)),
  79      fYParaFunc (Library->GetY (param, tname)),
  80      fIpParaFunc(Library->GetIp(param, tname)),
  81    fV2ParaFunc(Library->GetV2(param, tname)),
  82      fPtPara(0),
  83      fYPara(0),
  84    fV2Para(0),
  85    fdNdPhi(0),
  86      fParam(param),
  87      fdNdy0(0.),
  88      fYWgt(0.),
  89      fPtWgt(0.),
  90      fBias(0.),
  91      fTrials(0),
  92      fDeltaPt(0.01),
  93      fSelectAll(kFALSE),
  94      fDecayer(0)
  95 {
  96   // Constructor using number of particles parameterisation id and library
  97     fName = "Param";
  98     fTitle= "Particle Generator using pT and y parameterisation";
  99     fAnalog = kAnalog;
 100     SetForceDecay();
 101 }
 102 //____________________________________________________________
 103 AliGenParam::AliGenParam(Int_t npart, Int_t param, const char* tname, const char* name):
 104     AliGenMC(npart),
 105     fPtParaFunc(0),
 106     fYParaFunc (0),
 107     fIpParaFunc(0),
 108   fV2ParaFunc(0),
 109     fPtPara(0),
 110     fYPara(0),
 111   fV2Para(0),
 112   fdNdPhi(0),
 113     fParam(param),
 114     fdNdy0(0.),
 115     fYWgt(0.),
 116     fPtWgt(0.),
 117     fBias(0.),
 118     fTrials(0),
 119     fDeltaPt(0.01),
 120     fSelectAll(kFALSE),
 121     fDecayer(0)
 122 {
 123   // Constructor using parameterisation id and number of particles
 124   //
 125     fName = name;
 126     fTitle= "Particle Generator using pT and y parameterisation";
 127
 128     AliGenLib* pLibrary = new AliGenMUONlib();
 129     fPtParaFunc = pLibrary->GetPt(param, tname);
 130     fYParaFunc  = pLibrary->GetY (param, tname);
 131     fIpParaFunc = pLibrary->GetIp(param, tname);
 132   fV2ParaFunc = pLibrary->GetV2(param, tname);
 133
 134     fAnalog = kAnalog;
 135     fChildSelect.Set(5);
 136     for (Int_t i=0; i<5; i++) fChildSelect[i]=0;
 137     SetForceDecay();
 138     SetCutOnChild();
 139     SetChildMomentumRange();
 140     SetChildPtRange();
 141     SetChildPhiRange();
 142     SetChildThetaRange();
 143 }
 144 //____________________________________________________________
 145
 146 AliGenParam::AliGenParam(Int_t npart, Int_t param,
 147                          Double_t (*PtPara) (const Double_t*, const Double_t*),
 148                          Double_t (*YPara ) (const Double_t* ,const Double_t*),
 149                          Double_t (*V2Para) (const Double_t* ,const Double_t*),
 150                          Int_t    (*IpPara) (TRandom *))
 151     :AliGenMC(npart),
 152
 153      fPtParaFunc(PtPara),
 154      fYParaFunc(YPara),
 155      fIpParaFunc(IpPara),
 156    fV2ParaFunc(V2Para),
 157      fPtPara(0),
 158      fYPara(0),
 159    fV2Para(0),
 160    fdNdPhi(0),
 161      fParam(param),
 162      fdNdy0(0.),
 163      fYWgt(0.),
 164      fPtWgt(0.),
 165      fBias(0.),
 166      fTrials(0),
 167      fDeltaPt(0.01),
 168      fSelectAll(kFALSE),
 169      fDecayer(0)
 170 {
 171   // Constructor
 172   // Gines Martinez 1/10/99
 173     fName = "Param";
 174     fTitle= "Particle Generator using pT and y parameterisation";
 175
 176     fAnalog = kAnalog;
 177     fChildSelect.Set(5);
 178     for (Int_t i=0; i<5; i++) fChildSelect[i]=0;
 179     SetForceDecay();
 180     SetCutOnChild();
 181     SetChildMomentumRange();
 182     SetChildPtRange();
 183     SetChildPhiRange();
 184     SetChildThetaRange();
 185 }
 186
 187 //____________________________________________________________
 188 AliGenParam::~AliGenParam()
 189 {
 190   // Destructor
 191     delete  fPtPara;
 192     delete  fYPara;
 193   delete  fV2Para;
 194   delete  fdNdPhi;
 195 }
 196
 197 //____________________________________________________________
 198 void AliGenParam::Init()
 199 {
 200   // Initialisation
 201
 202     if (gMC) fDecayer = gMC->GetDecayer();
 203   //Begin_Html
 204   /*
 205     <img src="picts/AliGenParam.gif">
 206   */
 207   //End_Html
 208     char name[256];
 209     snprintf(name, 256, "pt-parameterisation for %s", GetName());
 210
 211     if (fPtPara) fPtPara->Delete();
 212     fPtPara = new TF1(name, fPtParaFunc, fPtMin, fPtMax,0);
 213     gROOT->GetListOfFunctions()->Remove(fPtPara);
 214   //  Set representation precision to 10 MeV
 215     Int_t npx= Int_t((fPtMax - fPtMin) / fDeltaPt);
 216
 217     fPtPara->SetNpx(npx);
 218
 219     snprintf(name, 256, "y-parameterisation  for %s", GetName());
 220     if (fYPara) fYPara->Delete();
 221     fYPara  = new TF1(name, fYParaFunc, fYMin, fYMax, 0);
 222     gROOT->GetListOfFunctions()->Remove(fYPara);
 223
 224   snprintf(name, 256, "v2-parameterisation for %s", GetName());
 225   if (fV2Para) fV2Para->Delete();
 226   fV2Para  = new TF1(name, fV2ParaFunc, fPtMin, fPtMax, 0);
 227   // fV2Para  = new TF1(name, "2*[0]/(1+TMath::Exp([1]*([2]-x)))-[0]", fPtMin, fPtMax);
 228   // fV2Para->SetParameter(0, 0.236910);
 229   // fV2Para->SetParameter(1, 1.71122);
 230   // fV2Para->SetParameter(2, 0.0827617);
 231   //gROOT->GetListOfFunctions()->Remove(fV2Para);  //TR: necessary?
 232
 233   snprintf(name, 256, "dNdPhi for %s", GetName());
 234   if (fdNdPhi) fdNdPhi->Delete();
 235   fdNdPhi  = new TF1(name, "1+2*[0]*TMath::Cos(2*(x-[1]))", fPhiMin, fPhiMax);
 236   //gROOT->GetListOfFunctions()->Remove(fdNdPhi);  //TR: necessary?
 237
 238     snprintf(name, 256, "pt-for-%s", GetName());
 239     TF1 ptPara(name ,fPtParaFunc, 0, 15, 0);
 240     snprintf(name, 256, "y-for-%s", GetName());
 241     TF1 yPara(name, fYParaFunc, -6, 6, 0);
 242
 243   //
 244   // dN/dy| y=0
 245     Double_t y1=0;
 246     Double_t y2=0;
 247
 248     fdNdy0=fYParaFunc(&y1,&y2);
 249   //
 250   // Integral over generation region
 251 #if ROOT_VERSION_CODE < ROOT_VERSION(5,99,0)
 252     Float_t intYS  = yPara.Integral(fYMin, fYMax,(Double_t*) 0x0,1.e-6);
 253     Float_t intPt0 = ptPara.Integral(0,15,(Double_t *) 0x0,1.e-6);
 254     Float_t intPtS = ptPara.Integral(fPtMin,fPtMax,(Double_t*) 0x0,1.e-6);
 255 #else
 256     Float_t intYS  = yPara.Integral(fYMin, fYMax,1.e-6);
 257     Float_t intPt0 = ptPara.Integral(0,15,1.e-6);
 258     Float_t intPtS = ptPara.Integral(fPtMin,fPtMax,1.e-6);
 259 #endif
 260   Float_t phiWgt=(fPhiMax-fPhiMin)/2./TMath::Pi();    //TR: should probably be done differently in case of anisotropic phi...
 261     if (fAnalog == kAnalog) {
 262         fYWgt  = intYS/fdNdy0;
 263         fPtWgt = intPtS/intPt0;
 264         fParentWeight = fYWgt*fPtWgt*phiWgt/fNpart;
 265     } else {
 266         fYWgt = intYS/fdNdy0;
 267         fPtWgt = (fPtMax-fPtMin)/intPt0;
 268         fParentWeight = fYWgt*fPtWgt*phiWgt/fNpart;
 269     }
 270   //
 271   // particle decay related initialization
 272     fDecayer->SetForceDecay(fForceDecay);
 273     fDecayer->Init();
 274
 275   //
 276     AliGenMC::Init();
 277 }
 278
 279 //____________________________________________________________
 280 void AliGenParam::Generate()
 281 {
 282   //
 283   // Generate 1 event (see Generate(Int_t ntimes) for details
 284   //
 285   GenerateN(1);
 286 }
 287 //____________________________________________________________
 288 void AliGenParam::GenerateN(Int_t ntimes)
 289 {
 290   //
 291   // Generate ntimes*'npart' light and heavy mesons (J/Psi, upsilon or phi, Pion,
 292   // Kaons, Etas, Omegas) and Baryons (proton, antiprotons, neutrons and
 293   // antineutrons in the the desired theta, phi and momentum windows;
 294   // Gaussian smearing on the vertex is done if selected.
 295   // The decay of heavy mesons is done using lujet,
 296   //    and the childern particle are tracked by GEANT
 297   // However, light mesons are directly tracked by GEANT
 298   // setting fForceDecay = nodecay (SetForceDecay(nodecay))
 299   //
 300   //
 301   //  Reinitialize decayer
 302   fDecayer->SetForceDecay(fForceDecay);
 303   fDecayer->Init();
 304
 305   //
 306   Float_t polar[3]= {0,0,0};  // Polarisation of the parent particle (for GEANT tracking)
 307   Float_t origin0[3];         // Origin of the generated parent particle (for GEANT tracking)
 308   Float_t time0;              // Time0 of the generated parent particle
 309   Float_t pt, pl, ptot;       // Transverse, logitudinal and total momenta of the parent particle
 310   Float_t phi, theta;         // Phi and theta spherical angles of the parent particle momentum
 311   Float_t p[3], pc[3],
 312           och[3];             // Momentum, polarisation and origin of the children particles from lujet
 313   Double_t ty, xmt;
 314   Int_t nt, i, j;
 315   Float_t  wgtp, wgtch;
 316   Double_t dummy;
 317   static TClonesArray *particles;
 318   //
 319   if(!particles) particles = new TClonesArray("TParticle",1000);
 320
 321   TDatabasePDG *pDataBase = TDatabasePDG::Instance();
 322   //
 323   Float_t random[6];
 324
 325   // Calculating vertex position per event
 326   for (j=0;j<3;j++) origin0[j]=fOrigin[j];
 327   time0 = fTimeOrigin;
 328   if(fVertexSmear==kPerEvent) {
 329       Vertex();
 330       for (j=0;j<3;j++) origin0[j]=fVertex[j];
 331       time0 = fTime;
 332   }
 333
 334   Int_t ipa=0;
 335
 336   // Generating fNpart particles
 337   fNprimaries = 0;
 338
 339   Int_t nGen = fNpart*ntimes;
 340   while (ipa<nGen) {
 341       while(1) {
 342       //
 343       // particle type
 344           Int_t iPart = fIpParaFunc(fRandom);
 345           fChildWeight=(fDecayer->GetPartialBranchingRatio(iPart))*fParentWeight;
 346           TParticlePDG *particle = pDataBase->GetParticle(iPart);
 347           Float_t am = particle->Mass();
 348
 349           Rndm(random,2);
 350       //
 351       // y
 352           ty = TMath::TanH(fYPara->GetRandom());
 353       //
 354       // pT
 355           if (fAnalog == kAnalog) {
 356               pt=fPtPara->GetRandom();
 357               wgtp=fParentWeight;
 358               wgtch=fChildWeight;
 359           } else {
 360               pt=fPtMin+random[1]*(fPtMax-fPtMin);
 361               Double_t ptd=pt;
 362               wgtp=fParentWeight*fPtParaFunc(& ptd, &dummy);
 363               wgtch=fChildWeight*fPtParaFunc(& ptd, &dummy);
 364           }
 365           xmt=sqrt(pt*pt+am*am);
 366           if (TMath::Abs(ty)==1.) {
 367               ty=0.;
 368               Fatal("AliGenParam",
 369                     "Division by 0: Please check you rapidity range !");
 370           }
 371       //
 372       // phi
 373       // if(!ipa)
 374       //        phi=fEvPlane; //align first particle of each event with event plane
 375       // else{
 376         double v2 = fV2Para->Eval(pt);
 377         fdNdPhi->SetParameter(0,v2);
 378         fdNdPhi->SetParameter(1,fEvPlane);
 379         phi=fdNdPhi->GetRandom();
 380       // }
 381
 382           pl=xmt*ty/sqrt((1.-ty)*(1.+ty));
 383           theta=TMath::ATan2(pt,pl);
 384       // Cut on theta
 385           if(theta<fThetaMin || theta>fThetaMax) continue;
 386           ptot=TMath::Sqrt(pt*pt+pl*pl);
 387       // Cut on momentum
 388           if(ptot<fPMin || ptot>fPMax) continue;
 389       //
 390           p[0]=pt*TMath::Cos(phi);
 391           p[1]=pt*TMath::Sin(phi);
 392           p[2]=pl;
 393           if(fVertexSmear==kPerTrack) {
 394               Rndm(random,6);
 395               for (j=0;j<3;j++) {
 396                   origin0[j]=
 397                       fOrigin[j]+fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
 398                       TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
 399               }
 400               Rndm(random,2);
 401               time0 = fTimeOrigin + fOsigma[2]/TMath::Ccgs()*
 402                 TMath::Cos(2*random[0]*TMath::Pi())*
 403                 TMath::Sqrt(-2*TMath::Log(random[1]));
 404           }
 405
 406       // Looking at fForceDecay :
 407       // if fForceDecay != none Primary particle decays using
 408       // AliPythia and children are tracked by GEANT
 409       //
 410       // if fForceDecay == none Primary particle is tracked by GEANT
 411       // (In the latest, make sure that GEANT actually does all the decays you want)
 412       //
 413           Bool_t decayed = kFALSE;
 414
 415
 416           if (fForceDecay != kNoDecay) {
 417         // Using lujet to decay particle
 418               Float_t energy=TMath::Sqrt(ptot*ptot+am*am);
 419               TLorentzVector pmom(p[0], p[1], p[2], energy);
 420               fDecayer->Decay(iPart,&pmom);
 421         //
 422         // select decay particles
 423               Int_t np=fDecayer->ImportParticles(particles);
 424
 425               //  Selecting  GeometryAcceptance for particles fPdgCodeParticleforAcceptanceCut;
 426               if (fGeometryAcceptance)
 427                 if (!CheckAcceptanceGeometry(np,particles)) continue;
 428               Int_t ncsel=0;
 429               Int_t* pFlag      = new Int_t[np];
 430               Int_t* pParent    = new Int_t[np];
 431               Int_t* pSelected  = new Int_t[np];
 432               Int_t* trackIt    = new Int_t[np];
 433
 434               for (i=0; i<np; i++) {
 435                   pFlag[i]     =  0;
 436                   pSelected[i] =  0;
 437                   pParent[i]   = -1;
 438               }
 439
 440               if (np >1) {
 441                   decayed = kTRUE;
 442                   TParticle* iparticle =  0;
 443                   Int_t ipF, ipL;
 444                   for (i = 1; i<np ; i++) {
 445                       trackIt[i] = 1;
 446                       iparticle = (TParticle *) particles->At(i);
 447                       Int_t kf = iparticle->GetPdgCode();
 448                       Int_t ks = iparticle->GetStatusCode();
 449             // flagged particle
 450
 451                       if (pFlag[i] == 1) {
 452                           ipF = iparticle->GetFirstDaughter();
 453                           ipL = iparticle->GetLastDaughter();
 454                           if (ipF > 0) for (j=ipF-1; j<ipL; j++) pFlag[j]=1;
 455                           continue;
 456                       }
 457
 458             // flag decay products of particles with long life-time (c tau > .3 mum)
 459
 460                       if (ks != 1) {
 461               //                          TParticlePDG *particle = pDataBase->GetParticle(kf);
 462
 463                           Double_t lifeTime = fDecayer->GetLifetime(kf);
 464               //                          Double_t mass     = particle->Mass();
 465               //                          Double_t width    = particle->Width();
 466                           if (lifeTime > (Double_t) fMaxLifeTime) {
 467                               ipF = iparticle->GetFirstDaughter();
 468                               ipL = iparticle->GetLastDaughter();
 469                               if (ipF > 0) for (j=ipF-1; j<ipL; j++) pFlag[j]=1;
 470                           } else{
 471                               trackIt[i]     = 0;
 472                               pSelected[i]   = 1;
 473                           }
 474                       } // ks==1 ?
 475             //
 476             // children
 477
 478                       if ((ChildSelected(TMath::Abs(kf)) || fForceDecay == kAll || fSelectAll) && trackIt[i])
 479                       {
 480                           if (fCutOnChild) {
 481                               pc[0]=iparticle->Px();
 482                               pc[1]=iparticle->Py();
 483                               pc[2]=iparticle->Pz();
 484                               Bool_t  childok = KinematicSelection(iparticle, 1);
 485                               if(childok) {
 486                                   pSelected[i]  = 1;
 487                                   ncsel++;
 488                               } else {
 489                                   ncsel=-1;
 490                                   break;
 491                               } // child kine cuts
 492                           } else {
 493                               pSelected[i]  = 1;
 494                               ncsel++;
 495                           } // if child selection
 496                       } // select muon
 497                   } // decay particle loop
 498               } // if decay products
 499
 500               Int_t iparent;
 501               if ((fCutOnChild && ncsel >0) || !fCutOnChild){
 502                   ipa++;
 503           //
 504           // Parent
 505
 506
 507                   PushTrack(0, -1, iPart, p, origin0, polar, time0, kPPrimary, nt, wgtp, ((decayed)? 11 : 1));
 508                   pParent[0] = nt;
 509                   KeepTrack(nt);
 510                   fNprimaries++;
 511
 512           //
 513           // Decay Products
 514           //
 515                   for (i = 1; i < np; i++) {
 516                       if (pSelected[i]) {
 517                           TParticle* iparticle = (TParticle *) particles->At(i);
 518                           Int_t kf   = iparticle->GetPdgCode();
 519                           Int_t ksc  = iparticle->GetStatusCode();
 520                           Int_t jpa  = iparticle->GetFirstMother()-1;
 521
 522                           och[0] = origin0[0]+iparticle->Vx();
 523                           och[1] = origin0[1]+iparticle->Vy();
 524                           och[2] = origin0[2]+iparticle->Vz();
 525                           pc[0]  = iparticle->Px();
 526                           pc[1]  = iparticle->Py();
 527                           pc[2]  = iparticle->Pz();
 528
 529                           if (jpa > -1) {
 530                               iparent = pParent[jpa];
 531                           } else {
 532                               iparent = -1;
 533                           }
 534
 535                           PushTrack(fTrackIt * trackIt[i], iparent, kf,
 536                                            pc, och, polar,
 537                                            time0 + iparticle->T(), kPDecay, nt, wgtch, ksc);
 538                           pParent[i] = nt;
 539                           KeepTrack(nt);
 540                           fNprimaries++;
 541                       } // Selected
 542                   } // Particle loop
 543               }  // Decays by Lujet
 544               particles->Clear();
 545               if (pFlag)      delete[] pFlag;
 546               if (pParent)    delete[] pParent;
 547               if (pSelected)  delete[] pSelected;
 548               if (trackIt)    delete[] trackIt;
 549           } // kinematic selection
 550           else  // nodecay option, so parent will be tracked by GEANT (pions, kaons, eta, omegas, baryons)
 551           {
 552             gAlice->GetMCApp()->
 553                 PushTrack(fTrackIt,-1,iPart,p,origin0,polar,time0,kPPrimary,nt,wgtp, 1);
 554             ipa++;
 555             fNprimaries++;
 556           }
 557           break;
 558     } // while
 559   } // event loop
 560
 561   SetHighWaterMark(nt);
 562
 563   AliGenEventHeader* header = new AliGenEventHeader("PARAM");
 564   header->SetPrimaryVertex(fVertex);
 565   header->SetInteractionTime(fTime);
 566   header->SetNProduced(fNprimaries);
 567   AddHeader(header);
 568 }
 569 //____________________________________________________________________________________
 570 Float_t AliGenParam::GetRelativeArea(Float_t ptMin, Float_t ptMax, Float_t yMin, Float_t yMax, Float_t phiMin, Float_t phiMax)
 571 {
 572   //
 573   // Normalisation for selected kinematic region
 574   //
 575 #if ROOT_VERSION_CODE < ROOT_VERSION(5,99,0)
 576   Float_t ratio =
 577     fPtPara->Integral(ptMin,ptMax,(Double_t *)0,1.e-6) / fPtPara->Integral( fPtPara->GetXmin(), fPtPara->GetXmax(),(Double_t *)0,1.e-6) *
 578     fYPara->Integral(yMin,yMax,(Double_t *)0,1.e-6)/fYPara->Integral(fYPara->GetXmin(),fYPara->GetXmax(),(Double_t *)0,1.e-6)   *
 579     (phiMax-phiMin)/360.;
 580 #else
 581   Float_t ratio =
 582     fPtPara->Integral(ptMin,ptMax,1.e-6) / fPtPara->Integral( fPtPara->GetXmin(), fPtPara->GetXmax(),1.e-6) *
 583     fYPara->Integral(yMin,yMax,1.e-6)/fYPara->Integral(fYPara->GetXmin(),fYPara->GetXmax(),1.e-6)   *
 584     (phiMax-phiMin)/360.;
 585 #endif
 586   return TMath::Abs(ratio);
 587 }
 588
 589 //____________________________________________________________________________________
 590
 591 void AliGenParam::Draw( const char * /*opt*/)
 592 {
 593     //
 594     // Draw the pT and y Distributions
 595     //
 596      TCanvas *c0 = new TCanvas("c0","Canvas 0",400,10,600,700);
 597      c0->Divide(2,1);
 598      c0->cd(1);
 599      fPtPara->Draw();
 600      fPtPara->GetHistogram()->SetXTitle("p_{T} (GeV)");
 601      c0->cd(2);
 602      fYPara->Draw();
 603      fYPara->GetHistogram()->SetXTitle("y");
 604 }
 605
 606
 607
 608