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