TDPMjet/AliGenDPMjet.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
  17 // Generator using DPMJET as an external generator
  18 // The main DPMJET options are accessable for the user through this interface.
  19 // Uses the TDPMjet implementation of TGenerator.
  20
  21 #include <TDPMjet.h>
  22 #include "DPMcommon.h"
  23 #include <TRandom.h>
  24 #include <TArrayI.h>
  25 #include <TParticle.h>
  26 #include <TGraph.h>
  27 #include <TDatabasePDG.h>
  28 #include <TParticlePDG.h>
  29 #include <TParticleClassPDG.h>
  30 #include <TPDGCode.h>
  31 #include <TLorentzVector.h>
  32 #include <TClonesArray.h>
  33 #include "AliRunLoader.h"
  34 #include "AliGenDPMjet.h"
  35 #include "AliGenDPMjetEventHeader.h"
  36 #include "AliRun.h"
  37 #include "AliDpmJetRndm.h"
  38 #include "AliIonPDGCodes.h"
  39 #include "AliHeader.h"
  40 #include "AliStack.h"
  41 #include "AliMC.h"
  42 #include "AliLog.h"
  43
  44 ClassImp(AliGenDPMjet)
  45
  46 //______________________________________________________________________________
  47 AliGenDPMjet::AliGenDPMjet()
  48     :AliGenMC(),
  49      fBeamEn(0.),
  50      fMinImpactParam(0.),
  51      fMaxImpactParam(100.),
  52      fICentr(0),
  53      fSelectAll(0),
  54      fFlavor(0),
  55      fTrials(0),
  56      fSpectators(1),
  57      fSpecn(0),
  58      fSpecp(0),
  59      fDPMjet(0),
  60      fNoGammas(0),
  61      fLHC(0),
  62      fPi0Decay(1),
  63      fDecayAll(0),
  64      fGenImpPar(0.),
  65      fProcess(kDpmMb),
  66      fTriggerMultiplicity(0),
  67      fTriggerMultiplicityEta(0),
  68      fTriggerMultiplicityPtMin(0),
  69      fkTuneForDiff(0),
  70      fProcDiff(0),
  71      fFragmentation(kFALSE),
  72      fHeader(AliGenDPMjetEventHeader("DPMJET"))
  73
  74 {
  75 // Constructor
  76     fEnergyCMS = 5500.;
  77     AliDpmJetRndm::SetDpmJetRandom(GetRandom());
  78 }
  79
  80
  81 //______________________________________________________________________________
  82 AliGenDPMjet::AliGenDPMjet(Int_t npart)
  83     :AliGenMC(npart),
  84      fBeamEn(0.),
  85      fMinImpactParam(0.),
  86      fMaxImpactParam(100.),
  87      fICentr(0),
  88      fSelectAll(0),
  89      fFlavor(0),
  90      fTrials(0),
  91      fSpectators(1),
  92      fSpecn(0),
  93      fSpecp(0),
  94      fDPMjet(0),
  95      fNoGammas(0),
  96      fLHC(0),
  97      fPi0Decay(1),
  98      fDecayAll(0),
  99      fGenImpPar(0.),
 100      fProcess(kDpmMb),
 101      fTriggerMultiplicity(0),
 102      fTriggerMultiplicityEta(0),
 103      fTriggerMultiplicityPtMin(0),
 104      fkTuneForDiff(0),
 105      fProcDiff(0),
 106      fFragmentation(kFALSE),
 107      fHeader(AliGenDPMjetEventHeader("DPMJET"))
 108
 109
 110 {
 111 // Default PbPb collisions at 5. 5 TeV
 112 //
 113     fEnergyCMS = 5500.;
 114     fName = "DPMJET";
 115     fTitle= "Particle Generator using DPMJET";
 116     SetTarget();
 117     SetProjectile();
 118     fVertex.Set(3);
 119     AliDpmJetRndm::SetDpmJetRandom(GetRandom());
 120 }
 121
 122 AliGenDPMjet::AliGenDPMjet(const AliGenDPMjet &/*Dpmjet*/)
 123     :AliGenMC(),
 124      fBeamEn(0.),
 125      fMinImpactParam(0.),
 126      fMaxImpactParam(100.),
 127      fICentr(0),
 128      fSelectAll(0),
 129      fFlavor(0),
 130      fTrials(0),
 131      fSpectators(1),
 132      fSpecn(0),
 133      fSpecp(0),
 134      fDPMjet(0),
 135      fNoGammas(0),
 136      fLHC(0),
 137      fPi0Decay(1),
 138      fDecayAll(0),
 139      fGenImpPar(0.),
 140      fProcess(kDpmMb),
 141      fTriggerMultiplicity(0),
 142      fTriggerMultiplicityEta(0),
 143      fTriggerMultiplicityPtMin(0),
 144      fkTuneForDiff(0),
 145      fProcDiff(0),
 146      fFragmentation(kFALSE),
 147      fHeader(0x0)
 148
 149
 150 {
 151     // Dummy copy constructor
 152     fEnergyCMS = 5500.;
 153 }
 154
 155 //______________________________________________________________________________
 156 AliGenDPMjet::~AliGenDPMjet()
 157 {
 158 // Destructor
 159 }
 160 //______________________________________________________________________________
 161 void AliGenDPMjet::Init()
 162 {
 163 // Initialization
 164
 165     if(fEnergyCMS>0. && fBeamEn<0.1) fBeamEn = fEnergyCMS/2;
 166     SetMC(new TDPMjet(fProcess, fAProjectile, fZProjectile, fATarget, fZTarget,
 167                       fBeamEn,fEnergyCMS));
 168
 169     fDPMjet=(TDPMjet*) fMCEvGen;
 170     if (fAProjectile == 1 && TMath::Abs(fZProjectile == 1)) fDPMjet->SetfIdp(1);
 171
 172     // **** Flag to force central production
 173     // fICentr=1. central production forced
 174     // fICentr<0 && fICentr>-100 -> bmin = fMinImpactParam, bmax = fMaxImpactParam
 175     // fICentr<-99 -> fraction of x-sec. = XSFRAC
 176     // fICentr=-1. -> evaporation/fzc suppressed
 177     // fICentr<-1. -> evaporation/fzc allowed
 178     fDPMjet->SetfFCentr(fICentr);
 179
 180     fDPMjet->SetbRange(fMinImpactParam, fMaxImpactParam);
 181     fDPMjet->SetPi0Decay(fPi0Decay);
 182     fDPMjet->SetDecayAll(fDecayAll);
 183     fDPMjet->SetFragmentProd(fFragmentation);
 184
 185 //
 186 //  Initialize DPMjet
 187 //
 188     fDPMjet->Initialize();
 189 }
 190
 191
 192 //______________________________________________________________________________
 193 void AliGenDPMjet::Generate()
 194 {
 195 // Generate one event
 196
 197   Double_t polar[3]    =   {0,0,0};
 198   Double_t origin[3]   =   {0,0,0};
 199   Double_t p[4]        =   {0};
 200   Float_t tof;
 201
 202 //  converts from mm/c to s
 203   const Float_t kconv = 0.001/2.999792458e8;
 204   Int_t nt  = 0;
 205   Int_t jev = 0;
 206   Int_t kf, ks, imo;
 207   kf = 0;
 208   fTrials = 0;
 209   //  Set collision vertex position
 210   if (fVertexSmear == kPerEvent) Vertex();
 211
 212   while(1)
 213   {
 214 //    Generate one event
 215 // --------------------------------------------------------------------------
 216       fSpecn = 0;
 217       fSpecp = 0;
 218 // --------------------------------------------------------------------------
 219       fDPMjet->GenerateEvent();
 220
 221       fTrials++;
 222
 223       fDPMjet->ImportParticles(&fParticles,"All");
 224       if (fLHC) Boost();
 225
 226       // Temporaneo
 227       fGenImpPar = fDPMjet->GetBImpac();
 228
 229       if(TMath::Abs(fXingAngleY) > 1.e-10) BeamCrossAngle();
 230
 231       Int_t np = fParticles.GetEntriesFast();
 232       //
 233       // Multiplicity Trigger
 234       if (fTriggerMultiplicity > 0) {
 235         Int_t multiplicity = 0;
 236         for (Int_t i = 0; i < np; i++) {
 237           TParticle *  iparticle = (TParticle *) fParticles.At(i);
 238
 239           Int_t statusCode = iparticle->GetStatusCode();
 240
 241           // Initial state particle
 242           if (statusCode != 1)
 243             continue;
 244           // eta cut
 245           if (fTriggerMultiplicityEta > 0 && TMath::Abs(iparticle->Eta()) > fTriggerMultiplicityEta)
 246             continue;
 247           // pt cut
 248           if (iparticle->Pt() < fTriggerMultiplicityPtMin)
 249             continue;
 250
 251           TParticlePDG* pdgPart = iparticle->GetPDG();
 252           if (pdgPart && pdgPart->Charge() == 0)
 253             continue;
 254           ++multiplicity;
 255         }
 256         //
 257         //
 258         if (multiplicity < fTriggerMultiplicity) continue;
 259         Printf("Triggered on event with multiplicity of %d >= %d", multiplicity, fTriggerMultiplicity);
 260       }
 261
 262
 263       if(fkTuneForDiff && ( (TMath::Abs(fEnergyCMS - 900) < 1) || (TMath::Abs(fEnergyCMS - 2760) < 1) || (TMath::Abs(fEnergyCMS - 7000) < 1)) ) {
 264         if(!CheckDiffraction() ) continue;
 265       }
 266
 267
 268       Int_t nc = 0;
 269       if (np == 0) continue;
 270
 271       Int_t i;
 272       Int_t* newPos     = new Int_t[np];
 273       Int_t* pSelected  = new Int_t[np];
 274
 275       for (i = 0; i<np; i++) {
 276           newPos[i]    = i;
 277           pSelected[i] = 0;
 278       }
 279
 280 //      First select parent particles
 281
 282       for (i = 0; i<np; i++) {
 283           TParticle *iparticle = (TParticle *) fParticles.At(i);
 284
 285 // Is this a parent particle ?
 286
 287           if (Stable(iparticle)) continue;
 288
 289           Bool_t  selected             =  kTRUE;
 290           Bool_t  hasSelectedDaughters =  kFALSE;
 291
 292           kf = iparticle->GetPdgCode();
 293           if (kf == 92 || kf == 99999) continue;
 294           ks = iparticle->GetStatusCode();
 295 // No initial state partons
 296           if (ks==21) continue;
 297           if (!fSelectAll) selected = KinematicSelection(iparticle, 0) &&
 298                                SelectFlavor(kf);
 299
 300
 301           hasSelectedDaughters = DaughtersSelection(iparticle);
 302
 303
 304 // Put particle on the stack if it is either selected or
 305 // it is the mother of at least one seleted particle
 306
 307           if (selected || hasSelectedDaughters) {
 308               nc++;
 309               pSelected[i] = 1;
 310           } // selected
 311       } // particle loop parents
 312
 313 // Now select the final state particles
 314
 315
 316       for (i=0; i<np; i++) {
 317           TParticle *iparticle = (TParticle *) fParticles.At(i);
 318
 319 // Is this a final state particle ?
 320
 321           if (!Stable(iparticle)) continue;
 322
 323           Bool_t  selected =  kTRUE;
 324           kf = iparticle->GetPdgCode();
 325           ks = iparticle->GetStatusCode();
 326
 327 // --------------------------------------------------------------------------
 328 // Count spectator neutrons and protons (ks == 13, 14)
 329           if(ks == 13 || ks == 14){
 330               if(kf == kNeutron) fSpecn += 1;
 331               if(kf == kProton)  fSpecp += 1;
 332           }
 333 // --------------------------------------------------------------------------
 334
 335           if (!fSelectAll) {
 336               selected = KinematicSelection(iparticle,0)&&SelectFlavor(kf);
 337               if (!fSpectators && selected) selected = (ks == 13 || ks == 14);
 338           }
 339
 340 // Put particle on the stack if selected
 341
 342           if (selected) {
 343               nc++;
 344               pSelected[i] = 1;
 345           } // selected
 346       } // particle loop final state
 347
 348 // Write particles to stack
 349
 350       for (i = 0; i<np; i++) {
 351           TParticle *  iparticle = (TParticle *) fParticles.At(i);
 352           Bool_t  hasMother   = (iparticle->GetFirstMother()>=0);
 353           if (pSelected[i]) {
 354
 355               kf   = iparticle->GetPdgCode();
 356               ks   = iparticle->GetStatusCode();
 357
 358               p[0] = iparticle->Px();
 359               p[1] = iparticle->Py();
 360               p[2] = iparticle->Pz();
 361               p[3] = iparticle->Energy();
 362               origin[0] = fVertex[0]+iparticle->Vx()/10; // [cm]
 363               origin[1] = fVertex[1]+iparticle->Vy()/10; // [cm]
 364               origin[2] = fVertex[2]+iparticle->Vz()/10; // [cm]
 365
 366               tof = fTime + kconv*iparticle->T();
 367
 368               imo = -1;
 369               TParticle* mother = 0;
 370               if (hasMother) {
 371                   imo = iparticle->GetFirstMother();
 372                   mother = (TParticle *) fParticles.At(imo);
 373                   imo = (mother->GetPdgCode() != 92 && mother->GetPdgCode() != 99999) ? newPos[imo] : -1;
 374               } // if has mother
 375
 376
 377
 378               Bool_t tFlag = (fTrackIt && (ks==1 || ks==-1));
 379               //printf(" AliGemDPMJet->PushTrack: kf %d  ks %d  flag %d\n",kf,ks,tFlag);
 380               if(kf>10000 && (ks==-1 || ks==1000 || ks==1001)) kf += 1000000000;
 381               PushTrack(tFlag, imo, kf,
 382                         p[0], p[1], p[2], p[3],
 383                         origin[0], origin[1], origin[2], tof,
 384                         polar[0], polar[1], polar[2],
 385                         kPNoProcess, nt, 1., ks);
 386               KeepTrack(nt);
 387               newPos[i] = nt;
 388           } // if selected
 389       } // particle loop
 390       delete[] newPos;
 391       delete[] pSelected;
 392       if (nc>0) {
 393           jev += nc;
 394           if (jev >= fNpart || fNpart == -1) {
 395               break;
 396           }
 397       }
 398   } // event loop
 399   MakeHeader();
 400   SetHighWaterMark(nt);
 401 }
 402
 403 //______________________________________________________________________________
 404 Bool_t AliGenDPMjet::DaughtersSelection(TParticle* iparticle)
 405 {
 406 //
 407 // Looks recursively if one of the daughters has been selected
 408 //
 409 //    printf("\n Consider daughters %d:",iparticle->GetPdgCode());
 410     Int_t imin = -1;
 411     Int_t imax = -1;
 412     Int_t i;
 413     Bool_t hasDaughters = (iparticle->GetFirstDaughter() >=0);
 414     Bool_t selected = kFALSE;
 415     if (hasDaughters) {
 416         imin = iparticle->GetFirstDaughter();
 417         imax = iparticle->GetLastDaughter();
 418         for (i = imin; i <= imax; i++){
 419             TParticle *  jparticle = (TParticle *) fParticles.At(i);
 420             Int_t ip = jparticle->GetPdgCode();
 421             if (KinematicSelection(jparticle,0)&&SelectFlavor(ip)) {
 422                 selected=kTRUE; break;
 423             }
 424             if (DaughtersSelection(jparticle)) {selected=kTRUE; break; }
 425         }
 426     } else {
 427         return kFALSE;
 428     }
 429     return selected;
 430 }
 431
 432
 433
 434 //______________________________________________________________________________
 435 Bool_t AliGenDPMjet::SelectFlavor(Int_t pid)
 436 {
 437 // Select flavor of particle
 438 // 0: all
 439 // 4: charm and beauty
 440 // 5: beauty
 441     Bool_t res = 0;
 442
 443     if (fFlavor == 0) {
 444         res = kTRUE;
 445     } else {
 446         Int_t ifl = TMath::Abs(pid/100);
 447         if (ifl > 10) ifl/=10;
 448         res = (fFlavor == ifl);
 449     }
 450 //
 451 //  This part if gamma writing is inhibited
 452     if (fNoGammas)
 453         res = res && (pid != kGamma && pid != kPi0);
 454 //
 455     return res;
 456 }
 457
 458 //______________________________________________________________________________
 459 Bool_t AliGenDPMjet::Stable(TParticle*  particle)
 460 {
 461 // Return true for a stable particle
 462 //
 463     int st = particle->GetStatusCode();
 464     if(st == 1 || st == -1) return kTRUE;
 465     else return kFALSE;
 466
 467 }
 468
 469 //______________________________________________________________________________
 470 void AliGenDPMjet::MakeHeader()
 471 {
 472 // Builds the event header, to be called after each event
 473     fHeader.SetNProduced(fDPMjet->GetNumStablePc());
 474     fHeader.SetImpactParameter(fDPMjet->GetBImpac());
 475     fHeader.SetTotalEnergy(fDPMjet->GetTotEnergy());
 476     fHeader.SetParticipants(fDPMjet->GetProjParticipants(),
 477                             fDPMjet->GetTargParticipants());
 478
 479     fHeader.SetCollisions(DTGLCP.ncp, DTGLCP.nct,
 480         fDPMjet->GetProjWounded(),fDPMjet->GetTargWounded());
 481
 482     if(fProcDiff>0)  fHeader.SetProcessType(fProcDiff);
 483     else fHeader.SetProcessType(fDPMjet->GetProcessCode());
 484
 485     // Bookkeeping for kinematic bias
 486     fHeader.SetTrials(fTrials);
 487     // Event Vertex
 488     fHeader.SetPrimaryVertex(fVertex);
 489     fHeader.SetInteractionTime(fTime);
 490     fHeader.SetNDiffractive(POEVT1.nsd1, POEVT1.nsd2, POEVT1.ndd);
 491 //    gAlice->SetGenEventHeader(fHeader);
 492     AddHeader(&fHeader);
 493     fCollisionGeometry = &fHeader;
 494 }
 495
 496 //______________________________________________________________________________
 497 /*void AliGenDPMjet::AddHeader(AliGenEventHeader* fHeader)
 498 {
 499     // Add fHeader to container or runloader
 500     if (fContainer) {
 501         fContainer->AddHeader(fHeader);
 502     } else {
 503         AliRunLoader::Instance()->GetHeader()->SetGenEventHeader(fHeader);
 504     }
 505 }*/
 506
 507
 508 //______________________________________________________________________________
 509 AliGenDPMjet& AliGenDPMjet::operator=(const  AliGenDPMjet& /*rhs*/)
 510 {
 511 // Assignment operator
 512     return *this;
 513 }
 514
 515
 516 //______________________________________________________________________________
 517 void AliGenDPMjet::FinishRun()
 518 {
 519     // Print run statistics
 520     fDPMjet->Dt_Dtuout();
 521 }
 522
 523
 524 //______________________________________________________________________________
 525 Bool_t AliGenDPMjet::CheckDiffraction()
 526 {
 527
 528   //  printf("AAA\n");
 529
 530    Int_t np = fParticles.GetEntriesFast();
 531
 532    Int_t iPart1=-1;
 533    Int_t iPart2=-1;
 534
 535    Double_t y1 = 1e10;
 536    Double_t y2 = -1e10;
 537
 538   const Int_t kNstable=20;
 539   const Int_t pdgStable[20] = {
 540     22,             // Photon
 541     11,             // Electron
 542     12,             // Electron Neutrino
 543     13,             // Muon
 544     14,             // Muon Neutrino
 545     15,             // Tau
 546     16,             // Tau Neutrino
 547     211,            // Pion
 548     321,            // Kaon
 549     311,            // K0
 550     130,            // K0s
 551     310,            // K0l
 552     2212,           // Proton
 553     2112,           // Neutron
 554     3122,           // Lambda_0
 555     3112,           // Sigma Minus
 556     3222,           // Sigma Plus
 557     3312,           // Xsi Minus
 558     3322,           // Xsi0
 559     3334            // Omega
 560   };
 561
 562      for (Int_t i = 0; i < np; i++) {
 563         TParticle *  part = (TParticle *) fParticles.At(i);
 564
 565         Int_t statusCode = part->GetStatusCode();
 566
 567         // Initial state particle
 568         if (statusCode != 1)
 569           continue;
 570
 571         Int_t pdg = TMath::Abs(part->GetPdgCode());
 572         Bool_t isStable = kFALSE;
 573         for (Int_t i1 = 0; i1 < kNstable; i1++) {
 574           if (pdg == pdgStable[i1]) {
 575             isStable = kTRUE;
 576             break;
 577           }
 578         }
 579         if(!isStable)
 580           continue;
 581
 582         Double_t y = part->Y();
 583
 584         if (y < y1)
 585           {
 586             y1 = y;
 587             iPart1 = i;
 588           }
 589         if (y > y2)
 590         {
 591           y2 = y;
 592           iPart2 = i;
 593         }
 594      }
 595
 596      if(iPart1<0 || iPart2<0) return kFALSE;
 597
 598      y1=TMath::Abs(y1);
 599      y2=TMath::Abs(y2);
 600
 601      TParticle *  part1 = (TParticle *) fParticles.At(iPart1);
 602      TParticle *  part2 = (TParticle *) fParticles.At(iPart2);
 603
 604      Int_t pdg1 = part1->GetPdgCode();
 605      Int_t pdg2 = part2->GetPdgCode();
 606
 607
 608      Int_t iPart = -1;
 609      if (pdg1 == 2212 && pdg2 == 2212)
 610        {
 611          if(y1 > y2)
 612            iPart = iPart1;
 613          else if(y1 < y2)
 614            iPart = iPart2;
 615          else {
 616            iPart = iPart1;
 617            if((AliDpmJetRndm::GetDpmJetRandom())->Uniform(0.,1.)>0.5) iPart = iPart2;
 618          }
 619        }
 620      else if (pdg1 == 2212)
 621        iPart = iPart1;
 622      else if (pdg2 == 2212)
 623        iPart = iPart2;
 624
 625
 626
 627
 628
 629      Double_t M=-1.;
 630      if(iPart>0) {
 631        TParticle *  part = (TParticle *) fParticles.At(iPart);
 632        Double_t E= part->Energy();
 633        Double_t P= part->P();
 634        Double_t M2 = (fEnergyCMS-E-P)*(fEnergyCMS-E+P);
 635        if(M2<0)  return kFALSE;
 636        M= TMath::Sqrt(M2);
 637      }
 638
 639      Double_t Mmin, Mmax, wSD, wDD, wND;
 640
 641      if(!GetWeightsDiffraction(M, Mmin, Mmax, wSD, wDD, wND)) return kFALSE;
 642      if(M>-1 && M<Mmin) return kFALSE;
 643      if(M>Mmax) M=-1;
 644
 645
 646  Int_t procType=fDPMjet->GetProcessCode();//fPythia->GetMSTI(1);
 647  if(procType== 4) return kFALSE;
 648  Int_t proc0=2;
 649  if(procType== 7) proc0=1;
 650  if(procType== 5 || procType== 6) proc0=0;
 651
 652
 653  Int_t proc=2;
 654  if(M>0) proc=0;
 655  else if(proc0==1) proc=1;
 656
 657  if(proc==1 && (AliDpmJetRndm::GetDpmJetRandom())->Uniform(0.,1.) > wDD) return kFALSE;
 658  if(proc==2 && (AliDpmJetRndm::GetDpmJetRandom())->Uniform(0.,1.) > wND) return kFALSE;
 659
 660
 661
 662     if(proc!=0) {
 663       if(proc0!=0) fProcDiff = procType;
 664       else       fProcDiff = 1;
 665       return kTRUE;
 666     }
 667
 668     if(wSD<0)  AliError("wSD<0 ! \n");
 669
 670     if((AliDpmJetRndm::GetDpmJetRandom())->Uniform(0.,1.)> wSD) return kFALSE;
 671
 672     //    printf("iPart = %d\n", iPart);
 673
 674     if(iPart==iPart1) fProcDiff=5;
 675     else if(iPart==iPart2) fProcDiff=6;
 676     else {
 677       printf("EROOR:  iPart!=iPart1 && iPart!=iPart2\n");
 678
 679     }
 680
 681     return kTRUE;
 682 }
 683
 684 // -------------------------------------------------------
 685 void AliGenDPMjet::SetIonPDGCodes()
 686 {
 687    // Defining PDG codes for the ions
 688    AliIonPDGCodes *pdgcodes = new AliIonPDGCodes();
 689    pdgcodes->AddParticlesToPdgDataBase();
 690 }
 691
 692 // -------------------------------------------------------
 693 Bool_t AliGenDPMjet::GetWeightsDiffraction(Double_t M, Double_t &Mmin, Double_t &Mmax,
 694                                                        Double_t &wSD, Double_t &wDD, Double_t &wND)
 695 {
 696
 697   if(TMath::Abs(fEnergyCMS-900)<1 ){
 698 const Int_t nbin=400;
 699 Double_t bin[]={
 700 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
 701 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
 702 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
 703 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
 704 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
 705 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
 706 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
 707 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
 708 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
 709 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
 710 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
 711 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
 712 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
 713 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
 714 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
 715 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
 716 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
 717 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
 718 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
 719 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
 720 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
 721 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
 722 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
 723 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
 724 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
 725 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
 726 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
 727 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
 728 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
 729 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
 730 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
 731 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
 732 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
 733 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
 734 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
 735 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
 736 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
 737 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
 738 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
 739 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
 740 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
 741 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
 742 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
 743 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
 744 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
 745 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
 746 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
 747 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
 748 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
 749 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
 750 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
 751 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
 752 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
 753 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
 754 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
 755 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
 756 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
 757 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
 758 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
 759 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
 760 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
 761 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
 762 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
 763 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
 764 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
 765 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
 766 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
 767 Double_t w[]={
 768 1.000000, 0.389079, 0.326612, 0.356363, 0.310726, 0.264037,
 769 0.253810, 0.224655, 0.207990, 0.198149, 0.186803, 0.180209,
 770 0.178928, 0.161772, 0.162314, 0.158409, 0.148379, 0.143557,
 771 0.140357, 0.136272, 0.136957, 0.136606, 0.129376, 0.127671,
 772 0.128406, 0.132239, 0.119311, 0.130227, 0.130233, 0.123312,
 773 0.115257, 0.120105, 0.119733, 0.117116, 0.110096, 0.117509,
 774 0.109149, 0.114372, 0.100322, 0.106227, 0.108696, 0.110352,
 775 0.106748, 0.101475, 0.101837, 0.098585, 0.094433, 0.100148,
 776 0.096505, 0.100155, 0.103526, 0.098161, 0.093647, 0.100121,
 777 0.097426, 0.093414, 0.090241, 0.097228, 0.098098, 0.098340,
 778 0.096779, 0.099265, 0.095101, 0.090342, 0.097469, 0.090663,
 779 0.093422, 0.093103, 0.089648, 0.096923, 0.088079, 0.089360,
 780 0.097772, 0.092932, 0.093554, 0.085784, 0.096395, 0.096304,
 781 0.090183, 0.089255, 0.090265, 0.086262, 0.087044, 0.088965,
 782 0.085825, 0.095073, 0.085933, 0.091005, 0.091707, 0.092428,
 783 0.091689, 0.091224, 0.095256, 0.083407, 0.087983, 0.090320,
 784 0.081580, 0.088077, 0.083478, 0.091309, 0.083734, 0.089906,
 785 0.100155, 0.092728, 0.086542, 0.078091, 0.085261, 0.094302,
 786 0.078063, 0.090070, 0.086566, 0.095020, 0.082124, 0.084791,
 787 0.090624, 0.090236, 0.086344, 0.085706, 0.085913, 0.083107,
 788 0.092557, 0.081144, 0.091254, 0.096139, 0.079458, 0.088124,
 789 0.088777, 0.076652, 0.085168, 0.082326, 0.080435, 0.085022,
 790 0.085693, 0.093957, 0.088134, 0.082347, 0.082139, 0.090980,
 791 0.078397, 0.093595, 0.091562, 0.081789, 0.085026, 0.080868,
 792 0.083455, 0.084590, 0.084124, 0.081486, 0.087242, 0.080928,
 793 0.096392, 0.078324, 0.092093, 0.074267, 0.082925, 0.090717,
 794 0.081467, 0.080435, 0.078825, 0.075843, 0.088940, 0.081451,
 795 0.096047, 0.088102, 0.091710, 0.095208, 0.086160, 0.091305,
 796 0.085406, 0.094432, 0.078227, 0.099870, 0.094140, 0.080130,
 797 0.090707, 0.083268, 0.082222, 0.088767, 0.084477, 0.088069,
 798 0.089382, 0.086164, 0.089123, 0.092799, 0.076710, 0.090727,
 799 0.077097, 0.099905, 0.089733, 0.088101, 0.093705, 0.089215,
 800 0.085110, 0.086032, 0.083719, 0.079693, 0.088116, 0.090519,
 801 0.091150, 0.090855, 0.090547, 0.077773, 0.081914, 0.080864,
 802 0.082935, 0.092952, 0.079390, 0.080255, 0.091123, 0.091331,
 803 0.095160, 0.089343, 0.083353, 0.087445, 0.094036, 0.084719,
 804 0.099665, 0.085104, 0.100912, 0.085958, 0.083972, 0.098284,
 805 0.078318, 0.082042, 0.088007, 0.085469, 0.089984, 0.081181,
 806 0.098850, 0.086409, 0.089070, 0.081579, 0.089622, 0.089396,
 807 0.093922, 0.089472, 0.090806, 0.075034, 0.090346, 0.083871,
 808 0.086931, 0.089207, 0.094425, 0.087830, 0.079537, 0.097316,
 809 0.096513, 0.092264, 0.082211, 0.083841, 0.081861, 0.092591,
 810 0.077785, 0.079646, 0.093721, 0.090735, 0.086910, 0.075837,
 811 0.090729, 0.092800, 0.098704, 0.092441, 0.086404, 0.089344,
 812 0.083650, 0.082569, 0.085753, 0.091430, 0.085460, 0.095210,
 813 0.083226, 0.094882, 0.084856, 0.093350, 0.084579, 0.096744,
 814 0.099211, 0.085142, 0.085499, 0.083395, 0.088352, 0.110418,
 815 0.093788, 0.096140, 0.082758, 0.097145, 0.089170, 0.090720,
 816 0.084708, 0.095927, 0.090330, 0.074239, 0.089532, 0.090514,
 817 0.089823, 0.089709, 0.101840, 0.082676, 0.087157, 0.079221,
 818 0.096460, 0.108192, 0.088904, 0.084510, 0.096624, 0.099242,
 819 0.094470, 0.092473, 0.091745, 0.090439, 0.095316, 0.087963,
 820 0.096400, 0.089479, 0.094880, 0.085964, 0.092775, 0.083200,
 821 0.095133, 0.090079, 0.088828, 0.087600, 0.094123, 0.089135,
 822 0.082617, 0.085109, 0.114091, 0.096331, 0.108465, 0.080318,
 823 0.106576, 0.089671, 0.092023, 0.090722, 0.086603, 0.091788,
 824 0.098375, 0.082712, 0.102681, 0.098869, 0.089051, 0.109972,
 825 0.080440, 0.097860, 0.093074, 0.097028, 0.107826, 0.117152,
 826 0.090968, 0.096790, 0.096725, 0.094641, 0.089535, 0.092561,
 827 0.095828, 0.084320, 0.089942, 0.105476, 0.087495, 0.089805,
 828 0.092238, 0.094141, 0.104390, 0.082958, 0.097449, 0.099594,
 829 0.092640, 0.097332, 0.093223, 0.100183, 0.092511, 0.087035,
 830 0.088741, 0.097856, 0.116682, 0.091732, 0.097753, 0.094283,
 831 0.094235, 0.082016, 0.098370, 0.085676, 0.104529, 0.087319,
 832 0.090289, 0.105314, 0.103634, 0.101261, 0.092764, 0.098217,
 833 0.098939, 0.096071, 0.096071, 0.094027, 0.092713, 0.089542,
 834 0.112293, 0.112293, 0.089531, 0.087752};
 835 wDD = 0.188264;
 836 wND = 0.102252;
 837 wSD = -1;
 838
 839     Mmin = bin[0];
 840     Mmax = bin[nbin];
 841     if(M<Mmin || M>Mmax) return kTRUE;
 842
 843     Int_t ibin=nbin-1;
 844     for(Int_t i=1; i<=nbin; i++)
 845       if(M<=bin[i]) {
 846         ibin=i-1;
 847         //      printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);
 848         break;
 849       }
 850     wSD=w[ibin];
 851     return kTRUE;
 852
 853   }
 854   else if(TMath::Abs(fEnergyCMS-2760)<1 ){
 855
 856 const Int_t nbin=400;
 857 Double_t bin[]={
 858 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
 859 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
 860 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
 861 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
 862 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
 863 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
 864 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
 865 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
 866 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
 867 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
 868 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
 869 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
 870 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
 871 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
 872 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
 873 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
 874 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
 875 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
 876 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
 877 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
 878 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
 879 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
 880 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
 881 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
 882 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
 883 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
 884 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
 885 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
 886 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
 887 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
 888 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
 889 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
 890 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
 891 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
 892 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
 893 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
 894 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
 895 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
 896 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
 897 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
 898 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
 899 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
 900 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
 901 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
 902 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
 903 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
 904 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
 905 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
 906 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
 907 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
 908 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
 909 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
 910 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
 911 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
 912 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
 913 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
 914 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
 915 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
 916 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
 917 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
 918 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
 919 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
 920 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
 921 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
 922 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
 923 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
 924 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
 925 Double_t w[]={
 926 1.000000, 0.402402, 0.347976, 0.386866, 0.304413, 0.275746,
 927 0.256941, 0.250439, 0.228486, 0.219225, 0.203368, 0.196124,
 928 0.180551, 0.169230, 0.159021, 0.157136, 0.154743, 0.160522,
 929 0.153993, 0.143448, 0.148062, 0.136337, 0.130853, 0.127798,
 930 0.137576, 0.137530, 0.124175, 0.120225, 0.120774, 0.123584,
 931 0.121883, 0.114522, 0.118338, 0.119587, 0.106017, 0.125802,
 932 0.106256, 0.101585, 0.097737, 0.100817, 0.102279, 0.110935,
 933 0.101429, 0.107629, 0.105802, 0.095688, 0.096617, 0.100985,
 934 0.101521, 0.100450, 0.111836, 0.088828, 0.102264, 0.101012,
 935 0.097673, 0.097634, 0.091720, 0.096167, 0.084537, 0.104443,
 936 0.091966, 0.097204, 0.089456, 0.088273, 0.098220, 0.095151,
 937 0.088201, 0.085570, 0.088431, 0.077625, 0.088314, 0.084582,
 938 0.091545, 0.083774, 0.086183, 0.089195, 0.086255, 0.077167,
 939 0.096118, 0.089881, 0.086206, 0.086968, 0.079827, 0.082002,
 940 0.083081, 0.076587, 0.081820, 0.092382, 0.082964, 0.080901,
 941 0.093512, 0.082907, 0.086691, 0.093639, 0.089010, 0.082857,
 942 0.090795, 0.082403, 0.078602, 0.087284, 0.079638, 0.082905,
 943 0.077189, 0.084329, 0.075900, 0.082559, 0.084210, 0.070053,
 944 0.068453, 0.083369, 0.077659, 0.073286, 0.075396, 0.073946,
 945 0.083293, 0.076235, 0.082184, 0.080104, 0.073755, 0.077369,
 946 0.068281, 0.084593, 0.074923, 0.069467, 0.072094, 0.077973,
 947 0.067034, 0.080314, 0.087441, 0.071728, 0.072302, 0.076783,
 948 0.079931, 0.075557, 0.076318, 0.079029, 0.075126, 0.076859,
 949 0.080253, 0.074344, 0.083387, 0.087553, 0.084437, 0.076322,
 950 0.077748, 0.087559, 0.075649, 0.076615, 0.077098, 0.074559,
 951 0.075540, 0.085883, 0.079269, 0.079866, 0.068922, 0.069770,
 952 0.068024, 0.077783, 0.086620, 0.086769, 0.068433, 0.078059,
 953 0.084463, 0.072790, 0.073889, 0.079090, 0.081759, 0.069576,
 954 0.069160, 0.071695, 0.078569, 0.075727, 0.088055, 0.082395,
 955 0.079915, 0.075150, 0.073580, 0.069968, 0.071141, 0.065823,
 956 0.077270, 0.076229, 0.071735, 0.079271, 0.078209, 0.069503,
 957 0.064732, 0.076312, 0.088579, 0.087271, 0.080566, 0.073527,
 958 0.071376, 0.068785, 0.071812, 0.072680, 0.067252, 0.059543,
 959 0.072818, 0.071856, 0.083427, 0.080664, 0.077315, 0.066060,
 960 0.078772, 0.065442, 0.067255, 0.077676, 0.068306, 0.071578,
 961 0.059728, 0.073526, 0.078932, 0.070316, 0.067620, 0.067323,
 962 0.078316, 0.080366, 0.091194, 0.066739, 0.076238, 0.070382,
 963 0.077118, 0.073822, 0.072830, 0.062947, 0.077375, 0.076284,
 964 0.070951, 0.086127, 0.071074, 0.076621, 0.070502, 0.070897,
 965 0.058764, 0.065322, 0.073620, 0.062051, 0.074548, 0.083770,
 966 0.072502, 0.066823, 0.074192, 0.081570, 0.067340, 0.081896,
 967 0.080680, 0.072835, 0.069675, 0.078881, 0.072826, 0.069616,
 968 0.082744, 0.069138, 0.073333, 0.080881, 0.069091, 0.064266,
 969 0.062849, 0.072846, 0.067462, 0.067229, 0.071319, 0.071099,
 970 0.069381, 0.073799, 0.076106, 0.075006, 0.073074, 0.092589,
 971 0.074277, 0.061324, 0.069788, 0.070310, 0.063160, 0.073077,
 972 0.075447, 0.081013, 0.074102, 0.076565, 0.060192, 0.058326,
 973 0.071508, 0.072974, 0.065098, 0.060432, 0.062077, 0.075644,
 974 0.075417, 0.066947, 0.066744, 0.065111, 0.075468, 0.078185,
 975 0.073216, 0.066502, 0.079582, 0.065405, 0.069828, 0.072165,
 976 0.071112, 0.064656, 0.080396, 0.066505, 0.063419, 0.066128,
 977 0.080616, 0.071463, 0.064867, 0.074782, 0.063103, 0.068179,
 978 0.064911, 0.073029, 0.075746, 0.062264, 0.099126, 0.067890,
 979 0.074940, 0.070090, 0.073583, 0.056904, 0.063745, 0.063780,
 980 0.076876, 0.066519, 0.067164, 0.068668, 0.075065, 0.070983,
 981 0.075751, 0.066317, 0.077520, 0.073139, 0.073961, 0.085492,
 982 0.074629, 0.059555, 0.081026, 0.059274, 0.055448, 0.084269,
 983 0.066222, 0.078874, 0.062021, 0.075156, 0.090478, 0.066265,
 984 0.067845, 0.076653, 0.086671, 0.082362, 0.067246, 0.075975,
 985 0.072594, 0.068570, 0.095034, 0.065644, 0.070945, 0.066190,
 986 0.075906, 0.055154, 0.075575, 0.058078, 0.071004, 0.070850,
 987 0.064374, 0.058251, 0.069430, 0.072293, 0.065484, 0.084197,
 988 0.090119, 0.091619, 0.067620, 0.079462, 0.063025, 0.068128,
 989 0.056927, 0.076351, 0.073869, 0.061597, 0.083741, 0.063762,
 990 0.064489, 0.074269, 0.068832, 0.058648, 0.069536, 0.074824,
 991 0.081009, 0.073389, 0.076792, 0.084855, 0.075280, 0.061473,
 992 0.063840, 0.062891, 0.071328, 0.065250};
 993 wDD = 0.077845;
 994 wND = 0.066355;
 995 wSD = -1;
 996
 997     Mmin = bin[0];
 998     Mmax = bin[nbin];
 999     if(M<Mmin || M>Mmax) return kTRUE;
1000
1001     Int_t ibin=nbin-1;
1002     for(Int_t i=1; i<=nbin; i++)
1003       if(M<=bin[i]) {
1004         ibin=i-1;
1005         //      printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);
1006         break;
1007       }
1008     wSD=w[ibin];
1009     return kTRUE;
1010
1011   }
1012   else if(TMath::Abs(fEnergyCMS-7000)<1 ){
1013
1014 const Int_t nbin=400;
1015 Double_t bin[]={
1016 1.080000, 1.577300, 2.074600, 2.571900, 3.069200, 3.566500,
1017 4.063800, 4.561100, 5.058400, 5.555700, 6.053000, 6.550300,
1018 7.047600, 7.544900, 8.042200, 8.539500, 9.036800, 9.534100,
1019 10.031400, 10.528700, 11.026000, 11.523300, 12.020600, 12.517900,
1020 13.015200, 13.512500, 14.009800, 14.507100, 15.004400, 15.501700,
1021 15.999000, 16.496300, 16.993600, 17.490900, 17.988200, 18.485500,
1022 18.982800, 19.480100, 19.977400, 20.474700, 20.972000, 21.469300,
1023 21.966600, 22.463900, 22.961200, 23.458500, 23.955800, 24.453100,
1024 24.950400, 25.447700, 25.945000, 26.442300, 26.939600, 27.436900,
1025 27.934200, 28.431500, 28.928800, 29.426100, 29.923400, 30.420700,
1026 30.918000, 31.415300, 31.912600, 32.409900, 32.907200, 33.404500,
1027 33.901800, 34.399100, 34.896400, 35.393700, 35.891000, 36.388300,
1028 36.885600, 37.382900, 37.880200, 38.377500, 38.874800, 39.372100,
1029 39.869400, 40.366700, 40.864000, 41.361300, 41.858600, 42.355900,
1030 42.853200, 43.350500, 43.847800, 44.345100, 44.842400, 45.339700,
1031 45.837000, 46.334300, 46.831600, 47.328900, 47.826200, 48.323500,
1032 48.820800, 49.318100, 49.815400, 50.312700, 50.810000, 51.307300,
1033 51.804600, 52.301900, 52.799200, 53.296500, 53.793800, 54.291100,
1034 54.788400, 55.285700, 55.783000, 56.280300, 56.777600, 57.274900,
1035 57.772200, 58.269500, 58.766800, 59.264100, 59.761400, 60.258700,
1036 60.756000, 61.253300, 61.750600, 62.247900, 62.745200, 63.242500,
1037 63.739800, 64.237100, 64.734400, 65.231700, 65.729000, 66.226300,
1038 66.723600, 67.220900, 67.718200, 68.215500, 68.712800, 69.210100,
1039 69.707400, 70.204700, 70.702000, 71.199300, 71.696600, 72.193900,
1040 72.691200, 73.188500, 73.685800, 74.183100, 74.680400, 75.177700,
1041 75.675000, 76.172300, 76.669600, 77.166900, 77.664200, 78.161500,
1042 78.658800, 79.156100, 79.653400, 80.150700, 80.648000, 81.145300,
1043 81.642600, 82.139900, 82.637200, 83.134500, 83.631800, 84.129100,
1044 84.626400, 85.123700, 85.621000, 86.118300, 86.615600, 87.112900,
1045 87.610200, 88.107500, 88.604800, 89.102100, 89.599400, 90.096700,
1046 90.594000, 91.091300, 91.588600, 92.085900, 92.583200, 93.080500,
1047 93.577800, 94.075100, 94.572400, 95.069700, 95.567000, 96.064300,
1048 96.561600, 97.058900, 97.556200, 98.053500, 98.550800, 99.048100,
1049 99.545400, 100.042700, 100.540000, 101.037300, 101.534600, 102.031900,
1050 102.529200, 103.026500, 103.523800, 104.021100, 104.518400, 105.015700,
1051 105.513000, 106.010300, 106.507600, 107.004900, 107.502200, 107.999500,
1052 108.496800, 108.994100, 109.491400, 109.988700, 110.486000, 110.983300,
1053 111.480600, 111.977900, 112.475200, 112.972500, 113.469800, 113.967100,
1054 114.464400, 114.961700, 115.459000, 115.956300, 116.453600, 116.950900,
1055 117.448200, 117.945500, 118.442800, 118.940100, 119.437400, 119.934700,
1056 120.432000, 120.929300, 121.426600, 121.923900, 122.421200, 122.918500,
1057 123.415800, 123.913100, 124.410400, 124.907700, 125.405000, 125.902300,
1058 126.399600, 126.896900, 127.394200, 127.891500, 128.388800, 128.886100,
1059 129.383400, 129.880700, 130.378000, 130.875300, 131.372600, 131.869900,
1060 132.367200, 132.864500, 133.361800, 133.859100, 134.356400, 134.853700,
1061 135.351000, 135.848300, 136.345600, 136.842900, 137.340200, 137.837500,
1062 138.334800, 138.832100, 139.329400, 139.826700, 140.324000, 140.821300,
1063 141.318600, 141.815900, 142.313200, 142.810500, 143.307800, 143.805100,
1064 144.302400, 144.799700, 145.297000, 145.794300, 146.291600, 146.788900,
1065 147.286200, 147.783500, 148.280800, 148.778100, 149.275400, 149.772700,
1066 150.270000, 150.767300, 151.264600, 151.761900, 152.259200, 152.756500,
1067 153.253800, 153.751100, 154.248400, 154.745700, 155.243000, 155.740300,
1068 156.237600, 156.734900, 157.232200, 157.729500, 158.226800, 158.724100,
1069 159.221400, 159.718700, 160.216000, 160.713300, 161.210600, 161.707900,
1070 162.205200, 162.702500, 163.199800, 163.697100, 164.194400, 164.691700,
1071 165.189000, 165.686300, 166.183600, 166.680900, 167.178200, 167.675500,
1072 168.172800, 168.670100, 169.167400, 169.664700, 170.162000, 170.659300,
1073 171.156600, 171.653900, 172.151200, 172.648500, 173.145800, 173.643100,
1074 174.140400, 174.637700, 175.135000, 175.632300, 176.129600, 176.626900,
1075 177.124200, 177.621500, 178.118800, 178.616100, 179.113400, 179.610700,
1076 180.108000, 180.605300, 181.102600, 181.599900, 182.097200, 182.594500,
1077 183.091800, 183.589100, 184.086400, 184.583700, 185.081000, 185.578300,
1078 186.075600, 186.572900, 187.070200, 187.567500, 188.064800, 188.562100,
1079 189.059400, 189.556700, 190.054000, 190.551300, 191.048600, 191.545900,
1080 192.043200, 192.540500, 193.037800, 193.535100, 194.032400, 194.529700,
1081 195.027000, 195.524300, 196.021600, 196.518900, 197.016200, 197.513500,
1082 198.010800, 198.508100, 199.005400, 199.502700, 200.000000};
1083 Double_t w[]={
1084 1.000000, 0.526293, 0.446686, 0.437789, 0.366854, 0.333320,
1085 0.291931, 0.266464, 0.253870, 0.248706, 0.232788, 0.220736,
1086 0.209886, 0.202741, 0.188617, 0.182767, 0.178748, 0.169039,
1087 0.175911, 0.169098, 0.171256, 0.146728, 0.144543, 0.159470,
1088 0.153171, 0.151883, 0.144693, 0.136307, 0.140226, 0.135388,
1089 0.141317, 0.151121, 0.131209, 0.144039, 0.124688, 0.128020,
1090 0.119122, 0.121868, 0.122332, 0.119561, 0.115744, 0.102504,
1091 0.114726, 0.109518, 0.115418, 0.114860, 0.110026, 0.107693,
1092 0.103005, 0.115985, 0.108629, 0.105937, 0.101056, 0.101228,
1093 0.113305, 0.110302, 0.104696, 0.107447, 0.099095, 0.107378,
1094 0.103090, 0.111384, 0.090821, 0.109806, 0.093639, 0.096655,
1095 0.110416, 0.104446, 0.098530, 0.095105, 0.106641, 0.101804,
1096 0.091798, 0.094775, 0.104130, 0.088436, 0.119692, 0.099160,
1097 0.086163, 0.089250, 0.101898, 0.091786, 0.087346, 0.095312,
1098 0.107186, 0.085671, 0.093283, 0.091992, 0.114654, 0.086172,
1099 0.084559, 0.097862, 0.079417, 0.094737, 0.089678, 0.105265,
1100 0.099092, 0.087901, 0.086828, 0.084756, 0.077769, 0.084254,
1101 0.083262, 0.092935, 0.088858, 0.093377, 0.083569, 0.084771,
1102 0.084957, 0.084033, 0.095007, 0.076690, 0.087726, 0.097913,
1103 0.084119, 0.076361, 0.084172, 0.091068, 0.094564, 0.091406,
1104 0.088405, 0.090355, 0.085388, 0.084586, 0.085752, 0.084298,
1105 0.091294, 0.081491, 0.075024, 0.082543, 0.088048, 0.074621,
1106 0.078514, 0.077249, 0.078378, 0.092945, 0.075859, 0.077602,
1107 0.074592, 0.078181, 0.081527, 0.080865, 0.078854, 0.078237,
1108 0.085506, 0.089823, 0.072738, 0.078430, 0.077162, 0.081495,
1109 0.090878, 0.099417, 0.100966, 0.087960, 0.085556, 0.081661,
1110 0.078066, 0.089059, 0.073564, 0.092704, 0.073148, 0.098057,
1111 0.067906, 0.079452, 0.090567, 0.082454, 0.077851, 0.079694,
1112 0.091272, 0.078628, 0.096906, 0.098779, 0.088906, 0.070174,
1113 0.083822, 0.084241, 0.093237, 0.071062, 0.075771, 0.096405,
1114 0.098441, 0.086007, 0.069599, 0.078400, 0.083481, 0.081054,
1115 0.085552, 0.069582, 0.071336, 0.073207, 0.087913, 0.080125,
1116 0.075189, 0.067217, 0.073509, 0.099694, 0.080781, 0.073943,
1117 0.062822, 0.076611, 0.087672, 0.077632, 0.063959, 0.077754,
1118 0.084651, 0.097348, 0.063909, 0.055053, 0.087616, 0.074428,
1119 0.101165, 0.078446, 0.070683, 0.071162, 0.091516, 0.069641,
1120 0.070969, 0.069778, 0.083318, 0.066703, 0.074027, 0.067589,
1121 0.070620, 0.081307, 0.094406, 0.076188, 0.061663, 0.077561,
1122 0.076159, 0.071851, 0.074417, 0.076136, 0.069073, 0.075450,
1123 0.070218, 0.089709, 0.079974, 0.082077, 0.076979, 0.070497,
1124 0.076296, 0.087347, 0.073508, 0.088073, 0.067186, 0.059898,
1125 0.085495, 0.100152, 0.079468, 0.093783, 0.082680, 0.077272,
1126 0.081995, 0.068240, 0.071881, 0.073737, 0.072346, 0.079034,
1127 0.078721, 0.067518, 0.068196, 0.081738, 0.082814, 0.082480,
1128 0.069072, 0.066853, 0.081523, 0.073588, 0.082308, 0.091641,
1129 0.072747, 0.084317, 0.072190, 0.079372, 0.067424, 0.077450,
1130 0.062343, 0.090708, 0.065470, 0.086588, 0.071172, 0.066705,
1131 0.070670, 0.070407, 0.096111, 0.066738, 0.081249, 0.072882,
1132 0.075144, 0.060331, 0.074589, 0.076968, 0.085913, 0.072561,
1133 0.064645, 0.078742, 0.075670, 0.065984, 0.080932, 0.069898,
1134 0.065303, 0.096856, 0.057690, 0.065720, 0.066545, 0.068521,
1135 0.068278, 0.069245, 0.086643, 0.063401, 0.070933, 0.070752,
1136 0.066978, 0.058891, 0.070073, 0.078031, 0.082691, 0.101414,
1137 0.075814, 0.072790, 0.057622, 0.093002, 0.084660, 0.079216,
1138 0.070371, 0.070141, 0.076944, 0.067285, 0.078016, 0.077807,
1139 0.066668, 0.066459, 0.059962, 0.062774, 0.083450, 0.064554,
1140 0.067887, 0.064165, 0.072782, 0.067285, 0.052710, 0.096824,
1141 0.071931, 0.064190, 0.074442, 0.082647, 0.055797, 0.078632,
1142 0.061116, 0.063092, 0.049131, 0.074517, 0.069915, 0.079021,
1143 0.088656, 0.101045, 0.090432, 0.076460, 0.067451, 0.071381,
1144 0.069790, 0.077330, 0.077115, 0.084403, 0.065138, 0.056436,
1145 0.088024, 0.069893, 0.055985, 0.089655, 0.062911, 0.075311,
1146 0.086679, 0.093645, 0.068632, 0.064498, 0.057677, 0.081391,
1147 0.062781, 0.075467, 0.061314, 0.073394, 0.084462, 0.068470,
1148 0.071267, 0.060556, 0.072487, 0.063785, 0.079164, 0.070406,
1149 0.073394, 0.063168, 0.066968, 0.064125, 0.063954, 0.072504,
1150 0.058948, 0.057740, 0.083383, 0.067262};
1151 wDD = 0.098024;
1152 wND = 0.053452;
1153 wSD = -1;
1154
1155     Mmin = bin[0];
1156     Mmax = bin[nbin];
1157     if(M<Mmin || M>Mmax) return kTRUE;
1158
1159     Int_t ibin=nbin-1;
1160     for(Int_t i=1; i<=nbin; i++)
1161       if(M<=bin[i]) {
1162         ibin=i-1;
1163         //      printf("Mi> %f && Mi< %f\n", bin[i-1], bin[i]);
1164         break;
1165       }
1166     wSD=w[ibin];
1167     return kTRUE;
1168
1169   }
1170
1171   return kFALSE;
1172 }
1173
1174
1175 //______________________________________________________________________________