Bug correction (energy recalculation when adding new incarnation of the particle).

diff --git a/HBTAN/AliHBTParticle.h b/HBTAN/AliHBTParticle.h
index 3e7d1cc22f699afa7255cf162e9f4abbe7380942..429f37eb0a46da9aa7bcae1397a668aec50a75cd 100644 (file)
--- a/HBTAN/AliHBTParticle.h
+++ b/HBTAN/AliHBTParticle.h
@@ -97,6 +97,7 @@ public:
                              {fVx=vx; fVy=vy; fVz=vz; fVt=t;}
   void           SetProductionVertex(const TLorentzVector& v)
                              {SetProductionVertex(v.X(),v.Y(),v.Z(),v.T());}
                              {fVx=vx; fVy=vy; fVz=vz; fVt=t;}
   void           SetProductionVertex(const TLorentzVector& v)
                              {SetProductionVertex(v.X(),v.Y(),v.Z(),v.T());}

+  void           SetCalcMass(Double_t mass) {fCalcMass = mass;}

   const Char_t*  GetName() const; 
   void           Print() const;
 
   const Char_t*  GetName() const; 
   void           Print() const;
 

diff --git a/HBTAN/AliHBTReaderInternal.cxx b/HBTAN/AliHBTReaderInternal.cxx
index 333828278e60e1f599d4d802632ab16fd3fc506d..5b25c1a2d45ff4370024abebab0c7eca087764e4 100644 (file)
--- a/HBTAN/AliHBTReaderInternal.cxx
+++ b/HBTAN/AliHBTReaderInternal.cxx
@@ -126,13 +126,14 @@ Int_t AliHBTReaderInternal::ReadNext()
 
     counter = 0;  
     if (fPartBranch && fTrackBranch)
 
     counter = 0;  
     if (fPartBranch && fTrackBranch)

-    {
+     {

        for(i = 0; i < fPartBuffer->GetEntries(); i++)
          {
            tpart = dynamic_cast<AliHBTParticle*>(fPartBuffer->At(i));
            ttrack =  dynamic_cast<AliHBTParticle*>(fTrackBuffer->At(i));
 
            if( tpart == 0x0 ) continue; //if returned pointer is NULL
        for(i = 0; i < fPartBuffer->GetEntries(); i++)
          {
            tpart = dynamic_cast<AliHBTParticle*>(fPartBuffer->At(i));
            ttrack =  dynamic_cast<AliHBTParticle*>(fTrackBuffer->At(i));
 
            if( tpart == 0x0 ) continue; //if returned pointer is NULL

+           if( ttrack == 0x0 ) continue; //if returned pointer is NULL

 
            if (ttrack->GetUID() != tpart->GetUID())
              {
 
            if (ttrack->GetUID() != tpart->GetUID())
              {


@@ -140,19 +141,39 @@ Int_t AliHBTReaderInternal::ReadNext()
                Error("ReadNext","They probobly do not correspond to each other.");
              }
 
                Error("ReadNext","They probobly do not correspond to each other.");
              }
 

-           for (Int_t s = 0; s < tpart->GetNumberOfPids(); s++)
+           for (Int_t s = 0; s < ttrack->GetNumberOfPids(); s++)

             {
             {

-              if( pdgdb->GetParticle(tpart->GetNthPid(s)) == 0x0 ) continue; //if particle has crazy PDG code (not known to our database)
-              if( Pass(tpart->GetNthPid(s)) ) continue; //check if we are intersted with particles of this type
+              if( pdgdb->GetParticle(ttrack->GetNthPid(s)) == 0x0 ) continue; //if particle has crazy PDG code (not known to our database)
+              if( Pass(ttrack->GetNthPid(s)) ) continue; //check if we are intersted with particles of this type

                                           //if not take next partilce
                                           //if not take next partilce

-              AliHBTParticle* part = new AliHBTParticle(*tpart);
-              part->SetPdgCode(tpart->GetNthPid(s),tpart->GetNthPidProb(s));
-              if( Pass(part) )
+              
+              AliHBTParticle* track = new AliHBTParticle(*ttrack);
+              
+              //apart of setting PDG code of an incarnation
+              //it is necessary tu recalculate energy on the basis of
+              //new PDG code (mass) hypothesis
+              TParticlePDG* pdgp = pdgdb->GetParticle(ttrack->GetNthPid(s));
+              if (pdgp == 0x0)//PDG part corresponding to new incarnation
+               {
+                 Error("ReadNext","Particle code unknown to PDG DB.");
+                 continue;
+               }
+              Double_t mass = pdgp->Mass();//mass of new incarnation
+              Double_t tEtot = TMath::Sqrt( ttrack->Px()*ttrack->Px() + 
+                                            ttrack->Py()*ttrack->Py() + 
+                                            ttrack->Pz()*ttrack->Pz() + 
+                                            mass*mass);//total energy of the new incarnation
+              //update Energy and Calc Mass 
+              track->SetMomentum(ttrack->Px(),ttrack->Py(),ttrack->Pz(),tEtot);
+              track->SetCalcMass(mass);
+              track->SetPdgCode(ttrack->GetNthPid(s),ttrack->GetNthPidProb(s));
+              
+              if( Pass(track) )

                 {
                 {

-                  delete part;
+                  delete track;

                   continue; 
                 }
                   continue; 
                 }

-              AliHBTParticle* track = new AliHBTParticle(*ttrack);
+              AliHBTParticle* part = new AliHBTParticle(*tpart);//particle has only 1 incarnation (real)

 
               fParticlesEvent->AddParticle(part);
               fTracksEvent->AddParticle(track);
 
               fParticlesEvent->AddParticle(part);
               fTracksEvent->AddParticle(track);


@@ -176,7 +197,6 @@ Int_t AliHBTReaderInternal::ReadNext()
             {
               if( pdgdb->GetParticle(tpart->GetNthPid(s)) == 0x0 ) continue; //if particle has crazy PDG code (not known to our database)
               if( Pass(tpart->GetNthPid(s)) ) continue; //check if we are intersted with particles of this type
             {
               if( pdgdb->GetParticle(tpart->GetNthPid(s)) == 0x0 ) continue; //if particle has crazy PDG code (not known to our database)
               if( Pass(tpart->GetNthPid(s)) ) continue; //check if we are intersted with particles of this type

-

               AliHBTParticle* part = new AliHBTParticle(*tpart);
               part->SetPdgCode(tpart->GetNthPid(s),tpart->GetNthPidProb(s));
               if( Pass(part) )
               AliHBTParticle* part = new AliHBTParticle(*tpart);
               part->SetPdgCode(tpart->GetNthPid(s),tpart->GetNthPidProb(s));
               if( Pass(part) )


@@ -196,21 +216,42 @@ Int_t AliHBTReaderInternal::ReadNext()
       {
         for(i = 0; i < fTrackBuffer->GetEntries(); i++)
          {
       {
         for(i = 0; i < fTrackBuffer->GetEntries(); i++)
          {

-           tpart = dynamic_cast<AliHBTParticle*>(fTrackBuffer->At(i));
-           if(tpart == 0x0) continue; //if returned pointer is NULL
-           for (Int_t s = 0; s < tpart->GetNumberOfPids(); s++)
-            {
-              if( pdgdb->GetParticle(tpart->GetNthPid(s)) == 0x0 ) continue; //if particle has crazy PDG code (not known to our database)
-              if( Pass(tpart->GetNthPid(s)) ) continue; //check if we are intersted with particles of this type
+           ttrack =  dynamic_cast<AliHBTParticle*>(fTrackBuffer->At(i));
+           if( ttrack == 0x0 ) continue; //if returned pointer is NULL

 
 

-              AliHBTParticle* part = new AliHBTParticle(*tpart);
-              part->SetPdgCode(tpart->GetNthPid(s),tpart->GetNthPidProb(s));
-              if( Pass(part) )
+           for (Int_t s = 0; s < ttrack->GetNumberOfPids(); s++)
+            {
+              if( pdgdb->GetParticle(ttrack->GetNthPid(s)) == 0x0 ) continue; //if particle has crazy PDG code (not known to our database)
+              if( Pass(ttrack->GetNthPid(s)) ) continue; //check if we are intersted with particles of this type
+                                          //if not take next partilce
+              AliHBTParticle* track = new AliHBTParticle(*ttrack);
+              
+              //apart of setting PDG code of an incarnation
+              //it is necessary tu recalculate energy on the basis of
+              //new PDG code (mass) hypothesis
+              TParticlePDG* pdgp = pdgdb->GetParticle(ttrack->GetNthPid(s));
+              if (pdgp == 0x0)//PDG part corresponding to new incarnation
+               {
+                 Error("ReadNext","Particle code unknown to PDG DB.");
+                 continue;
+               }
+              Double_t mass = pdgp->Mass();//mass of new incarnation
+              Double_t tEtot = TMath::Sqrt( ttrack->Px()*ttrack->Px() + 
+                                            ttrack->Py()*ttrack->Py() + 
+                                            ttrack->Pz()*ttrack->Pz() + 
+                                            mass*mass);//total energy of the new incarnation
+              //update Energy and Calc Mass 
+              track->SetMomentum(ttrack->Px(),ttrack->Py(),ttrack->Pz(),tEtot);
+              track->SetCalcMass(mass);
+              track->SetPdgCode(ttrack->GetNthPid(s),ttrack->GetNthPidProb(s));
+              
+              if( Pass(track) )

                 {
                 {

-                  delete part;
+                  delete track;

                   continue; 
                 }
                   continue; 
                 }

-              fTracksEvent->AddParticle(part);
+              fTracksEvent->AddParticle(track);
+

               counter++;
             }
          }
               counter++;
             }
          }