Correction needed for decayer.
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 24 May 2012 10:59:22 +0000 (10:59 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 24 May 2012 10:59:22 +0000 (10:59 +0000)
PYTHIA8/AliTPythia8.cxx

index 776bec0..01b48e0 100644 (file)
@@ -164,19 +164,24 @@ Int_t AliTPythia8::ImportParticles(TClonesArray *particles, Option_t *option)
    clonesParticles.Clear();
    Int_t nparts=0;
    Int_t i;
-   fNumberOfParticles  = fPythia->event.size() - 1;
+   Int_t ioff = 0;
+   fNumberOfParticles  = fPythia->event.size();
+   if (fPythia->event[0].id() == 90) {
+     fNumberOfParticles--; 
+     ioff = -1;
+   }
 
    if (!strcmp(option,"") || !strcmp(option,"Final")) {
-      for (i = 0; i <= fNumberOfParticles; i++) {
+      for (i = 0; i < fNumberOfParticles; i++) {
        if (fPythia->event[i].id() == 90) continue;
          if (fPythia->event[i].isFinal()) {
             new(clonesParticles[nparts]) TParticle(
                 fPythia->event[i].id(),
                 fPythia->event[i].isFinal(),
-                fPythia->event[i].mother1() - 1,
-                fPythia->event[i].mother2() - 1,
-                fPythia->event[i].daughter1() - 1, 
-                fPythia->event[i].daughter2() - 1,
+                fPythia->event[i].mother1() + ioff,
+                fPythia->event[i].mother2() + ioff,
+                fPythia->event[i].daughter1() + ioff, 
+                fPythia->event[i].daughter2() + ioff,
                 fPythia->event[i].px(),     // [GeV/c]
                 fPythia->event[i].py(),     // [GeV/c]
                 fPythia->event[i].pz(),     // [GeV/c]
@@ -189,15 +194,15 @@ Int_t AliTPythia8::ImportParticles(TClonesArray *particles, Option_t *option)
            } // final state partice
        } // particle loop
     } else if (!strcmp(option,"All")) {
-       for (i = 0; i <= fNumberOfParticles; i++) {
+       for (i = 0; i < fNumberOfParticles; i++) {
          if (fPythia->event[i].id() == 90) continue;
            new(clonesParticles[nparts]) TParticle(
                fPythia->event[i].id(),
                fPythia->event[i].isFinal(),
-               fPythia->event[i].mother1() - 1,
-               fPythia->event[i].mother2() - 1,
-               fPythia->event[i].daughter1() - 1,
-               fPythia->event[i].daughter2() - 1,
+               fPythia->event[i].mother1() + ioff,
+               fPythia->event[i].mother2() + ioff,
+               fPythia->event[i].daughter1() + ioff,
+               fPythia->event[i].daughter2() + ioff,
                fPythia->event[i].px(),       // [GeV/c]
                fPythia->event[i].py(),       // [GeV/c]
                fPythia->event[i].pz(),       // [GeV/c]
@@ -217,16 +222,23 @@ TObjArray* AliTPythia8::ImportParticles(Option_t* /* option */)
 {
    // Import particles from Pythia stack
    fParticles->Clear();
-   Int_t numpart   = fPythia->event.size() - 1;
+   Int_t ioff = 0;
+   Int_t numpart   = fPythia->event.size();
+   if (fPythia->event[0].id() == 90) {
+     numpart--; 
+     ioff = -1;
+   }
+
+
    TClonesArray &a = *((TClonesArray*)fParticles);
    for (Int_t i = 1; i <= numpart; i++) {
       new(a[i]) TParticle(
          fPythia->event[i].id(),
          fPythia->event[i].isFinal(),
-         fPythia->event[i].mother1()  - 1,
-         fPythia->event[i].mother2()  - 1,
-         fPythia->event[i].daughter1() - 1,
-         fPythia->event[i].daughter2() - 1,
+         fPythia->event[i].mother1()   + ioff,
+         fPythia->event[i].mother2()   + ioff,
+         fPythia->event[i].daughter1() + ioff,
+         fPythia->event[i].daughter2() + ioff,
          fPythia->event[i].px(),       // [GeV/c]
          fPythia->event[i].py(),       // [GeV/c]
          fPythia->event[i].pz(),       // [GeV/c]