From 7836ee9457a76c09d25c8dee306f8e337f4d1fb8 Mon Sep 17 00:00:00 2001 From: skowron Date: Thu, 10 Jan 2002 18:37:33 +0000 Subject: [PATCH] Bug correction: files closed to early --- HBTAN/.rootrc | 2 +- HBTAN/AliHBTAnalysis.cxx | 10 ++++++ HBTAN/AliHBTFunction.cxx | 71 +++++++++++++++++++++++++++++++++++-- HBTAN/AliHBTFunction.h | 50 +++++++++++++++++++++++++- HBTAN/AliHBTPair.cxx | 8 ++--- HBTAN/AliHBTPair.h | 36 +++++++++++++++++++ HBTAN/AliHBTReaderITSv2.cxx | 4 +-- HBTAN/HBTAnalysisLinkDef.h | 5 +++ HBTAN/Makefile | 3 +- HBTAN/libHBTAN.pkg | 3 +- 10 files changed, 179 insertions(+), 13 deletions(-) diff --git a/HBTAN/.rootrc b/HBTAN/.rootrc index 413ae4c2a45..c85cc3332c8 100644 --- a/HBTAN/.rootrc +++ b/HBTAN/.rootrc @@ -1,7 +1,7 @@ Unix.*.Root.MacroPath: .:$(ALICE_ROOT)/macros:$(ALICE_ROOT) -Root.Html.OutputDir: $(ALICE_ROOT)/html/roothtml/HBTANALYSIS +Root.Html.OutputDir: $(ALICE_ROOT)/html/roothtml/HBTAN Root.Html.SourceDir: ./ Root.Html.Author: Piotr Krzysztof Skowronski Root.Html.Root: http://root.cern.ch/root/html diff --git a/HBTAN/AliHBTAnalysis.cxx b/HBTAN/AliHBTAnalysis.cxx index eba0546685c..56e4bc76c90 100644 --- a/HBTAN/AliHBTAnalysis.cxx +++ b/HBTAN/AliHBTAnalysis.cxx @@ -113,12 +113,22 @@ void AliHBTAnalysis::Process(Option_t* option) if(oT) { + if (fReader->GetNumberOfTrackEvents() <1) + { + Error("Process","There is no data to analyze."); + return; + } ProcessTracks(); return; } if(oP) { + if (fReader->GetNumberOfPartEvents() <1) + { + Error("Process","There is no data to analyze."); + return; + } ProcessParticles(); return; } diff --git a/HBTAN/AliHBTFunction.cxx b/HBTAN/AliHBTFunction.cxx index 1107ab56174..04e2806a84c 100644 --- a/HBTAN/AliHBTFunction.cxx +++ b/HBTAN/AliHBTFunction.cxx @@ -32,8 +32,13 @@ ClassImp( AliHBTFunction ) AliHBTFunction::AliHBTFunction() { - - fPairCut = new AliHBTEmptyPairCut(); //dummy cut + fPairCut = new AliHBTEmptyPairCut(); //dummy cut +} +/******************************************************************/ +AliHBTFunction::AliHBTFunction(const char* name,const char* title) +{ + fPairCut = new AliHBTEmptyPairCut(); //dummy cut + Rename(name,title); } /******************************************************************/ @@ -350,15 +355,71 @@ AliHBTOnePairFctn3D(Int_t nXbins, Double_t maxXval, Double_t minXval, fDenominator->Sumw2(); } - +/******************************************************************/ AliHBTOnePairFctn3D::~AliHBTOnePairFctn3D() { delete fNumerator; delete fDenominator; } +/******************************************************************/ +/******************************************************************/ +/******************************************************************/ +/******************************************************************/ +ClassImp( AliHBTTwoPairFctn1D) + +AliHBTTwoPairFctn1D:: +AliHBTTwoPairFctn1D(Int_t nbins, Double_t maxval, Double_t minval) + { + TString numstr = fName + " Numerator"; //title and name of the + //numerator histogram + TString denstr = fName + " Denominator";//title and name of the + //denominator histogram + + fNumerator = new TH1D(numstr.Data(),numstr.Data(), + nbins,minval,maxval); + + fDenominator = new TH1D(denstr.Data(),denstr.Data(), + nbins,minval,maxval); + + fNumerator->Sumw2(); + fDenominator->Sumw2(); + } + +/******************************************************************/ +AliHBTTwoPairFctn1D::~AliHBTTwoPairFctn1D() +{ + delete fNumerator; + delete fDenominator; +} +void AliHBTTwoPairFctn1D:: +ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair) +{ + partpair = CheckPair(partpair); + trackpair = CheckPair(trackpair); + if( partpair && trackpair) + { + Double_t x = GetValue(trackpair,partpair); + fNumerator->Fill(x); + } +} +/******************************************************************/ + +void AliHBTTwoPairFctn1D:: +ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair) +{ + partpair = CheckPair(partpair); + trackpair = CheckPair(trackpair); + if( partpair && trackpair) + { + Double_t x = GetValue(trackpair,partpair); + fDenominator->Fill(x); + } + +} + /******************************************************************/ /******************************************************************/ /******************************************************************/ @@ -419,3 +480,7 @@ ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair) } +/******************************************************************/ +/******************************************************************/ +/******************************************************************/ +ClassImp(AliHBTTwoPairFctn3D) diff --git a/HBTAN/AliHBTFunction.h b/HBTAN/AliHBTFunction.h index 9fc75fd8231..44140565d45 100644 --- a/HBTAN/AliHBTFunction.h +++ b/HBTAN/AliHBTFunction.h @@ -18,6 +18,7 @@ class AliHBTFunction: public TNamed { public: AliHBTFunction(); + AliHBTFunction(const char* name,const char* title); virtual ~AliHBTFunction(); virtual TH1* GetNumerator() =0; @@ -179,6 +180,29 @@ class AliHBTOnePairFctn3D: public AliHBTOnePairFctn /******************************************************************/ /******************************************************************/ +class AliHBTTwoPairFctn1D: public AliHBTTwoPairFctn +{ + public: + AliHBTTwoPairFctn1D(Int_t nbins = 200, Double_t maxval = 1.5, Double_t minval = 0.0); + AliHBTTwoPairFctn1D(const char*,const char*, + Int_t nbins = 200, Double_t maxval = 1.5, Double_t minval = 0.0); + ~AliHBTTwoPairFctn1D(); + + TH1* GetNumerator(){return fNumerator;} + TH1* GetDenominator(){return fDenominator;} + + void ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair); + void ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair); + + protected: + virtual Double_t GetValue(AliHBTPair* trackpair, AliHBTPair* partpair) = 0; + + TH1D* fNumerator; + TH1D* fDenominator; + + public: + ClassDef(AliHBTTwoPairFctn1D,1) +}; /******************************************************************/ @@ -189,7 +213,7 @@ class AliHBTTwoPairFctn2D: public AliHBTTwoPairFctn public: AliHBTTwoPairFctn2D(Int_t nXbins = 200, Double_t maxXval = 1.5, Double_t minXval = 0.0, Int_t nYbins = 200, Double_t maxYval = .15, Double_t minYval =-0.15); - ~AliHBTTwoPairFctn2D(); + virtual ~AliHBTTwoPairFctn2D(); TH1* GetNumerator(){return fNumerator;} TH1* GetDenominator(){return fDenominator;} @@ -212,6 +236,30 @@ class AliHBTTwoPairFctn2D: public AliHBTTwoPairFctn /******************************************************************/ /******************************************************************/ /******************************************************************/ +class AliHBTTwoPairFctn3D: public AliHBTTwoPairFctn +{ + public: + AliHBTTwoPairFctn3D(Int_t nXbins = 200, Double_t maxXval = 1.5, Double_t minXval = 0.0, + Int_t nYbins = 200, Double_t maxYval = .15, Double_t minYval =-0.15, + Int_t nZbins = 200, Double_t maxZval = .15, Double_t minZval =-0.15){} + virtual ~AliHBTTwoPairFctn3D(){} + + TH1* GetNumerator(){return fNumerator;} + TH1* GetDenominator(){return fDenominator;} + + void ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair){} + void ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair){} + + + protected: + virtual void GetValues(AliHBTPair*,AliHBTPair*, Double_t&, Double_t&,Double_t&) = 0; + + TH3D* fNumerator; + TH3D* fDenominator; + + public: + ClassDef(AliHBTTwoPairFctn3D,1) +}; /******************************************************************/ /******************************************************************/ diff --git a/HBTAN/AliHBTPair.cxx b/HBTAN/AliHBTPair.cxx index 732f004692b..6d05d09fe67 100644 --- a/HBTAN/AliHBTPair.cxx +++ b/HBTAN/AliHBTPair.cxx @@ -98,7 +98,7 @@ Double_t AliHBTPair::GetKStar() CalculateQInvL(); Q = TMath::Sqrt( Q*Q - fQInvL); - fKStar = Q/2; + fKStar = Q/2.; fKStarNotCalc = kFALSE; } return fKStar; @@ -120,19 +120,19 @@ Double_t AliHBTPair::GetQInv() Double_t AliHBTPair::GetQSide() { //returns Q side - return -1; + return fPart1->Py() - fPart1->Py(); } Double_t AliHBTPair::GetQLong() { //returns Q long - return -1; + return fPart1->Pz() - fPart1->Pz(); } Double_t AliHBTPair::GetQOut() { //returns Q out - return -1; + return fPart1->Px() - fPart1->Px(); } /************************************************************************/ diff --git a/HBTAN/AliHBTPair.h b/HBTAN/AliHBTPair.h index 9195f2fcadc..10bedc6d7c0 100644 --- a/HBTAN/AliHBTPair.h +++ b/HBTAN/AliHBTPair.h @@ -37,6 +37,10 @@ class AliHBTPair: public TObject Double_t GetQLong(); //returns Q long Double_t GetQOut(); //returns Q out + Double_t GetDeltaP(); //return difference of momenta + Double_t GetDeltaPx(); + Double_t GetDeltaPy(); + Double_t GetDeltaPz(); protected: AliHBTParticle* fPart1; //pointer to first particle @@ -157,6 +161,7 @@ void AliHBTPair::CalculateInvMassSqr() inline void AliHBTPair::CalculateQInvL() { + //Calculates square root of Qinv if (fQInvLNotCalc) { CalculateDiffs(); @@ -191,4 +196,35 @@ void AliHBTPair::CalculateDiffs() } } +/****************************************************************/ +inline +Double_t AliHBTPair::GetDeltaP() //return difference of momenta +{ + CalculateDiffs(); + return TMath::Sqrt(fPxDiff*fPxDiff + fPyDiff*fPyDiff + fPzDiff*fPzDiff); +} +/****************************************************************/ +inline +Double_t AliHBTPair::GetDeltaPx() + { + CalculateDiffs(); + return fPxDiff; + } +/****************************************************************/ +inline +Double_t AliHBTPair::GetDeltaPy() + { + CalculateDiffs(); + return fPyDiff; + } + +/****************************************************************/ +inline +Double_t AliHBTPair::GetDeltaPz() + { + CalculateDiffs(); + return fPzDiff; + } + + #endif diff --git a/HBTAN/AliHBTReaderITSv2.cxx b/HBTAN/AliHBTReaderITSv2.cxx index 7867873a5ae..92687b3003a 100644 --- a/HBTAN/AliHBTReaderITSv2.cxx +++ b/HBTAN/AliHBTReaderITSv2.cxx @@ -295,11 +295,11 @@ Int_t AliHBTReaderITSv2::Read(AliHBTRun* particles, AliHBTRun *tracks) // delete tracker; totalNevents++; - CloseFiles(aTracksFile,aClustersFile,aGAliceFile); cout<<"all: "<