As part of the efforts of the combined PID subgroup we are continuing to develop...

[u/mrichter/AliRoot.git] / PYTHIA6 / AliDecayerPythia.cxx

diff --git a/PYTHIA6/AliDecayerPythia.cxx b/PYTHIA6/AliDecayerPythia.cxx
index 8f30e291a11e6afb8997b3071dd63aed3bfb0281..3e628a71b349e9116f637329e0082a54340306ff 100644 (file)
--- a/PYTHIA6/AliDecayerPythia.cxx
+++ b/PYTHIA6/AliDecayerPythia.cxx
@@ -75,6 +75,7 @@ AliDecayerPythia::AliDecayerPythia(const AliDecayerPythia &decayer):
 {
     // Copy Constructor
     decayer.Copy(*this);
+    for (Int_t i = 0; i < 501; i++) fBraPart[i] = 0.;
 }
 
 void AliDecayerPythia::Init()
@@ -122,7 +123,7 @@ void AliDecayerPythia::Init()
     fPythia->SetMDCY(fPythia->Pycomp(310) ,1, isw);
     fPythia->SetMDCY(fPythia->Pycomp(3122),1, isw);
     fPythia->SetMDCY(fPythia->Pycomp(3112),1, isw);
-    fPythia->SetMDCY(fPythia->Pycomp(3212),1, isw);
+//    fPythia->SetMDCY(fPythia->Pycomp(3212),1, isw); // Sigma0 decays elem.
     fPythia->SetMDCY(fPythia->Pycomp(3222),1, isw);
     fPythia->SetMDCY(fPythia->Pycomp(3312),1, isw);
     fPythia->SetMDCY(fPythia->Pycomp(3322),1, isw);
@@ -160,7 +161,6 @@ void AliDecayerPythia::ForceDecay()
 //
     Decay_t decay=fDecay;
     fPythia->SetMSTJ(21,2);
-    if (decay == kNoDecayHeavy) return;
 
 //
 // select mode    
@@ -250,6 +250,9 @@ void AliDecayerPythia::ForceDecay()
        ForceParticleDecay(100553,13,2);// Upsilon'
        ForceParticleDecay(200553,13,2);// Upsilon''
        break;
+    case kJpsiDiMuon:
+       ForceParticleDecay(  443,13,2); // J/Psi
+       break;
     case kBSemiElectronic:
        ForceParticleDecay(  511,11,1); // B0     
        ForceParticleDecay(  521,11,1); // B+/-     
@@ -286,6 +289,14 @@ void AliDecayerPythia::ForceDecay()
        ForceParticleDecay(100553,11,2);// Upsilon'
        ForceParticleDecay(200553,11,2);// Upsilon''
        break;
+    case kPsiPrimeJpsiDiElectron:
+       products[0] =    443;
+        products[1] =    211;
+        mult[0] = 1;
+        mult[1] = 2;
+        ForceParticleDecay( 100443, products, mult, 2, 1);
+        ForceParticleDecay(    443,11,2);             
+        break; 
     case kBJpsiDiMuon:
 
        products[0] =    443;
@@ -370,8 +381,10 @@ void AliDecayerPythia::ForceDecay()
        break;
     case kOmega:
        ForceOmega();
+       break;
     case kLambda:
        ForceLambda();
+       break;
     case kAll:
        break;
     case kNoDecay:
@@ -380,6 +393,33 @@ void AliDecayerPythia::ForceDecay()
     case kNoDecayHeavy:
     case kNeutralPion:
        break;
+    case kNoDecayBeauty:                                                    
+        SwitchOffBDecay(); 
+       break;
+    case kElectronEM:
+        ForceParticleDecay(  111,11,1); // pi^0     
+       ForceParticleDecay(  221,11,1); // eta     
+       ForceParticleDecay(  113,11,1); // rho     
+       ForceParticleDecay(  223,11,1); // omega     
+       ForceParticleDecay(  331,11,1); // etaprime     
+       ForceParticleDecay(  333,11,1); // phi     
+       break;
+    case kDiElectronEM:
+       ForceParticleDecay(  111,11,2); // pi^0     
+       ForceParticleDecay(  221,11,2); // eta     
+       ForceParticleDecay(  113,11,2); // rho     
+       ForceParticleDecay(  223,11,2); // omega     
+       ForceParticleDecay(  331,11,2); // etaprime     
+       ForceParticleDecay(  333,11,2); // phi     
+       break;
+    case kGammaEM:
+       ForceParticleDecay(  111,22,1); // pi^0     
+       ForceParticleDecay(  221,22,1); // eta     
+       ForceParticleDecay(  113,22,1); // rho     
+       ForceParticleDecay(  223,22,1); // omega     
+       ForceParticleDecay(  331,22,1); // etaprime     
+       ForceParticleDecay(  333,22,1); // phi     
+       break;
     }
 }
 
@@ -610,7 +650,7 @@ void AliDecayerPythia::ForceParticleDecay(Int_t particle, Int_t product, Int_t m
     if (norm > 0.) fBraPart[kc] /= norm;
 }
 
-void AliDecayerPythia::ForceParticleDecay(Int_t particle, Int_t* products, Int_t* mult, Int_t npart)
+void AliDecayerPythia::ForceParticleDecay(Int_t particle, Int_t* products, Int_t* mult, Int_t npart, Bool_t flag)
 {
 //
 //  Force decay of particle into products with multiplicity mult
@@ -634,9 +674,9 @@ void AliDecayerPythia::ForceParticleDecay(Int_t particle, Int_t* products, Int_t
        for (Int_t i = 0; i < npart; i++) {
            nprod += (CountProducts(channel, products[i]) >= mult[i]);
        }
-       if (nprod) {
+       if ((nprod && !flag) || ((nprod == npart) && flag)) {
            fPythia->SetMDME(channel,1,1);
-       } else {
+       } else {                // 
            fPythia->SetMDME(channel,1,0);
            fBraPart[kc] -= fPythia->GetBRAT(channel);
        }
@@ -747,6 +787,14 @@ void  AliDecayerPythia::ForceLambda()
 }
 
 
+void AliDecayerPythia::SwitchOffBDecay(){
+  Int_t heavyB[]={511,521,531,5122,5132,5232,5332};
+  for(int i=0;i<4;i++)
+    {
+      fPythia->SetMDCY(fPythia->Pycomp(heavyB[i]),1,0);
+    }
+}
+
 Float_t  AliDecayerPythia::GetPartialBranchingRatio(Int_t kf)
 {
 // Get branching ratio