AliRICHhit *pHit=(AliRICHhit*)Hits()->At(iHitN);
TVector2 x2 = Param()->ShiftToWirePos(C(pHit->C())->Glob2Loc(pHit->OutX3()));
Int_t iTotQdc=Param()->TotQdc(x2,pHit->Eloss());
-
+ if(iTotQdc==0) continue;
Int_t iPadXmin,iPadXmax,iPadYmin,iPadYmax;
Param()->Loc2Area(x2,iPadXmin,iPadYmin,iPadXmax,iPadYmax);//determine affected pads
if(GetDebug()) Info("Hits2SDigits","left-down=(%i,%i) right-up=(%i,%i)",iPadXmin,iPadYmin,iPadXmax,iPadYmax);
#include <AliDigit.h>
#include "AliRICHDigitizer.h"
#include "AliRICHParam.h"
+#include <Riostream.h>
//__________________AliRICHhit______________________________________________________________________
class AliRICHhit : public AliHit
{TClonesArray &tmp=*fHits;new(tmp[fNhits++])AliRICHhit(chamber,tid,iX3,oX3,eloss);}
inline void AddSDigit(Int_t c,Int_t x,Int_t y,Double_t q,Int_t pid,Int_t tid);
void AddDigit(int c,int x,int y,int q,int cpid,int *tid){TClonesArray &tmp=*((TClonesArray*)fDigitsNew->At(c-1));new(tmp[fNdigitsNew[c-1]++])AliRICHdigit(c,x,y,q,cpid,tid[0],tid[1],tid[2]);}
- void AddCluster(AliRICHcluster &cl) {TClonesArray &tmp=*((TClonesArray*)fClusters->At(cl.C()-1));new(tmp[fNclusters[cl.C()-1]++])AliRICHcluster(cl);}
+ void AddCluster(AliRICHcluster &cl) {Int_t c=cl.C()-1;cout<<c<<endl;TClonesArray &tmp=*((TClonesArray*)fClusters->At(c));new(tmp[fNclusters[c]++])AliRICHcluster(cl);}
void AddReco(Int_t tid,Double_t thetaCherenkov,Int_t nPhotons) {TClonesArray &tmp=*(TClonesArray*)fRecos;new(tmp[fNrecos++])AliRICHreco(tid,thetaCherenkov,nPhotons);}
// void ResetHits() {AliDetector::ResetHits();} //virtual
void ResetSDigits() {fNsdigits=0; if(fSdigits) fSdigits ->Clear();}
void AliRICHClusterFinder::WriteRawCluster()
{
// out the current raw cluster
-// Info("WriteRawCluster","Start.");
+ Info("WriteRawCluster","Start.");
FindClusterContribs(&fRawCluster);
+ fRawCluster.Dump();
Rich()->AddCluster(fRawCluster);
// fRawCluster.Print();
}//WriteRawCluster()
void AliRICHClusterFinder::WriteResolvedCluster()
{
// out the current resolved cluster
-// Info("WriteResolvedCluster","Start.");
+ Info("WriteResolvedCluster","Start.");
// FindClusterContribs(&fResolvedCluster);
Rich()->AddCluster(fResolvedCluster);
WriteResolvedCluster();
}
if(fNlocals==5) Info("CoG","Stop.");
-}
+}//FitCoG()
//__________________________________________________________________________________________________
void RICHMinMathieson(Int_t &npar, Double_t *, Double_t &chi2, Double_t *par, Int_t )
-{// Minimization function of Mathieson
+{
+// Mathieson minimization function
AliRICHcluster *pRawCluster = ((AliRICHClusterFinder*)gMinuit->GetObjectFit())->GetRawCluster();
ClassImp(AliRICHDigitizer)
//__________________________________________________________________________________________________
-AliRICHDigitizer::AliRICHDigitizer()
-{//default constructor
-}//default ctor
-//__________________________________________________________________________________________________
-AliRICHDigitizer::AliRICHDigitizer(AliRunDigitizer *pManager)
- :AliDigitizer(pManager)
+Bool_t AliRICHDigitizer::Init()
{
-//main ctor which should be used
- Info("main ctor","Start.");
+//This methode is called from AliRunDigitizer after the corresponding file is open
+ if(GetDebug())Info("Init","Start.");
fRich=(AliRICH*)gAlice->GetDetector("RICH");
Rich()->Param()->GenSigmaThMap();
-}//main ctor
-//__________________________________________________________________________________________________
-AliRICHDigitizer::~AliRICHDigitizer()
-{
-//dtor
- Info("dtor","Start.");
-}//dtor
+ return kTRUE;
+}//Init()
//__________________________________________________________________________________________________
void AliRICHDigitizer::Exec(Option_t*)
{
- Info("Exec","\n\n\n\n");
- Info("Exec","Start with %i input(s) for event %i",fManager->GetNinputs(),fManager->GetOutputEventNr());
+ if(GetDebug())Info("Exec","Start with %i input(s) for event %i",fManager->GetNinputs(),fManager->GetOutputEventNr());
AliRunLoader *pInAL=0, *pOutAL;//in and out Run loaders
AliLoader *pInRL=0, *pOutRL;//in and out RICH loaders
pOutAL = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
pOutRL = pOutAL->GetLoader("RICHLoader");
pOutRL->MakeTree("D"); Rich()->MakeBranch("D"); //create TreeD with RICH branches in output stream
-
for(Int_t inFileN=0;inFileN<fManager->GetNinputs();inFileN++){//files loop
pInAL = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(inFileN)); pInRL = pInAL->GetLoader("RICHLoader");
- pInAL->GetEvent(fManager->GetOutputEventNr()); pInRL->LoadSDigits(); pInRL->TreeS()->GetEntry(0);
- Info("Exec","Run Loader %p RICH Loader %p RICH %p",pInAL,pInRL,Rich());
- Rich()->SDigits()->Print();
+ pInRL->LoadSDigits(); pInAL->GetEvent(fManager->GetOutputEventNr()); pInRL->TreeS()->GetEntry(0);
pInRL->UnloadSDigits();
}//files loop
-// Rich()->SDigits()->Sort(); //sort them according to Id() methode
-//
-// Int_t combiPid=0,chamber=0,x=0,y=0,tid[3],id=0; Double_t q=0;
-// Int_t iNdigitsPerPad=0;//how many sdigits for a given pad
-// for(Int_t i=0;i<Rich()->SDigits()->GetEntries();i++){//sdigits loop (sorted)
-// AliRICHdigit *pSdig=(AliRICHdigit*)Rich()->SDigits()->At(i);
-// if(pSdig->Id()==id){//still the same pad
-// iNdigitsPerPad++;
-// q+=pSdig->Q();
-// combiPid+=pSdig->CombiPid();
-// if(iNdigitsPerPad<=3)
-// tid[iNdigitsPerPad-1]=pSdig->Tid(0);
-// else
-// Warning("SDigits2Digits","More then 3 sdigits for the given pad");
-// }else{//new pad, add the pevious one
-// if(id!=kBad&&Rich()->Param()->IsOverTh(chamber,x,y,q)) {
-// Rich()->AddDigit(chamber,x,y,(Int_t)q,combiPid,tid); //add newly created digit to the list of digits
-// }
-// combiPid=pSdig->CombiPid();chamber=pSdig->C();id=pSdig->Id();
-// x=pSdig->X();y=pSdig->Y();
-// q=pSdig->Q();
-// tid[0]=pSdig->Tid(0);
-// iNdigitsPerPad=1;tid[1]=tid[2]=kBad;
-// }
-// }//sdigits loop (sorted)
-//
-// if(Rich()->SDigits()->GetEntries()&&Rich()->Param()->IsOverTh(chamber,x,y,q))
-// Rich()->AddDigit(chamber,x,y,(Int_t)q,combiPid,tid);//add the last digit
-//
+ Rich()->SDigits()->Sort(); //sort them according to Id() methode
+ Int_t combiPid=0,chamber=0,x=0,y=0,tid[3],id=0; Double_t q=0;
+ Int_t iNdigitsPerPad=0;//how many sdigits for a given pad
+ for(Int_t i=0;i<Rich()->SDigits()->GetEntries();i++){//sdigits loop (sorted)
+ AliRICHdigit *pSdig=(AliRICHdigit*)Rich()->SDigits()->At(i);
+ if(pSdig->Id()==id){//still the same pad
+ iNdigitsPerPad++; q+=pSdig->Q(); combiPid+=pSdig->CombiPid();//sum up charge and cfm
+ if(iNdigitsPerPad<=3) tid[iNdigitsPerPad-1]=pSdig->Tid(0);
+ else Warning("SDigits2Digits","More then 3 sdigits for the given pad");
+ }else{//new pad, add the pevious one
+ if(id!=kBad&&Rich()->Param()->IsOverTh(chamber,x,y,q)) Rich()->AddDigit(chamber,x,y,(Int_t)q,combiPid,tid); //add newly created dig
+ combiPid=pSdig->CombiPid(); chamber=pSdig->C(); id=pSdig->Id(); x=pSdig->X(); y=pSdig->Y(); q=pSdig->Q(); //init all values by current sdig
+ iNdigitsPerPad=1; tid[0]=pSdig->Tid(0); tid[1]=tid[2]=kBad;
+ }
+ }//sdigits loop (sorted)
+ if(Rich()->SDigits()->GetEntries()&&Rich()->Param()->IsOverTh(chamber,x,y,q)) Rich()->AddDigit(chamber,x,y,(Int_t)q,combiPid,tid);//add the last dig
pOutRL->TreeD()->Fill(); //fill the tree with the list of digits
pOutRL->WriteDigits("OVERWRITE"); //serialize them to file
- Info("Exec","Stop\n\n\n\n");
+ if(GetDebug())Info("Exec","Stop.");
}//Exec()
//__________________________________________________________________________________________________
class AliRICHDigitizer : public AliDigitizer
{
public:
- AliRICHDigitizer();
- AliRICHDigitizer(AliRunDigitizer * manager);
- virtual ~AliRICHDigitizer();
+ AliRICHDigitizer() {;}
+ AliRICHDigitizer(AliRunDigitizer * manager):AliDigitizer(manager) {if(GetDebug())Info("main ctor","Start.");}
+ virtual ~AliRICHDigitizer() {if(GetDebug())Info("dtor","Start.");}
+
void Exec(Option_t* option=0); //virtual
+ Bool_t Init(); //virtual
Bool_t GetDebug() const {return gAlice->GetDebug();}
AliRICH* Rich() const {return fRich;}
protected:
// Calculates the total charge produced by the eloss in point x2 (Chamber RS).
// Returns this change parametrised in QDC channels.
// eloss=0 means photons which provided for only 1 electron
-// eloss > 0 for Mip
+// eloss > 0 for Mip
+ if(Sector(x2)==kBad) return 0; //hit in the dead zone
Int_t iNelectrons=Int_t(eloss/IonisationPotential()); if(iNelectrons==0) iNelectrons=1;
Double_t qdc=0;
for(Int_t i=1;i<=iNelectrons;i++) qdc+=-Gain(x2)*TMath::Log(gRandom->Rndm());
-pointer to the stack:
-AliRunLoader::Stack()
total amount of particles in stack for a given event:
AliStack::GetNtrack() or AliRun::GetEvent()
total amount of primiry particles in stack for a given event:
AliStack::GetNprimary() or TreeH()::GetEntries()
-to open session
-AliRunLoader::Open(
+How to open session?
+ use static methode AliRunLoader::Open(
-Hits a stored on primiry by primiry basis. To retrieve all hits one needs to do:
+How to retrive hits:
+ Hits a stored on primiry by primiry basis. To retrieve all hits one needs to do:
+ initialise the root tree and containers: AliLoader::LoadHits()
+ read number of primiries in current event:
+ loop on the list of primiries:
-initialise the root tree and containers: AliLoader::LoadHits()
-read number of primiries in current event:
-loop on the list of primiries:
-
-How to retrive sdigits? Sdigits stored in tree S with the branch of TClonesArray, all sdigits in a single TClonesArray
- So the tree has only one entry.
- One needs to say:
- pRich->GetLoader()->LoadSDigits(); this one open file, get the tree and invoke AliRICH::SetTreeAddress()
+How to retrive sdigits?
+ Sdigits stored in tree S with the branch of TClonesArray, all sdigits in a single TClonesArray
+ So the tree has only one entry.
+ One needs to say:
+ pRich->GetLoader()->LoadSDigits(); this one open file, get the tree and invoke AliRICH::SetTreeAddress()
AliRun::SetDebug()
+How to get info for tid number?
+ Header and Kinematics trees must be loaded, then possible to retrive pointer to Stack of particles
+ Int_t AliRunLoader::LoadHeader(); Int_t AliRunLoader::LoadKinematics()
+ AliStack *AliRunLoader::Stack()
+ TParticle *AliStack::Particle(tid)
+ TParticle::Print()
+How to deal with AliRunDigitizer?
+ AliRunDigitizer::Exec() just call AliRunDigitizer::Digitize()
-pointer to gAlice
+How to avoid using the pointer to gAlice?
-//__________________________________________________________________________________________________
-void H_SD()
-{
- Info("H_SD","Start.");
-
- for(Int_t iEventN=0;iEventN<a->GetEventsPerRun();iEventN++){//events loop
- al->GetEvent(iEventN);
-
- if(!rl->TreeH()) rl->LoadHits(); al->LoadHeader(); al->LoadKinematics();//from
- if(!rl->TreeS()) rl->MakeTree("S"); r->MakeBranch("S");//to
-
- for(Int_t iPrimN=0;iPrimN<rl->TreeH()->GetEntries();iPrimN++){//prims loop
- rl->TreeH()->GetEntry(iPrimN);
- for(Int_t iHitN=0;iHitN<r->Hits()->GetEntries();iHitN++){//hits loop ???
- AliRICHhit *pHit=r->Hits()->At(iHitN);
-
- TVector2 x2 = r->Param()->ShiftToWirePos(r->C(pHit->C())->Glob2Loc(pHit->OutX3()));
-
- Int_t iTotQdc=r->Param()->TotQdc(x2,pHit->Eloss());
-
- Int_t iPadXmin,iPadXmax,iPadYmin,iPadYmax;
- r->Param()->Loc2Area(x2,iPadXmin,iPadYmin,iPadXmax,iPadYmax);
- cout<<"left-down=("<<iPadXmin<<","<<iPadYmin<<") right-up=("<<iPadXmax<<','<<iPadYmax<<')'<<endl;
- for(Int_t iPadY=iPadYmin;iPadY<=iPadYmax;iPadY++)
- for(Int_t iPadX=iPadXmin;iPadX<=iPadXmax;iPadX++){
- Double_t padQdc=iTotQdc*r->Param()->FracQdc(x2,iPadX,iPadY);
- cout<<padQdc<<endl;
- if(padQdc>0.1) r->AddSDigit(pHit->C(),iPadX,iPadY,padQdc,al->Stack()->Particle(pHit->GetTrack())->GetPdgCode(),pHit->GetTrack());
- }
- }//hits loop
- }//prims loop
- rl->TreeS()->Fill();
- rl->WriteSDigits("OVERWRITE");
- }//events loop
-
- rl->UnloadHits(); al->UnloadHeader(); al->UnloadKinematics();
- rl->UnloadSDigits();
- Info("H_SD","Stop.");
-}//H_SD()
-//__________________________________________________________________________________________________
Int_t countContrib[7][3];
Info("sd","totally %i digits",iTotalDigits);
rl->UnloadDigits();
}
-
+//__________________________________________________________________________________________________
void sc()
{
if(rl->LoadRecPoints()) return;
for(int i=1;i<=7;i++) r->Clusters(i)->Print();
rl->UnloadRecPoints();
}
+//__________________________________________________________________________________________________
+void sp(int tid)
+{
+ al->LoadHeader(); al->LoadKinematics();
+ PrintParticleInfo(tid);
+ al->UnloadKinematics(); al->UnloadHeader();
+}
+//__________________________________________________________________________________________________
+void PrintParticleInfo(int tid)
+{
+ TParticle *p=al->Stack()->Particle(tid);
+ cout<<p->GetName();
+ if(p->GetMother(0)!=-1){cout<<" from "; PrintParticleInfo(p->GetMother(0));}
+ else {cout<<endl;}
+}
+//__________________________________________________________________________________________________
Double_t r2d = TMath::RadToDeg();
Double_t d2r = TMath::DegToRad();
cout << "Info in Digits->Clusters: Time used: ";sw.Print();
}
//__________________________________________________________________________________________________
-void OLD_S_SD()
-{
- Info("OLD_S_SD","Start.");
- rl->LoadHits();
- for(Int_t iEventN=0;iEventN<a->GetEventsPerRun();iEventN++){//events loop
- al->GetEvent(iEventN); Info("OLD_S_SD","Processing event %i",iEventN);
-
- rl->MakeTree("S"); r->MakeBranch("S");
- r->ResetSDigits(); r->ResetSpecialsOld();
- for(Int_t iPrimN=0;iPrimN<rl->TreeH()->GetEntries();iPrimN++){//prims loop
- rl->TreeH()->GetEntry(iPrimN);
- for(Int_t i=0;i<r->Specials()->GetEntries();i++){//specials loop
- Int_t padx=1+ ((AliRICHSDigit*)r->Specials()->At(i))->PadX()+r->Param()->NpadsX()/2;
- Int_t pady=1+ ((AliRICHSDigit*)r->Specials()->At(i))->PadY()+r->Param()->NpadsY()/2;
- Double_t q= ((AliRICHSDigit*)r->Specials()->At(i))->QPad();
-
- Int_t hitN= ((AliRICHSDigit*)r->Specials()->At(i))->HitNumber()-1;//!!! important -1
- Int_t chamber=((AliRICHhit*)r->Hits()->At(hitN))->C();
- Int_t tid=((AliRICHhit*)r->Hits()->At(hitN))->GetTrack();
- Int_t pid=((AliRICHhit*)r->Hits()->At(hitN))->Pid();
- if(padx<1 || padx>r->Param()->NpadsX() ||pady<1 || pady>r->Param()->NpadsY())
- Warning("OLD_S_SD","pad is out of valid range padx= %i pady=%i event %i",padx,pady,iEventN);
- else
- r->AddSDigit(chamber,padx,pady,q,pid,tid);
- }//specials loop
- }//prims loop
- rl->TreeS()->Fill();
- rl->WriteSDigits("OVERWRITE");
- }//events loop
- rl->UnloadHits(); rl->UnloadSDigits();
- Info("OLD_S_SD","Stop.");
-}//OLD_S_SD()
-//__________________________________________________________________________________________________
+AliRICH * Rich() {return r;}
void SD_D()
{
Info("SD_D","Start.");
extern Int_t kBad;
- r->Param()->GenSigmaThMap();
+ Rich()->Param()->GenSigmaThMap();
rl->LoadSDigits();
for(Int_t iEventN=0;iEventN<a->GetEventsPerRun();iEventN++){//events loop
al->GetEvent(iEventN); cout<<"Event "<<iEventN<<endl;
- rl->MakeTree("D");r->MakeBranch("D"); //create TreeD with RICH branches
- r->ResetSDigits();r->ResetDigits();//reset lists of sdigits and digits
- rl->TreeS()->GetEntry(0);
- r->SDigits()->Sort();
+ rl->MakeTree("D"); Rich()->MakeBranch("D"); //create TreeD with RICH branches
+ Rich()->ResetSDigits(); Rich()->ResetDigits(); //reset lists of sdigits and digits
+ rl->TreeS()->GetEntry(0); //get sdigits to memory container
+ Rich()->SDigits()->Sort();
- Int_t combiPid,chamber,x,y,tid[3],id; Double_t q;
+ Int_t combiPid=0,chamber=0,x=0,y=0,tid[3],id=0; Double_t q=0;
Int_t iNdigitsPerPad;//how many sdigits for a given pad
const int kBad=-101;//??? to be removed in code
- for(Int_t i=0;i<r->SDigits()->GetEntries();i++){//sdigits loop (sorted)
- AliRICHdigit *pSdig=(AliRICHdigit*)r->SDigits()->At(i);
+ for(Int_t i=0;i<Rich()->SDigits()->GetEntries();i++){//sdigits loop (sorted)
+ AliRICHdigit *pSdig=(AliRICHdigit*)Rich()->SDigits()->At(i);
if(pSdig->Id()==id){//still the same pad
- iNdigitsPerPad++;
- q+=pSdig->Q();
- combiPid+=pSdig->CombiPid();
+ iNdigitsPerPad++; q+=pSdig->Q(); combiPid+=pSdig->CombiPid();
if(iNdigitsPerPad<=3)
tid[iNdigitsPerPad-1]=pSdig->Tid(0);
else
Warning("SDigits2Digits","More then 3 sdigits for the given pad");
}else{//new pad, add the pevious one
- if(id!=kBad&&r->Param()->IsOverTh(chamber,x,y,q)) {
- r->AddDigit(chamber,x,y,q,combiPid,tid);
+ if(id!=kBad&&Rich()->Param()->IsOverTh(chamber,x,y,q)) {
+ Rich()->AddDigit(chamber,x,y,q,combiPid,tid);
}
combiPid=pSdig->CombiPid();chamber=pSdig->C();id=pSdig->Id();
x=pSdig->X();y=pSdig->Y();
}
}//sdigits loop (sorted)
- if(r->SDigits()->GetEntries()&&r->Param()->IsOverTh(chamber,x,y,q))
- r->AddDigit(chamber,x,y,q,combiPid,tid);//add the last digit
+ if(Rich()->SDigits()->GetEntries() && Rich()->Param()->IsOverTh(chamber,x,y,q))
+ Rich()->AddDigit(chamber,x,y,q,combiPid,tid);//add the last digit
rl->TreeD()->Fill();
rl->WriteDigits("OVERWRITE");
}//events loop
rl->UnloadSDigits(); rl->UnloadDigits();
- r->ResetSDigits();r->ResetDigits();//reset lists of sdigits and digits
+ Rich()->ResetSDigits(); Rich()->ResetDigits();//reset lists of sdigits and digits
Info("SD_D","Stop.");
}//SD_D()
//__________________________________________________________________________________________________
iTotalHits+=r->Hits()->GetEntries();
TParticle *pPrim=al->Stack()->Particle(iPrimN);
Info("Show","Evt %4i prim %4i has %4i hits from %s (,%7.2f,%7.2f)",
- iEventN,
- iPrimN,
- r->Hits()->GetEntries(),
- pPrim->GetName(),
- pPrim->Theta()*r2d,pPrim->Phi()*r2d);
+ iEventN,iPrimN, r->Hits()->GetEntries(), pPrim->GetName(), pPrim->Theta()*r2d,pPrim->Phi()*r2d);
}//prims loop
Info("Show-HITS","Evt %i total: %i particles %i primaries %i hits",
iEventN, iNparticles, iNprims, iTotalHits);
for(int i=1;i<=7;i++)
Info("Show-DIGITS","Evt %i chamber %i contains %5i digits",
iEventN, i, r->Digits(i)->GetEntries());
- }
+ }else
+ Info("Show-DIGITS","There is no digits for this event");
if(isClusters){
rl->TreeR()->GetEntry(0);
for(int i=1;i<=7;i++)
if(ReadAlice()){//it's from file, reconstruct
pMenu->AddButton("hits->sdigits->digits->clusters","MainTrank()","Convert");
- pMenu->AddButton("hits->sdigits" ,"H_SD()" ,"AliRICH::Hits2SDigits");
pMenu->AddButton("sdigits->digits" ,"SD_D()" ,"AliRICHDigitizer");
pMenu->AddButton("digits->clusters" ,"D_C()" ,"AliRICHClusterFinder");
pMenu->AddButton("clusters->recos" ,"C_R()" ,"AliRICHRecon");