3D added in Lednicky alg.(L.Malinina)
authorskowron <skowron@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 21 Mar 2003 13:34:25 +0000 (13:34 +0000)
committerskowron <skowron@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 21 Mar 2003 13:34:25 +0000 (13:34 +0000)
HBTAN/AliHBTLLWeightFctn.cxx
HBTAN/AliHBTLLWeightFctn.h
HBTAN/AliHBTLLWeights.cxx
HBTAN/AliHBTLLWeights.h

index d2dc5a4..a17d366 100644 (file)
@@ -1,16 +1,7 @@
-/* $Id$ */
-
-//This class allows to obtain Q_inv correlation function with weights
-//calculated by Lednicky's alghorithm.
-//Numerator is filled with weighted events. Weights are attributed to reconstructed tracks.
-//Weights are calculated with corresponding simulated particles momenta.
-//Denominator is filled with mixing unweighted reconstructed tracks.
-//One needs both pairs 
-//(simulated and recontructed), thus function is of class AliHBTTwoPairFctn1D.
 //Author: Ludmila Malinina, JINR (malinina@sunhe.jinr.ru)
-
 #include "AliHBTLLWeightFctn.h"
 #include "AliHBTLLWeights.h"
+#include "AliHBTLLWeightsPID.h"
 
 //--for test--AliHBTLLWeightQInvFctn* yyy= new AliHBTLLWeightQInvFctn();
 
@@ -32,14 +23,13 @@ void  AliHBTLLWeightQInvFctn::ProcessSameEventParticles(AliHBTPair* trackpair, A
     Double_t weightPID=1.;
     Double_t weightHBT=AliHBTLLWeights::Instance()->GetWeight(partpair);
     Double_t weight=weightHBT*weightPID;
-    if(TMath::Abs(weight)<=10.) fNumerator->Fill(trackpair->GetQInv(),weight);
+    if(TMath::Abs(weight)<=10.)fNumerator->Fill(trackpair->GetQInv(),weight);
   }
 } 
 /****************************************************************/
 
 void  AliHBTLLWeightQInvFctn::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
 {
-  // Fills the denominator using mixed pairs
   trackpair = CheckPair(trackpair);
   partpair  = CheckPair(partpair);
   if ( trackpair && partpair)  
@@ -55,3 +45,156 @@ TH1* AliHBTLLWeightQInvFctn::GetResult()
  return GetRatio(Scale());                                                  
 }                    
                                                               
+/**************************************************************************************/
+/**************************************************************************************/
+/**************************************************************************************/
+/**************************************************************************************/
+
+ClassImp(AliHBTLLWeightQOutFctn)
+/****************************************************************/
+void AliHBTLLWeightQOutFctn::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
+{
+  trackpair = CheckPair(trackpair);
+  partpair  = CheckPair(partpair);
+  if ( trackpair && partpair)     
+  {
+    Double_t weightPID=1.;
+    Double_t weightHBT=AliHBTLLWeights::Instance()->GetWeight(partpair);
+    Double_t weight=weightHBT*weightPID;
+    if(TMath::Abs(weight)<=10.) fNumerator->Fill(trackpair->GetQOutCMSLC(),weight);
+  }
+} 
+/****************************************************************/
+
+void AliHBTLLWeightQOutFctn::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
+{
+  trackpair = CheckPair(trackpair);
+  partpair  = CheckPair(partpair);
+  if ( trackpair && partpair)  
+  {
+     fDenominator->Fill(trackpair->GetQOutCMSLC());
+  }
+}
+/**************************************************************/
+TH1* AliHBTLLWeightQOutFctn::GetResult() 
+                                                                               
+{ 
+//returns ratio of numerator and denominator                                    
+ return GetRatio(Scale());                                                  
+}                    
+                                                              
+/*************************************************************************************/ 
+/*************************************************************************************/ 
+/*************************************************************************************/ 
+
+ClassImp(AliHBTLLWeightQLongFctn)
+/****************************************************************/
+void  AliHBTLLWeightQLongFctn::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
+{
+//Processes Particles and tracks Same different even
+  trackpair = CheckPair(trackpair);
+  partpair  = CheckPair(partpair);
+  if ( trackpair && partpair)     
+  {
+    Double_t weightPID=1.;
+    Double_t weightHBT=AliHBTLLWeights::Instance()->GetWeight(partpair);
+    Double_t weight=weightHBT*weightPID;
+    if(TMath::Abs(weight)<=10.) fNumerator->Fill(trackpair->GetQLongCMSLC(),weight);
+  }
+} 
+/****************************************************************/
+
+void  AliHBTLLWeightQLongFctn::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
+{
+  trackpair = CheckPair(trackpair);
+  partpair  = CheckPair(partpair);
+  if ( trackpair && partpair)  
+  {
+     fDenominator->Fill(trackpair->GetQLongCMSLC());
+  }
+}
+/**************************************************************/
+TH1* AliHBTLLWeightQLongFctn::GetResult() 
+                                                                               
+{ 
+//returns ratio of numerator and denominator                                    
+ return GetRatio(Scale());                                                  
+}                    
+                                                              
+/*************************************************************************************/ 
+/*************************************************************************************/ 
+/*************************************************************************************/ 
+
+ClassImp(AliHBTLLWeightQSideFctn)
+/****************************************************************/
+void  AliHBTLLWeightQSideFctn::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
+{
+//Processes Particles and tracks Same different even
+  trackpair = CheckPair(trackpair);
+  partpair  = CheckPair(partpair);
+  if ( trackpair && partpair)     
+  {
+    Double_t weightPID=1.;
+    Double_t weightHBT=AliHBTLLWeights::Instance()->GetWeight(partpair);
+    Double_t weight=weightHBT*weightPID;
+    if(TMath::Abs(weight)<=10.) fNumerator->Fill(trackpair->GetQSideCMSLC(),weight);
+  }
+} 
+/****************************************************************/
+
+void  AliHBTLLWeightQSideFctn::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
+{
+  trackpair = CheckPair(trackpair);
+  partpair  = CheckPair(partpair);
+  if ( trackpair && partpair)  
+  {
+     fDenominator->Fill(trackpair->GetQSideCMSLC());
+  }
+}
+/**************************************************************/
+TH1* AliHBTLLWeightQSideFctn::GetResult() 
+                                                                               
+{ 
+//returns ratio of numerator and denominator                                    
+ return GetRatio(Scale());                                                  
+}                    
+                                                              
+/*************************************************************************************/ 
+/*************************************************************************************/ 
+
+ClassImp(AliHBTLLWeightTwoKStarFctn)
+/****************************************************************/
+void  AliHBTLLWeightTwoKStarFctn::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
+{
+//Processes Particles and tracks Same different even
+  trackpair = CheckPair(trackpair);
+  partpair  = CheckPair(partpair);
+  if ( trackpair && partpair)     
+  {
+    Double_t weightPID=1.;
+    Double_t weightHBT=AliHBTLLWeights::Instance()->GetWeight(partpair);
+    Double_t weight=weightHBT*weightPID;
+    if(TMath::Abs(weight)<=10.) fNumerator->Fill(2.0*(trackpair->GetKStar()),weight);
+  }
+} 
+/****************************************************************/
+
+void  AliHBTLLWeightTwoKStarFctn::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
+{
+  trackpair = CheckPair(trackpair);
+  partpair  = CheckPair(partpair);
+  if ( trackpair && partpair)  
+  {
+     fDenominator->Fill(2.0*(trackpair->GetKStar()));
+  }
+}
+/**************************************************************/
+TH1* AliHBTLLWeightTwoKStarFctn::GetResult() 
+                                                                               
+{ 
+//returns ratio of numerator and denominator                                    
+ return GetRatio(Scale());                                                  
+}                    
+                                                              
+/*************************************************************************************/ 
+
index 3200e91..885104b 100644 (file)
@@ -1,9 +1,13 @@
-#ifndef ALIHBTLLWEIGHTQINVFCTN_H
-#define ALIHBTLLWEIGHTQINVFCTN_H
-
-/* $Id$ */
-
+//This function allows to obtain Q_inv correlation function with weights
+//calculated by Lednicky's alghorithm.
+//Numerator is filled with weighted events. Weights are attributed to reconstructed tracks.
+//Weights are calculated with corresponding simulated particles momenta.
+//Denominator is filled with mixing unweighted reconstructed tracks.
+//One needs both pairs 
+//(simulated and recontructed), thus function is of class AliHBTTwoPairFctn1D.
 
+#ifndef ALIHBTLLWEIGHTFCTN_H
+#define ALIHBTLLWEIGHTFCTN_H
 #include "AliHBTFunction.h"
 
 
@@ -11,7 +15,7 @@ class AliHBTLLWeights;
 
 class AliHBTLLWeightQInvFctn: public AliHBTTwoPairFctn1D
 {
-  friend class AliHBTOnePairFctn1D;
+//  friend class AliHBTOnePairFctn1D;
 
   public:
       AliHBTLLWeightQInvFctn(Int_t nbins = 100, Double_t maxXval = 0.15, Double_t minXval = 0.0);
@@ -21,11 +25,98 @@ class AliHBTLLWeightQInvFctn: public AliHBTTwoPairFctn1D
       void   ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
       void   ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
       
-      Double_t GetValue(AliHBTPair* trackpair, AliHBTPair* partpair) const
-       { return trackpair->GetQInv()-partpair->GetQInv();} //isn't use                                                                    
-  private:
+      Double_t GetValue(AliHBTPair* trackpair, AliHBTPair* partpair)         
+           { return trackpair->GetQInv()-partpair->GetQInv();} //isn't use                                                                    
+       
+
+  protected:
 
+  private:
+  public:
      ClassDef(AliHBTLLWeightQInvFctn,1)
 };
+/*************************************************************************************/ 
+/*************************************************************************************/ 
+/*************************************************************************************/ 
+/*************************************************************************************/ 
+
+class AliHBTLLWeightQOutFctn: public AliHBTTwoPairFctn1D
+{
+
+ friend class AliHBTOnePairFctn1D;
+ public:
+   AliHBTLLWeightQOutFctn(Int_t nbins = 100, Double_t maxXval = 0.15, Double_t minXval = 0.0):
+                        AliHBTTwoPairFctn1D(nbins,maxXval,minXval){}
+   virtual ~AliHBTLLWeightQOutFctn(){};
+   TH1* GetResult();
+ protected:
+   void   ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
+   void   ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
+      
+    Double_t GetValue(AliHBTPair* trackpair, AliHBTPair* partpair)         
+           { return trackpair->GetQOutCMSLC()-partpair->GetQOutCMSLC();} //isn't use                                                                    
+  public:
+    ClassDef(AliHBTLLWeightQOutFctn,1)
+};
+/*************************************************************************************/ 
+class AliHBTLLWeightQLongFctn: public AliHBTTwoPairFctn1D
+{
+
+ friend class AliHBTOnePairFctn1D;
+ public:
+   AliHBTLLWeightQLongFctn(Int_t nbins = 100, Double_t maxXval = 0.15, Double_t minXval = 0.0):
+                        AliHBTTwoPairFctn1D(nbins,maxXval,minXval){}
+   virtual ~AliHBTLLWeightQLongFctn(){};
+   TH1* GetResult();
+ protected:
+   void   ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
+   void   ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
+      
+    Double_t GetValue(AliHBTPair* trackpair, AliHBTPair* partpair)         
+           { return trackpair->GetQLongCMSLC()-partpair->GetQLongCMSLC();} //isn't use                                                                    
+  public:
+    ClassDef(AliHBTLLWeightQLongFctn,1)
+};
+/*************************************************************************************/ 
+class AliHBTLLWeightQSideFctn: public AliHBTTwoPairFctn1D
+{
+
+ friend class AliHBTOnePairFctn1D;
+ public:
+   AliHBTLLWeightQSideFctn(Int_t nbins = 100, Double_t maxXval = 0.15, Double_t minXval = 0.0):
+                        AliHBTTwoPairFctn1D(nbins,maxXval,minXval){}
+   virtual ~AliHBTLLWeightQSideFctn(){};
+   TH1* GetResult();
+ protected:
+   void   ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
+   void   ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
+      
+    Double_t GetValue(AliHBTPair* trackpair, AliHBTPair* partpair)         
+           { return trackpair->GetQLongCMSLC()-partpair->GetQLongCMSLC();} //isn't use                                                                    
+  public:
+    ClassDef(AliHBTLLWeightQSideFctn,1) 
+};
+/*************************************************************************************/ 
+class AliHBTLLWeightTwoKStarFctn: public AliHBTTwoPairFctn1D
+{
+
+ friend class AliHBTOnePairFctn1D;
+ public:
+   AliHBTLLWeightTwoKStarFctn(Int_t nbins = 100, Double_t maxXval = 0.15, Double_t minXval = 0.0):
+                        AliHBTTwoPairFctn1D(nbins,maxXval,minXval){}
+   virtual ~AliHBTLLWeightTwoKStarFctn(){};
+   TH1* GetResult();
+ protected:
+   void   ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
+   void   ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair);
+      
+    Double_t GetValue(AliHBTPair* trackpair, AliHBTPair* partpair)         
+           { return trackpair->GetKStar()-partpair->GetKStar();} //isn't use                                                                    
+  public:
+    ClassDef(AliHBTLLWeightTwoKStarFctn,1) 
+};
+
   
 #endif
index 752f5fc..67adf0e 100644 (file)
@@ -1,23 +1,11 @@
-/* $Id$ */
-
-//------------------------------------------------------------------
-// This class introduces the weight's calculation 
-// according to the Lednicky's algorithm.
-// The detailed description of the algorithm can be found 
-// in comments to fortran code:
-// fsiw.f, fsiini.f
-// Author:
-//------------------------------------------------------------------
-
-#include <TMath.h>
-#include <TPDGCode.h>
-#include <TRandom.h>
-
 #include "AliHBTLLWeights.h"
+
 #include "AliHBTPair.h"
 #include "AliHBTParticle.h"
-#include "WLedCOMMONS.h"
-
+#include <TList.h>
+#include <TRandom.h>                                                                     
+#include <TMath.h>
+#include <TPDGCode.h>
 /*******************************************************************/
 /******      ROUTINES    USED    FOR     COMMUNUCATION      ********/
 /********************     WITH      FORTRAN     ********************/
@@ -48,189 +36,190 @@ AliHBTLLWeights* AliHBTLLWeights::fgLLWeights=NULL;
 
 AliHBTLLWeights::AliHBTLLWeights()
 {
-  // Default Constructor 
-  fPID1 = 0;
-  fPID2 = 0;
-  SetRandomPosition();
-  SetColWithResidNuclSwitch();
-  SetStrongInterSwitch();
-  SetQuantumStatistics();
-  SetColoumb();
-  SetTest();
+// Default Constructor 
+    fPID1 = 0;
+    fPID2 = 0;
+    SetRandomPosition();
+    SetColWithResidNuclSwitch();
+    SetStrongInterSwitch();
+    SetQuantumStatistics();
+    SetColoumb();
+    SetTest();
   
 }
 
 
-AliHBTLLWeights* AliHBTLLWeights::Instance()
-{
-  // Instantiates new object or returns a pointer to already exitsing one
-  if (fgLLWeights) {
-    return fgLLWeights;
-  } else {
-    fgLLWeights = new AliHBTLLWeights();
-    return fgLLWeights;
-  }
-}
+ AliHBTLLWeights* AliHBTLLWeights::Instance()
+{                                                                                             
+    if (fgLLWeights) {                                                                        
+      return fgLLWeights;                                                                   
+   } else {                                                                                  
+   fgLLWeights = new AliHBTLLWeights();                                                        
+      return fgLLWeights;                                                                   
+  }                                                                                         
+}                                                                                             
                                             
 
 Double_t AliHBTLLWeights::GetWeight(const AliHBTPair* partpair)
 {
-  // Returns the weignt of the pair "partpair"
-  AliHBTParticle *part1 = partpair->Particle1();
-  AliHBTParticle *part2 = partpair->Particle2();
+    AliHBTParticle *part1 = partpair->Particle1();
+    AliHBTParticle *part2 = partpair->Particle2();
 
-  if ( (part1 == 0x0) || (part2 == 0x0))
-    {
+    if ( (part1 == 0x0) || (part2 == 0x0))
+      {
       Error("GetWeight","Null particle pointer");
       return 0.0;
-    }
-  
+      }
+
    
-  Double_t part1Momentum[]={part1->Px(),part1->Py(),part1->Pz()};
-  Double_t part2Momentum[]={part2->Px(),part2->Py(),part2->Pz()};
-  
-  if ( (part1->Px() == part2->Px()) && 
-       (part1->Py() == part2->Py()) && 
-       (part1->Pz() == part2->Pz()) )
+    Double_t part1Momentum[]={part1->Px(),part1->Py(),part1->Pz()};
+    Double_t part2Momentum[]={part2->Px(),part2->Py(),part2->Pz()};
+              
+    if ( (part1->Px() == part2->Px()) && 
+         (part1->Py() == part2->Py()) && 
+         (part1->Pz() == part2->Pz()) )
     {
-      return 0.0;
+    return 0.0;
     }
-  
+
              
-  if ((!fRandomPosition) && 
-      (part1->Vx()  == part2->Vx()) && (part1->Vy()  == part2->Vy())
-      && (part1->Vz()  == part2->Vz()) )
+    if ((!fRandomPosition) && 
+    (part1->Vx()  == part2->Vx()) && (part1->Vy()  == part2->Vy())
+    && (part1->Vz()  == part2->Vz()) )
     {        
-      return 0.0;
+    return 0.0;
     }
-  
-  
-  
-  FSI_MOM.P1X=part1Momentum[0];
-  FSI_MOM.P1Y=part1Momentum[1];
-  FSI_MOM.P1Z=part1Momentum[2];
-  
-  FSI_MOM.P2X=part2Momentum[0];
-  FSI_MOM.P2Y=part2Momentum[1];
-  FSI_MOM.P2Z=part2Momentum[2];
-  
-  if (fRandomPosition){
-    
-    Double_t rxcm = fsigma*gRandom->Gaus();
-    Double_t rycm = fsigma*gRandom->Gaus();
-    Double_t rzcm = fsigma*gRandom->Gaus();   
-    
-    FSI_PRF.X=rxcm*fwcons;
-    FSI_PRF.Y=rycm*fwcons;
-    FSI_PRF.Z=rzcm*fwcons;
-    FSI_PRF.T=0.;
-    
-    Double_t rps=rxcm*rxcm+rycm*rycm+rzcm*rzcm;
-    Double_t rp=TMath::Sqrt(rps);
-    FSI_PRF.RP=rp;
-    FSI_PRF.RPS=rps;
-    
-  }        
-  
-  ltran12();
-  fsiw();
-  
-  return LEDWEIGHT.WEIN;
-  
-}
 
+
+
+        FSI_MOM.P1X=part1Momentum[0];
+        FSI_MOM.P1Y=part1Momentum[1];
+        FSI_MOM.P1Z=part1Momentum[2];
+      
+        FSI_MOM.P2X=part2Momentum[0];
+        FSI_MOM.P2Y=part2Momentum[1];
+        FSI_MOM.P2Z=part2Momentum[2];
+                
+        if (fRandomPosition){
+
+        Double_t rxcm = fsigma*gRandom->Gaus();
+        Double_t rycm = fsigma*gRandom->Gaus();
+        Double_t rzcm = fsigma*gRandom->Gaus();   
+        
+        FSI_PRF.X=rxcm*fwcons;
+        FSI_PRF.Y=rycm*fwcons;
+        FSI_PRF.Z=rzcm*fwcons;
+        FSI_PRF.T=0.;
+
+        Double_t rps=rxcm*rxcm+rycm*rycm+rzcm*rzcm;
+        Double_t rp=TMath::Sqrt(rps);
+        FSI_PRF.RP=rp;
+        FSI_PRF.RPS=rps;
+
+        }        
+              
+        ltran12();
+        fsiw();
+
+        if(flambda<1){
+       if(gRandom->Rndm()<(1-flambda))LEDWEIGHT.WEIN=1.;}
+       
+        return LEDWEIGHT.WEIN;
+       
+ }
 /************************************************************/
 void AliHBTLLWeights::Init()
-{
-  //---------------------------------------------------------------------  
-  //initial parameters of model
-  
-  FSI_NS.NS = fApproximationModel;      
-  
-  if(!ftest){LEDWEIGHT.ITEST=0;}
-  
-  if(ftest){
-    LEDWEIGHT.ITEST=1;
-    if(fColoumbSwitch){FSI_NS.ICH =1;}
-    else{FSI_NS.ICH=0;}
-    if(fStrongInterSwitch){FSI_NS.ISI=1;}
-    else{FSI_NS.ISI=0;}
-    if(fQuantStatSwitch){FSI_NS.IQS=1;}
-    else{FSI_NS.IQS=0;}
-    if(fColWithResidNuclSwitch){FSI_NS.I3C=1;}
-    else{FSI_NS.I3C=0;}
-  }
-  
-  if(fRandomPosition){LEDWEIGHT.IRANPOS=1;}
-  else{LEDWEIGHT.IRANPOS=0;}
-  
-  
-  if ( (fPID1 == 0) || (fPID2 == 0) )
-    {
-      Fatal("Init","Particles types are not set");
-      return;//pro forma
-    }
-  
-  FSI_NS.LL = GetPairCode(fPID1,fPID2);
-  
-  if (FSI_NS.LL == 0) 
-    {
-      Fatal("Init","Particles types are not supported");
-      return;//pro forma
-    }
-  
-  
-  TParticlePDG* tpart1 = TDatabasePDG::Instance()->GetParticle(fPID1);
-  if (tpart1 == 0x0)
-    {
-      Fatal("init","We can not find particle with ID=%d in our DataBase",fPID1);
-      return;
-    }
-  
-  FSI_POC.AM1=tpart1->Mass();
-  FSI_POC.C1=tpart1->Charge(); 
-  
-  TParticlePDG* tpart2 = TDatabasePDG::Instance()->GetParticle(fPID2);
-  //mlv
-  
-  
-  
-  if (tpart2 == 0x0)
-    {
-      Fatal("init","We can not find particle with ID=%d in our DataBase",fPID2);
-      return;
-    }
-  
-  FSI_POC.AM2=tpart2->Mass();
-  FSI_POC.C1=tpart2->Charge();
-  
-  led_bldata();
-  fsiini();
-  
-  
-  //constants for radii simulation 
-  
-  if(fRandomPosition){
-    fsigma =TMath::Sqrt(2.)*fRadius;     
-    fwcons =FSI_CONS.W;
-  } 
+ {
+//---------------------------------------------------------------------
+
+//initial parameters of model
+
+      FSI_NS.NS = approximationModel;      
+      
+      if(!ftest){LEDWEIGHT.ITEST=0;}
+    
+      if(ftest){
+      LEDWEIGHT.ITEST=1;
+      if(fColoumbSwitch){FSI_NS.ICH =1;}
+      else{FSI_NS.ICH=0;}
+      if(fStrongInterSwitch){FSI_NS.ISI=1;}
+      else{FSI_NS.ISI=0;}
+      if(fQuantStatSwitch){FSI_NS.IQS=1;}
+      else{FSI_NS.IQS=0;}
+      if(fColWithResidNuclSwitch){FSI_NS.I3C=1;}
+      else{FSI_NS.I3C=0;}
+      }
+      
+      if(fRandomPosition){LEDWEIGHT.IRANPOS=1;}
+      else{LEDWEIGHT.IRANPOS=0;}
+
+
+     if ( (fPID1 == 0) || (fPID2 == 0) )
+      {
+        Fatal("Init","Particles types are not set");
+        return;//pro forma
+      }
+
+      FSI_NS.LL = GetPairCode(fPID1,fPID2);
+       
+      if (FSI_NS.LL == 0) 
+       {
+         Fatal("Init","Particles types are not supported");
+         return;//pro forma
+       }
+
+
+      TParticlePDG* tpart1 = TDatabasePDG::Instance()->GetParticle(fPID1);
+      if (tpart1 == 0x0)
+       {
+         Fatal("init","We can not find particle with ID=%d in our DataBase",fPID1);
+         return;
+       }
+      
+      FSI_POC.AM1=tpart1->Mass();
+      FSI_POC.C1=tpart1->Charge(); 
+
+      TParticlePDG* tpart2 = TDatabasePDG::Instance()->GetParticle(fPID2);
+//mlv
+      
+      
+
+      if (tpart2 == 0x0)
+       {
+         Fatal("init","We can not find particle with ID=%d in our DataBase",fPID2);
+         return;
+       }
+
+      FSI_POC.AM2=tpart2->Mass();
+      FSI_POC.C1=tpart2->Charge();
+
+      led_bldata();
+      fsiini();
+
+
+//constants for radii simulation 
+
+      if(fRandomPosition){
+       fsigma =TMath::Sqrt(2.)*fRadius;     
+       fwcons =FSI_CONS.W;
+      } 
 } 
 
 Int_t AliHBTLLWeights::GetPairCode(const AliHBTPair* partpair)
 {
-  // Return the code of the pair "partpair"
-  return GetPairCode(partpair->Particle1()->GetPdgCode(),partpair->Particle2()->GetPdgCode());
+ return GetPairCode(partpair->Particle1()->GetPdgCode(),partpair->Particle2()->GetPdgCode());
 }
 
 Int_t AliHBTLLWeights::GetPairCode(Int_t pid1,Int_t pid2)
 {
-  //   pairCode   1  2  3  4   5    6   7  8  9 10  11  12  13  14 15 16 17 18  19  20   21   22  23 24 25 26    27     28
-  //   hpid:      n  p  n alfa pi+ pi0 pi+ n  p pi+ pi+ pi+ pi- K+ K+ K+ K-  d  d    t   t    K0  K0  d p  p      p      n
-  //   lpid:      n  p  p alfa pi- pi0 pi+ d  d  K-  K+  p   p  K- K+ p  p   d alfa  t  alfa  K0  K0b t t alfa lambda lambda
-  //   NS=1 y/n:  +  +  +  +   +    -   -  -  -  -   -   -   -  -  -  -  -   -  -    -    -    -  -   - -  -      -      -
-  
-  //alphas, deuterons and tyts are NOT supported here
-  
+//   pairCode   1  2  3  4   5    6   7  8  9 10  11  12  13  14 15 16 17 18  19  20   21   22  23 24 25 26    27     28
+//   hpid:      n  p  n alfa pi+ pi0 pi+ n  p pi+ pi+ pi+ pi- K+ K+ K+ K-  d  d    t   t    K0  K0  d p  p      p      n
+//   lpid:      n  p  p alfa pi- pi0 pi+ d  d  K-  K+  p   p  K- K+ p  p   d alfa  t  alfa  K0  K0b t t alfa lambda lambda
+//   NS=1 y/n:  +  +  +  +   +    -   -  -  -  -   -   -   -  -  -  -  -   -  -    -    -    -  -   - -  -      -      -
+
+//alphas, deuterons and tyts are NOT supported here
+
   Int_t chargefactor = 1;
   Int_t hpid; //pid in higher row
   Int_t lpid; //pid in lower row
@@ -238,166 +227,166 @@ Int_t AliHBTLLWeights::GetPairCode(Int_t pid1,Int_t pid2)
   
   Bool_t swap;
   
-  //determine the order of selcetion in switch  
+//determine the order of selcetion in switch  
   if (TMath::Abs(pid1) < TMath::Abs(pid2) ) 
-    {
-      if (pid1<0) chargefactor=-1;
-      hpid=pid2*chargefactor;
-      lpid=pid1*chargefactor;
-      swap = kFALSE;
-    } 
+   {
+    if (pid1<0) chargefactor=-1;
+    hpid=pid2*chargefactor;
+    lpid=pid1*chargefactor;
+    swap = kFALSE;
+   } 
   else 
-    {
-      if (pid2<0) chargefactor=-1;
-      hpid=pid1*chargefactor;
-      lpid=pid2*chargefactor;
-      swap = kTRUE;
-    }
-  
-  //mlv
-  hpid=pid1;
-  lpid=pid2;
-  
-  
-  //Determine the pair code
+   {
+    if (pid2<0) chargefactor=-1;
+    hpid=pid1*chargefactor;
+    lpid=pid2*chargefactor;
+    swap = kTRUE;
+   }
+
+//mlv
+   hpid=pid1;
+   lpid=pid2;
+
+
+//Determine the pair code
   switch (hpid) //switch on first  particle id
-    {
-    case kNeutron:
+   {
+     case kNeutron:
       switch (lpid)
-       {
-       case kNeutron: 
-         code = 1;  //neutron neutron
-         break;
-         
-       case kProton: 
-         code = 3;  //neutron proton
-         break;
-         
-       case kLambda0: 
-         code = 28;  //neutron lambda
-         break;
-         
-       default: 
-         return 0; //given pair not supported
-         break;
-       }
+       {
+         case kNeutron: 
+           code = 1;  //neutron neutron
+           break;
+        
+         case kProton: 
+           code = 3;  //neutron proton
+           break;
+           
+         case kLambda0: 
+           code = 28;  //neutron lambda
+           break;
+           
+         default: 
+           return 0; //given pair not supported
+           break;
+       }
       break;
-      
-    case kProton:
+
+     case kProton:
       switch (lpid)
-       {
-       case kProton:
-         code = 2; //proton proton
-         break;
-         
-       case kLambda0: 
-         code = 27;//proton lambda
-         break;
-         
-       default: 
-         return 0; //given pair not supported
-         break;
-         
-       }
+       {
+         case kProton:
+           code = 2; //proton proton
+           break;
+           
+         case kLambda0: 
+           code = 27;//proton lambda
+           break;
+           
+         default: 
+           return 0; //given pair not supported
+           break;
+           
+       }
       break;
-      
-    case kPiPlus:
-      
+
+     case kPiPlus:
+     
       switch (lpid)
-       {
-       case kPiPlus:
-         code = 7; //piplus piplus
-         break;
-         
-       case kPiMinus:
-         code = 5; //piplus piminus
-         break;
-         
-       case kKMinus:
-         code = 10; //piplus Kminus
-         break;
-         
-       case kKPlus:
-         code = 11; //piplus Kplus
-         break;
-         
-       case kProton:
-         code = 12; //piplus proton
-         chargefactor*=-1;
-         break;
-         
-       default: 
-         return 0; //given pair not supported
-         break;
-       }
+       {
+         case kPiPlus:
+           code = 7; //piplus piplus
+           break;
+
+         case kPiMinus:
+           code = 5; //piplus piminus
+           break;
+        
+         case kKMinus:
+           code = 10; //piplus Kminus
+           break;
+
+         case kKPlus:
+           code = 11; //piplus Kplus
+           break;
+
+         case kProton:
+           code = 12; //piplus proton
+           chargefactor*=-1;
+           break;
+
+         default: 
+           return 0; //given pair not supported
+           break;
+       }
       break;
-    case kPi0:
+     case kPi0:
       switch (lpid)
-       {
-       case kPi0:
-         code = 6;
-         break;
-         
-       default: 
-         return 0; //given pair not supported
-         break;
-       }
+       {
+         case kPi0:
+           code = 6;
+           break;
+           
+         default: 
+           return 0; //given pair not supported
+           break;
+       }
       break;
       
-    case kKPlus:
+     case kKPlus:
       switch (lpid)
-       {
-       case kKMinus:
-         code = 14; //Kplus Kminus
-         break;
-         
-       case kKPlus:
-         code = 15; //Kplus Kplus
-         break;
-         
-       case kProton:
-         code = 16; //Kplus proton
-         break;
-         
-       default: 
-         return 0; //given pair not supported
-         break;
-       }
+       {
+         case kKMinus:
+           code = 14; //Kplus Kminus
+           break;
+
+         case kKPlus:
+           code = 15; //Kplus Kplus
+           break;
+
+         case kProton:
+           code = 16; //Kplus proton
+           break;
+           
+         default: 
+           return 0; //given pair not supported
+           break;
+       }
       break;
       
-    case kKMinus:
+     case kKMinus:
       switch (lpid)
-       {
-       case kProton:
-         code = 17; //Kminus proton
-         chargefactor*=1;
-         break;
-         
-       default: 
-         return 0; //given pair not supported
-         break;
-       }
+       {
+         case kProton:
+           code = 17; //Kminus proton
+           chargefactor*=1;
+           break;
+           
+         default: 
+           return 0; //given pair not supported
+           break;
+       }
       break;
       
-    case kK0:
+     case kK0:
       switch (lpid)
-       {
-       case kK0:
-         code = 2; //Kzero Kzero
-         break;
-         
-       case kK0Bar:
-         code = 17; //Kzero KzeroBar
-         break;
-         
-       default: 
-         return 0; //given pair not supported
-         break;
-       }
+       {
+         case kK0:
+           code = 2; //Kzero Kzero
+           break;
+         
+         case kK0Bar:
+           code = 17; //Kzero KzeroBar
+           break;
+
+         default: 
+           return 0; //given pair not supported
+           break;
+       }
       break;
-      
-    default: return 0;
-    }
+
+     default: return 0;
+   }
   return code;
 }
 
index a533dec..d949502 100644 (file)
-/* $Id$ */
-
-// This class introduces the weight's calculation 
-// according to the Lednicky's algorithm.
-// The detailed description of the algorithm can be found 
-// in comments to fortran code:
-// fsiw.f, fsiini.f  
+//This class introduce the weights calculation according with Lednicky's algorithm.
+//The detailed description of the algorithm can be found in comments to fortran code:
+//fsiw.f, fsiini.f  
 #ifndef ALIHBTLLWEIGHTS_H
 #define ALIHBTLLWEIGHTS_H
 
 #include <TObject.h>
+#include "WLedCOMMONS.h"
 
 class AliHBTPair;
 class AliHBTLLWeights: public TObject
-{
- public:
-  virtual ~AliHBTLLWeights(){;}
-  AliHBTLLWeights(const AliHBTLLWeights &source) {
-    //Copy ctor needed by the coding conventions but not used
-    Fatal("AliHBTLLWeights","copy ctor not implemented");
-  }
-  AliHBTLLWeights & operator=(const AliHBTLLWeights &source) {
-    //Assignment operator needed by the coding conventions but not used
-    Fatal("AliHBTLLWeights","assignment operator not implemented");
-    return * this;
-  }
-  static AliHBTLLWeights* Instance();
-  
-  void Init(); //put the initial values in fortran commons fsiini, led_bldata
-  Double_t GetWeight(const AliHBTPair* partpair); //get weight calculated by Lednicky's algorithm
-
-  void SetTest(Bool_t rtest = kTRUE){ftest = rtest;} //if ftest=0: 
-  //physical values of the following  parameters are put automatically                       
-  //            in FSIINI (their values are not required)          
-  // ftest=1: any values of the following parameters are allowed,    
-  //the following parameters are required:                           
-
-  void SetColoumb(Bool_t col = kTRUE){
-    //(ICH in fortran code) Coulomb interaction between 
-    //the two particles ON (OFF)
-    fColoumbSwitch = col;
-  }
-  void SetQuantumStatistics(Bool_t qss = kTRUE){
-    //IQS: quantum statistics for the two particles ON (OFF)
-    //if non-identical particles automatically off
-    fQuantStatSwitch = qss;
-  }
-  void SetStrongInterSwitch(Bool_t sis = kTRUE){
-    //ISI: strong interaction between the two particles ON (OFF)
-    fStrongInterSwitch = sis;
-  }
-  void SetColWithResidNuclSwitch(Bool_t crn = kTRUE){
-    //I3C: Coulomb interaction with residual nucleus ON (OFF)
-    fColWithResidNuclSwitch = crn;
-  }
-  void SetApproxModel(Int_t ap){
-    //NS in Fortran code,
-    fApproximationModel=ap;
-  }
-  //   NS=1  Square well potential,
-  //   NS=3  not used
-  //   NS=4  scattered wave approximated by the spherical wave,
-  //   NS=2  same as NS=4 but the approx. of equal emission times in PRF
-  //         not required (t=0 approx. used in all other cases).    
+ {
+   public:
+     virtual ~AliHBTLLWeights(){;}
+     static AliHBTLLWeights* Instance();
+     
+     void Init(); //put the initial values in fortran commons fsiini, led_bldata
+     Double_t GetWeight(const AliHBTPair* partpair); //get weight calculated by Lednicky's algorithm
+
+     void SetTest(Bool_t rtest = kTRUE){ftest = rtest;} //if ftest=0: 
+     //physical values of the following  parameters are put automatically                       
+     //            in FSIINI (their values are not required)          
+     // ftest=1: any values of the following parameters are allowed,    
+     //the following parameters are required:                           
+
+     void SetColoumb(Bool_t col = kTRUE){fColoumbSwitch = col;}//: (ICH in fortran code) Coulomb interaction between the two particles ON (OFF)
+     void SetQuantumStatistics(Bool_t qss = kTRUE){fQuantStatSwitch = qss;}//IQS: quantum statistics for the two particles ON (OFF) //if non-identical particles automatically off
+     void SetStrongInterSwitch(Bool_t sis = kTRUE){fStrongInterSwitch = sis;}//ISI: strong interaction between the two particles ON (OFF)
+     void SetColWithResidNuclSwitch(Bool_t crn = kTRUE){fColWithResidNuclSwitch = crn;}//I3C: Coulomb interaction with residual nucleus ON (OFF)  
+     void SetApproxModel(Int_t ap){approximationModel=ap;}//NS in Fortran code, 
+     //   NS=1  Square well potential,                                                             
+     //   NS=3  not used                                                                           
+     //   NS=4  scattered wave approximated by the spherical wave,                                 
+     //   NS=2  same as NS=4 but the approx. of equal emission times in PRF                        
+     //         not required (t=0 approx. used in all other cases).      
+
+     
+     void SetRandomPosition(Bool_t rp = kTRUE){fRandomPosition = rp;} //ON=kTRUE(OFF=kFALSE)
+     // ON -- calculation of the Gauss source radii if the generator don't allows the source generation (for example MeVSim)
+     //if ON the following parameters are requested:
+     void SetR1dw(Double_t R){fRadius=R;}   //spherical source model radii                                                                                                                                                                                             
+     void SetLambdaw(Double_t la){flambda=la;}  //lambda=haoticity                                                   
 
      
-  void SetRandomPosition(Bool_t rp = kTRUE){
-    //ON=kTRUE(OFF=kFALSE)
-    fRandomPosition = rp;
-  } 
-  // ON -- calculation of the Gauss source radii if the generator 
-  //       don't allows the source generation (for example MeVSim)
-  //if ON the following parameters are requested:
-  void SetR1dw(Double_t R){fRadius=R;}   //spherical source model radii
-  void SetLambdaw(Double_t la){flambda=la;}  //lambda=haoticity
-  void SetParticlesTypes(Int_t pid1, Int_t pid2){
-    //set AliRoot particles types   
-    fPID1 = pid1; fPID2 = pid2;
-  }
+     void SetParticlesTypes(Int_t pid1, Int_t pid2){fPID1 = pid1; fPID2 = pid2;} //set AliRoot particles types   
     
-  void SetNucleusCharge(Double_t ch){
-    // not used now  (see comments in fortran code)
-    fNuclCharge=ch;
-  }
-  void SetNucleusMass(Double_t mass){
-    // (see comments in fortran code)
-    fNuclMass=mass;
-  }
-  
-  
- protected:
-  
-  Bool_t ftest; // Switch for automatic setting of all parameters
-  Bool_t fColoumbSwitch; // Swith for Couloumb interaction in the pair
-  Bool_t fQuantStatSwitch; //Switch for quantum statistics
-  Bool_t fStrongInterSwitch;//Switches strong interactions TRUE=ON
-  Bool_t fColWithResidNuclSwitch;//Switches couloumb interaction 
-  //with residual nucleus TRUE=ON          
-  Double_t fNuclMass; //mass 
-  Double_t fNuclCharge; //charge       
-  
-  Bool_t  fRandomPosition; // Radius of Gaussian source
-  Double_t fRadius; // Raduis of spheric source
-  Double_t flambda; // Chaoticity
-  
-  
-  //  Double_t fWein;
-  
-  Int_t fApproximationModel; //approximation used to calculate 
-  //  Bethe-Salpeter amplitude
-  //   ==1  Square well potential,
-  //   ==3  not used
-  //   ==4  scattered wave approximated by the spherical wave,
-  //   ==2  same as NS=4 but the approx. of equal emission times in PRF
-  //         not required (t=0 approx. used in all other cases).
-  //  Note: if ==2,4, the B-S amplitude diverges at zero distance r* in
-  //         the two-particle c.m.s.; user can specify a cutoff AA in
-  //         SUBROUTINE FSIINI, for example:
-  //         IF(NS.EQ.2.OR.NS.EQ.4)AA=5.D0 !! in 1/GeV --> AA=1. fm
-  
-  Int_t fPID1; // Type of the first particle
-  Int_t fPID2; // Type of the second particle
-  
-  static  AliHBTLLWeights *fgLLWeights;// pointer to wrapper of Fortran Lednicky code
-  
-  
-  static Int_t GetPairCode(Int_t pid1,Int_t pid2);
-  static Int_t GetPairCode(const AliHBTPair* partpair);//calculate automatically internal FSIW
-  //----------------------------------------------------------------------
-  // LL        1  2  3  4   5    6   7  8  9 10  11  12  13  14 15 16 17
-  // part. 1:  n  p  n alfa pi+ pi0 pi+ n  p pi+ pi+ pi+ pi- K+ K+ K+ K-
-  // part. 2:  n  p  p alfa pi- pi0 pi+ d  d  K-  K+  p   p  K- K+ p  p
-  // NS=1 y/n: +  +  +  +   +    -   -  -  -  -   -   -   -  -  -  -  -
-  //----------------------------------------------------------------------
-  // LL       18  19 20  21  22 23 24 25 26    27     28
-  // part. 1:  d  d   t  t   K0 K0  d p  p      p      n
-  // part. 2:  d alfa t alfa K0 K0b t t alfa lambda lambda
-  // NS=1 y/n:  -  -   -  -   -  -   - -  -      +      +
-  //----------------------------------------------------------------------
-  
-  Double_t fsigma; //constants for spherical source model 
-  Double_t fwcons; //weight of final state interaction?
-  
- private:
-  AliHBTLLWeights();
-  
-  ClassDef(AliHBTLLWeights,1)
+     void SetNucleusCharge(Double_t ch){fNuclCharge=ch;} // not used now  (see comments in fortran code)
+     void SetNucleusMass(Double_t mass){fNuclMass=mass;} // (see comments in fortran code)
+
+
+   protected:
+     
+     Bool_t ftest; 
+     Bool_t fColoumbSwitch; 
+     Bool_t fQuantStatSwitch; 
+     Bool_t fStrongInterSwitch;//Switches strong interactions TRUE=ON
+     Bool_t fColWithResidNuclSwitch;//Switches couloumb interaction 
+                                    //with residual nucleus TRUE=ON          
+     Double_t fNuclMass; //mass 
+     Double_t fNuclCharge; //charge    
+     
+     Bool_t  fRandomPosition;
+     Double_t fRadius;
+     Double_t flambda;
+
+     
+     Double_t wein;
+
+     Int_t approximationModel; //approximation used to calculate Bethe-Salpeter amplitude
+                             //   ==1  Square well potential,
+                             //   ==3  not used
+                             //   ==4  scattered wave approximated by the spherical wave,
+                             //   ==2  same as NS=4 but the approx. of equal emission times in PRF
+                             //         not required (t=0 approx. used in all other cases).
+                             //  Note: if ==2,4, the B-S amplitude diverges at zero distance r* in
+                             //         the two-particle c.m.s.; user can specify a cutoff AA in
+                             //         SUBROUTINE FSIINI, for example:
+                             //         IF(NS.EQ.2.OR.NS.EQ.4)AA=5.D0 !! in 1/GeV --> AA=1. fm
+     
+     Int_t fPID1;
+     Int_t fPID2;
+
+     static  AliHBTLLWeights *fgLLWeights;// pointer to wrapper of Fortran Lednicky code
+
+
+     static Int_t GetPairCode(Int_t pid1,Int_t pid2);
+     static Int_t GetPairCode(const AliHBTPair* partpair);//calculate automatically internal FSIW 
+         //     C----------------------------------------------------------------------               
+         //     C-   LL       1  2  3  4   5    6   7  8  9 10  11  12  13  14 15 16 17               
+         //     C-   part. 1: n  p  n alfa pi+ pi0 pi+ n  p pi+ pi+ pi+ pi- K+ K+ K+ K-               
+         //     C-   part. 2: n  p  p alfa pi- pi0 pi+ d  d  K-  K+  p   p  K- K+ p  p                
+         //     C   NS=1 y/n: +  +  +  +   +    -   -  -  -  -   -   -   -  -  -  -  -                
+         //     C----------------------------------------------------------------------               
+         //     C-   LL       18  19 20  21  22 23 24 25 26    27     28                              
+         //     C-   part. 1:  d  d   t  t   K0 K0  d p  p      p      n                              
+         //     C-   part. 2:  d alfa t alfa K0 K0b t t alfa lambda lambda                            
+         //     C   NS=1 y/n:  -  -   -  -   -  -   - -  -      +      +                              
+         //     C----------------------------------------------------------------------               
+                                                                                  
+
+     Double_t fsigma; //constants for spherical source model 
+     Double_t fwcons; //
+     
+   private:
+   AliHBTLLWeights();
+
+   public:
+     ClassDef(AliHBTLLWeights,1)
  };
 
 #endif