remove compiler warnings
authorsnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 19 Oct 2010 10:21:41 +0000 (10:21 +0000)
committersnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 19 Oct 2010 10:21:41 +0000 (10:21 +0000)
PWG2/FLOW/AliFlowTools/glauberMC/AliGlauberMC.cxx
PWG2/FLOW/AliFlowTools/glauberMC/AliGlauberMC.h
PWG2/FLOW/AliFlowTools/glauberMC/AliGlauberNucleus.cxx
PWG2/FLOW/AliFlowTools/glauberMC/AliGlauberNucleus.h

index 212d0c0..1e5bb93 100644 (file)
@@ -43,11 +43,14 @@ ClassImp(AliGlauberMC)
 
 //______________________________________________________________________________
 AliGlauberMC::AliGlauberMC(Option_t* NA, Option_t* NB, Double_t xsect) :
+  TNamed(),
   fANucleus(NA),
   fBNucleus(NB),
-  fXSect(0),
+  fXSect(xsect),
   fNucleonsA(0),
   fNucleonsB(0),
+  fAN(0),
+  fBN(0),
   fnt(0),
   fMeanX2(0),
   fMeanY2(0),
@@ -68,20 +71,20 @@ AliGlauberMC::AliGlauberMC(Option_t* NA, Option_t* NB, Double_t xsect) :
   fB_MC(0),
   fEvents(0),
   fTotalEvents(0),
-  fBMin(0),
-  fBMax(0),
+  fBMin(0.),
+  fBMax(20.),
   fMaxNpartFound(0),
   fNpart(0),
   fNcoll(0),
-  fSx2(0),
-  fSy2(0),
-  fSxy(0),
+  fSx2(0.),
+  fSy2(0.),
+  fSxy(0.),
+  fSx2Coll(0.),
+  fSy2Coll(0.),
+  fSxyColl(0.),
   fX(0.13),
   fNpp(8.)
 {
-  fBMin = 0;
-  fBMax = 20;
-  fXSect = xsect;
   fdNdEtaParam[0] = 8.0;
   fdNdEtaParam[1] = 0.13;
   fdNdEtaGBWParam[0] = 0.79;
@@ -90,17 +93,124 @@ AliGlauberMC::AliGlauberMC(Option_t* NA, Option_t* NB, Double_t xsect) :
   fdNdEtaNBDParam[0] = 3.0;
   fdNdEtaNBDParam[1] = 4.0;
   fdNdEtaNBDParam[2] = 0.13;
-  fdNdEtaTwoNBDParam[0] = 0.4;
-  fdNdEtaTwoNBDParam[1] = 3.0;
-  fdNdEtaTwoNBDParam[2] = 4.0;
-  fdNdEtaTwoNBDParam[3] = 2.0;
-  fdNdEtaTwoNBDParam[4] = 11.0;
+  fdNdEtaTwoNBDParam[0] = 3.0;
+  fdNdEtaTwoNBDParam[1] = 4.0;
+  fdNdEtaTwoNBDParam[2] = 2.0;
+  fdNdEtaTwoNBDParam[3] = 11.0;
+  fdNdEtaTwoNBDParam[4] = 0.4;
   fdNdEtaTwoNBDParam[5] = 0.13;
 
-  TString name(Form("Glauber_%s_%s",fANucleus.GetName(),fBNucleus.GetName()));
-  TString title(Form("Glauber %s+%s Version",fANucleus.GetName(),fBNucleus.GetName()));
-  SetName(name);
-  SetTitle(title);
+  SetName(Form("Glauber_%s_%s",fANucleus.GetName(),fBNucleus.GetName()));
+  SetTitle(Form("Glauber %s+%s Version",fANucleus.GetName(),fBNucleus.GetName()));
+}
+
+//______________________________________________________________________________
+AliGlauberMC::~AliGlauberMC()
+{
+  //dtor
+  delete fnt;
+}
+
+//______________________________________________________________________________
+AliGlauberMC::AliGlauberMC(const AliGlauberMC& in):
+  TNamed(in),
+  fANucleus(in.fANucleus),
+  fBNucleus(in.fBNucleus),
+  fXSect(in.fXSect),
+  fNucleonsA(in.fNucleonsA),
+  fNucleonsB(in.fNucleonsB),
+  fAN(in.fAN),
+  fBN(in.fBN),
+  fnt(in.fnt),
+  fMeanX2(in.fMeanX2),
+  fMeanY2(in.fMeanY2),
+  fMeanXY(in.fMeanXY),
+  fMeanXParts(in.fMeanXParts),
+  fMeanYParts(in.fMeanYParts),
+  fMeanXColl(in.fMeanXColl),
+  fMeanYColl(in.fMeanYColl),
+  fMeanX2Coll(in.fMeanX2Coll),
+  fMeanY2Coll(in.fMeanY2Coll),
+  fMeanXYColl(in.fMeanXYColl),
+  fMeanXSystem(in.fMeanXSystem),
+  fMeanYSystem(in.fMeanYSystem),
+  fMeanX_A(in.fMeanX_A),
+  fMeanY_A(in.fMeanY_A),
+  fMeanX_B(in.fMeanX_B),
+  fMeanY_B(in.fMeanY_B),
+  fB_MC(in.fB_MC),
+  fEvents(in.fEvents),
+  fTotalEvents(in.fTotalEvents),
+  fBMin(in.fBMin),
+  fBMax(in.fBMax),
+  fMaxNpartFound(in.fMaxNpartFound),
+  fNpart(in.fNpart),
+  fNcoll(in.fNcoll),
+  fSx2(in.fSx2),
+  fSy2(in.fSy2),
+  fSxy(in.fSxy),
+  fSx2Coll(in.fSx2Coll),
+  fSy2Coll(in.fSy2Coll),
+  fSxyColl(in.fSxyColl),
+  fX(in.fX),
+  fNpp(in.fNpp)
+{
+  //copy ctor
+  memcpy(fdNdEtaParam,in.fdNdEtaParam,2*sizeof(Double_t));
+  memcpy(fdNdEtaGBWParam,in.fdNdEtaGBWParam,3*sizeof(Double_t));
+  memcpy(fdNdEtaNBDParam,in.fdNdEtaNBDParam,3*sizeof(Double_t));
+  memcpy(fdNdEtaTwoNBDParam,in.fdNdEtaTwoNBDParam,6*sizeof(Double_t));
+}
+
+//______________________________________________________________________________
+AliGlauberMC& AliGlauberMC::operator=(const AliGlauberMC& in)
+{
+  //assignment
+  fANucleus=in.fANucleus; 
+  fBNucleus=in.fBNucleus; 
+  fXSect=in.fXSect;    
+  fNucleonsA=in.fNucleonsA;
+  fNucleonsB=in.fNucleonsB;
+  fAN=in.fAN;       
+  fBN=in.fBN;       
+  fnt=in.fnt;       
+  fMeanX2=in.fMeanX2;   
+  fMeanY2=in.fMeanY2;   
+  fMeanXY=in.fMeanXY;   
+  fMeanXParts=in.fMeanXParts;
+  fMeanYParts=in.fMeanYParts;
+  fMeanXColl=in.fMeanXColl;
+  fMeanYColl=in.fMeanYColl;
+  fMeanX2Coll=in.fMeanX2Coll;
+  fMeanY2Coll=in.fMeanY2Coll;
+  fMeanXYColl=in.fMeanXYColl;
+  fMeanXSystem=in.fMeanXSystem;
+  fMeanYSystem=in.fMeanYSystem;
+  fMeanX_A=in.fMeanX_A;  
+  fMeanY_A=in.fMeanY_A;  
+  fMeanX_B=in.fMeanX_B;  
+  fMeanY_B=in.fMeanY_B;  
+  fB_MC=in.fB_MC;     
+  fEvents=in.fEvents;   
+  fTotalEvents=in.fTotalEvents;
+  fBMin=in.fBMin;     
+  fBMax=in.fBMax;     
+  memcpy(fdNdEtaParam,in.fdNdEtaParam,2*sizeof(Double_t));
+  memcpy(fdNdEtaGBWParam,in.fdNdEtaGBWParam,3*sizeof(Double_t));
+  memcpy(fdNdEtaNBDParam,in.fdNdEtaNBDParam,3*sizeof(Double_t));
+  memcpy(fdNdEtaTwoNBDParam,in.fdNdEtaTwoNBDParam,6*sizeof(Double_t));
+  fMaxNpartFound=in.fMaxNpartFound;
+  fNpart=in.fNpart;   
+  fNcoll=in.fNcoll;    
+  fSx2=in.fSx2;      
+  fSy2=in.fSy2;      
+  fSxy=in.fSxy;      
+  fSx2Coll=in.fSx2Coll;  
+  fSy2Coll=in.fSy2Coll;  
+  fSxyColl=in.fSxyColl;  
+  fX=in.fX;      
+  fNpp=in.fNpp;       
+  return *this;
 }
 
 //______________________________________________________________________________
@@ -419,14 +529,19 @@ void AliGlauberMC::SetdNdEtaNBDParam(Double_t k, Double_t nmean, Double_t beta)
   fdNdEtaNBDParam[2]=beta;
 }
 //______________________________________________________________________________
-void AliGlauberMC::SetdNdEtaTwoNBDParam(Double_t alpha, Double_t k1, Double_t nmean1, Double_t k2, Double_t nmean2, Double_t beta)
+void AliGlauberMC::SetdNdEtaTwoNBDParam( Double_t k1, 
+                                         Double_t nmean1, 
+                                         Double_t k2, 
+                                         Double_t nmean2, 
+                                         Double_t alpha,
+                                         Double_t beta)
 {
   // set parameters used for calculating multiplicity see GetdNdEtaTwoNBD() for commments
-  fdNdEtaTwoNBDParam[0]=alpha;
-  fdNdEtaTwoNBDParam[1]=k1;
-  fdNdEtaTwoNBDParam[2]=nmean1;
-  fdNdEtaTwoNBDParam[3]=k2;
-  fdNdEtaTwoNBDParam[4]=nmean2;
+  fdNdEtaTwoNBDParam[0]=k1;
+  fdNdEtaTwoNBDParam[1]=nmean1;
+  fdNdEtaTwoNBDParam[2]=k2;
+  fdNdEtaTwoNBDParam[3]=nmean2;
+  fdNdEtaTwoNBDParam[4]=alpha;
   fdNdEtaTwoNBDParam[5]=beta;
 }
 //______________________________________________________________________________
@@ -472,13 +587,12 @@ Double_t AliGlauberMC::GetdNdEtaNBD ( Int_t k, Double_t nmean, Double_t beta)
 }
 //______________________________________________________________________________
 
-Double_t AliGlauberMC::GetdNdEtaTwoNBD (
-  Int_t k1,
-  Double_t nmean1,
-  Int_t k2,
-  Double_t nmean2,
-  Double_t alpha,
-  Double_t beta )
+Double_t AliGlauberMC::GetdNdEtaTwoNBD ( Int_t k1,
+                                         Double_t nmean1,
+                                         Int_t k2,
+                                         Double_t nmean2,
+                                         Double_t alpha,
+                                         Double_t beta )
 {
   //Get particle density per unit of rapidity
   //using random numbers from two negative binomial distributions
@@ -617,10 +731,10 @@ void AliGlauberMC::runAndSaveNtuple( Int_t n,
                                      Double_t mind,
                                      const char *fname)
 {
-  AliGlauberMC *mcg=new AliGlauberMC(sysA,sysB,signn);
-  mcg->SetMinDistance(mind);
-  mcg->Run(n);
-  TNtuple  *nt=mcg->GetNtuple();
+  AliGlauberMC mcg(sysA,sysB,signn);
+  mcg.SetMinDistance(mind);
+  mcg.Run(n);
+  TNtuple  *nt=mcg.GetNtuple();
   TFile out(fname,"recreate",fname,9);
   if(nt) nt->Write();
   out.Close();
@@ -635,8 +749,8 @@ void AliGlauberMC::runAndSaveNucleons( Int_t n,
                                        Bool_t verbose,
                                        const char *fname)
 {
-  AliGlauberMC *mcg=new AliGlauberMC(sysA,sysB,signn);
-  mcg->SetMinDistance(mind);
+  AliGlauberMC mcg(sysA,sysB,signn);
+  mcg.SetMinDistance(mind);
   TFile *out=0;
   if(fname)
     out=new TFile(fname,"recreate",fname,9);
@@ -645,10 +759,10 @@ void AliGlauberMC::runAndSaveNucleons( Int_t n,
   {
 
     //get an event with at least one collision
-    while(!mcg->NextEvent()) {}
+    while(!mcg.NextEvent()) {}
 
     //access, save and (if wanted) print out nucleons
-    TObjArray* nucleons=mcg->GetNucleons();
+    TObjArray* nucleons=mcg.GetNucleons();
     if(!nucleons) continue;
     if(out)
       nucleons->Write(Form("nucleonarray%d",ievent),TObject::kSingleKey);
@@ -656,7 +770,7 @@ void AliGlauberMC::runAndSaveNucleons( Int_t n,
     if(verbose)
     {
       cout<<endl<<endl<<"EVENT NO: "<<ievent<<endl;
-      cout<<"B = "<<mcg->GetB()<<"  Npart = "<<mcg->GetNpart()<<endl<<endl;
+      cout<<"B = "<<mcg.GetB()<<"  Npart = "<<mcg.GetNpart()<<endl<<endl;
       printf("Nucleus\t X\t Y\t Z\tNcoll\n");
       Int_t nNucls=nucleons->GetEntries();
       for(Int_t iNucl=0; iNucl<nNucls; iNucl++)
@@ -675,3 +789,9 @@ void AliGlauberMC::runAndSaveNucleons( Int_t n,
   if(out) delete out;
 }
 
+//---------------------------------------------------------------------------------
+void AliGlauberMC::Reset()
+{
+  //delete the ntuple
+  delete fnt; fnt=NULL;
+}
index 5ef95da..5898b67 100644 (file)
@@ -52,10 +52,10 @@ private:
    Int_t        fTotalEvents;    //All events within selected impact parameter range
    Double_t     fBMin;           //Minimum impact parameter to be generated
    Double_t     fBMax;           //Maximum impact parameter to be generated
-   Double_t    fdNdEtaParam[2];           //Parameters: nnp, x
-   Double_t fdNdEtaGBWParam[3];  //Parameters: delta, lambda, snn
-   Double_t fdNdEtaNBDParam[3];       //Parameters:  k, nmean, beta
-   Double_t fdNdEtaTwoNBDParam[6];    //Parameters: alpha, k1, nmean1, k2, nmean2, beta
+   Double_t        fdNdEtaParam[2];       //Parameters: nnp, x
+   Double_t     fdNdEtaGBWParam[3];  //Parameters: delta, lambda, snn
+   Double_t     fdNdEtaNBDParam[3];       //Parameters:  k, nmean, beta
+   Double_t     fdNdEtaTwoNBDParam[6];    //Parameters: k1, nmean1, k2, nmean2, alpha, beta
    Int_t        fMaxNpartFound;  //Largest value of Npart obtained
    Int_t        fNpart;          //Number of wounded (participating) nucleons in current event
    Int_t        fNcoll;          //Number of binary collisions in current event
@@ -71,7 +71,9 @@ private:
 
 public:
    AliGlauberMC(Option_t* NA = "Pb", Option_t* NB = "Pb", Double_t xsect = 72);
-   virtual     ~AliGlauberMC() {Reset();}
+   virtual     ~AliGlauberMC();
+   AliGlauberMC(const AliGlauberMC& in);
+   AliGlauberMC& operator=(const AliGlauberMC& in);
    void         Draw(Option_t* option);
 
    void         Run(Int_t nevents);
@@ -108,7 +110,7 @@ public:
    TObjArray   *GetNucleons();
    Double_t     GetTotXSect()        const;
    Double_t     GetTotXSectErr()     const;
-   void         Reset()                    {delete fnt; fnt=0; }
+   void         Reset();
    static Double_t     NegativeBinomialDistribution(Int_t x, Int_t k, Double_t nmean);
    Int_t NegativeBinomialRandom(Int_t k, Double_t nmean);
    Int_t DoubleNegativeBinomialRandom(Int_t k1, Double_t nmean1, Int_t k2, Double_t nmean2, Double_t alpha);
index d6d2c8b..0c25c3d 100644 (file)
@@ -37,9 +37,16 @@ ClassImp(AliGlauberNucleus)
 
 //______________________________________________________________________________
 AliGlauberNucleus::AliGlauberNucleus(Option_t* iname, Int_t iN, Double_t iR, Double_t ia, Double_t iw, TF1* ifunc) : 
-   fN(iN),fR(iR),fA(ia),fW(iw),fMinDist(-1),
-   fF(0),fTrials(0),fFunction(ifunc),
-   fNucleons(0)
+  TNamed(iname,""),
+  fN(iN),
+  fR(iR),
+  fA(ia),
+  fW(iw),
+  fMinDist(-1),
+  fF(0),
+  fTrials(0),
+  fFunction(ifunc),
+  fNucleons(NULL)
 {
    if (fN==0) {
       cout << "Setting up nucleus " << iname << endl;
@@ -57,6 +64,42 @@ AliGlauberNucleus::~AliGlauberNucleus()
 }
 
 //______________________________________________________________________________
+AliGlauberNucleus::AliGlauberNucleus(const AliGlauberNucleus& in):
+  TNamed(in),
+  fN(in.fN),
+  fR(in.fR),
+  fA(in.fA),
+  fW(in.fW),
+  fMinDist(in.fMinDist),
+  fF(in.fF),
+  fTrials(in.fTrials),
+  fFunction(in.fFunction),
+  fNucleons(NULL)
+{
+  //copy ctor
+  if (in.fNucleons)
+    fNucleons=static_cast<TObjArray*>((in.fNucleons)->Clone());
+}
+
+//______________________________________________________________________________
+AliGlauberNucleus& AliGlauberNucleus::operator=(const AliGlauberNucleus& in)
+{
+  //assignment
+  fN=in.fN;
+  fR=in.fR;
+  fA=in.fA;
+  fW=in.fW;
+  fMinDist=in.fMinDist;
+  fF=in.fF;
+  fTrials=in.fTrials;
+  fFunction=in.fFunction;
+  delete fNucleons;
+  fNucleons=static_cast<TObjArray*>((in.fNucleons)->Clone());
+  fNucleons->SetOwner();
+  return *this;
+}
+
+//______________________________________________________________________________
 void AliGlauberNucleus::Draw(Double_t xs, Int_t col)
 {
    Double_t r = 0.5*sqrt(xs/TMath::Pi()/10.);
index 473556e..09ae913 100644 (file)
@@ -36,6 +36,8 @@ private:
 public:
    AliGlauberNucleus(Option_t* iname="Au", Int_t iN=0, Double_t iR=0, Double_t ia=0, Double_t iw=0, TF1* ifunc=0);
    virtual ~AliGlauberNucleus();
+   AliGlauberNucleus(const AliGlauberNucleus& in);
+   AliGlauberNucleus& operator=(const AliGlauberNucleus& in);
 
    using      TObject::Draw;
    void       Draw(Double_t xs, Int_t col);