NcInvmass.cxx

Go to the documentation of this file.

 
#include "NcInvmass.h"
#include "Riostream.h"
 
ClassImp(NcInvmass); // Class implementation to enable ROOT I/O
 
NcInvmass::NcInvmass()
{
 
 fPi=acos(-1.);
 fMode=2;
 fBkg=0;
 fNewtheta=1;
 fNewphi=1;
 fMinv=0;
 fMbkg=0;
}

NcInvmass::~NcInvmass()
{
 
 if (fMinv)
 {
  delete fMinv;
  fMinv=0;
 }
 
 if (fMbkg)
 {
  delete fMbkg;
  fMbkg=0;
 }
}

void NcInvmass::SetStorageMode(Int_t m)
{
 
 fMode=2;
 if (m==1) fMode=1;
}

void NcInvmass::SetThetaSwitch(Int_t i)
{
 
 fNewtheta=1;
 if (i==0) fNewtheta=0;
}

void NcInvmass::SetPhiSwitch(Int_t i)
{
 
 fNewphi=1;
 if (i==0) fNewphi=0;
}

Int_t NcInvmass::GetStorageMode() const
{
 
 return fMode;
}

Int_t NcInvmass::GetThetaSwitch() const
{
 
 return fNewtheta;
}

Int_t NcInvmass::GetPhiSwitch() const
{
 
 return fNewphi;
}

TObjArray* NcInvmass::Invmass(TObjArray* a1,TObjArray* a2)
{
 
 fBkg=0;
 Combine(a1,a2);
 return fMinv;
}

TObjArray* NcInvmass::CombBkg(TObjArray* a1,TObjArray* a2)
{
 
 fBkg=1;
 Combine(a1,a2);
 if (fMode != 1)
 {
  return fMbkg;
 }
 else
 {
  return fMinv;
 }
}

void NcInvmass::Combine(TObjArray* a1,TObjArray* a2)
{
 
 if ((!fBkg || fMode==1) && fMinv)
 {
  delete fMinv;
  fMinv=0;
 }
 
 if (fBkg && (fMode !=1) && fMbkg)
 {
  delete fMbkg;
  fMbkg=0;
 }
 
 Int_t isame; // Indicates whether both lists are identical
 isame=0;
 if (a1==a2) isame=1;
 
 // Index i must loop over the shortest of a1 and a2
 TObjArray* listi=a1;
 TObjArray* listj=a2;
 Int_t ni=a1->GetEntries();
 Int_t nj=a2->GetEntries();
 if (nj < ni)
 {
  ni=a2->GetEntries();
  nj=a1->GetEntries();
  listi=a2;
  listj=a1;
 }
 
 NcTrack* p1=0;  
 NcTrack* p2=0;
 NcTrack* px=0;
 Nc4Vector ptot;
 NcTrack* t=0;
 Double_t v2[4],vx[4];
 Float_t q1,q2;
 
 Int_t jmin; // Start index for list j
 Int_t jx;   // Index for randomly picked particle for comb. bkg. reconstruction 
 
 for (Int_t i=0; i<ni; i++) // Select first a particle from list i
 {
  p1=(NcTrack*)listi->At(i);
  p2=0;
 
  if (!p1) continue;
 
  jmin=0;
  if (isame) jmin=i+1;
  for (Int_t j=jmin; j<nj; j++) // Select also a particle from list j
  {
   p2=(NcTrack*)listj->At(j);
   if (p1==p2) p2=0; // Don't combine particle with itself
 
   if (!p2) continue;
 
   p2->GetVector(v2,"sph");
 
   // Take theta and phi from randomly chosen other list j particle for bkg. reconstr.
   if (fBkg)
   {
    px=0; 
    if ((!isame && nj>1) || (isame && nj>2))
    {
     jx=int(fRndm.Uniform(0,float(nj)));
     px=(NcTrack*)listj->At(jx);
     
     while (!px || px==p2 || px==p1)
     {
      jx++;
      if (jx >= nj) jx=0;
      px=(NcTrack*)listj->At(jx);
     }
 
     px->GetVector(vx,"sph");
     if (fNewtheta) v2[2]=vx[2]; // Replace the theta angle in the v2 vector
     if (fNewphi) v2[3]=vx[3]; // Replace the phi angle in the v2 vector
    }
   }
 
   if ((!fBkg && p2) || (fBkg && px))
   {
    // Store the data of this two-particle combination
    ptot.SetVector(v2,"sph");
    ptot=(Nc4Vector)(ptot+(*p1));
    q1=p1->GetCharge();
    q2=p2->GetCharge();
    t=new NcTrack;
    t->Set4Momentum(ptot);
    t->SetCharge(q1+q2);
    if (!fBkg || fMode==1) 
    {
     if (!fMinv)
     {
      fMinv=new TObjArray();
      fMinv->SetOwner();
     }
     fMinv->Add(t);
    }
    else
    {
     if (!fMbkg)
     {
      fMbkg=new TObjArray();
      fMbkg->SetOwner();
     }
     fMbkg->Add(t);
    }
   }
 
  } // End of second particle loop
 
 } // End of first particle loop
}