///////////////////////////////////////////////////////////////////////////////
-//
-// An overview of the basic philosophy of the ITS code development
-// and analysis is show in the figure below.
-//Begin_Html
-/*
+// //
+// An overview of the basic philosophy of the ITS code development //
+// and analysis is show in the figure below. //
+//Begin_Html //
+/*
<img src="picts/ITS/ITS_Analysis_schema.gif">
</pre>
<br clear=left>
#include "AliITSClusterFinderSDD.h"
#include "AliITSClusterFinderSPD.h"
#include "AliITSClusterFinderSSD.h"
+#include "AliITSClusterFinderV2SDD.h"
+#include "AliITSClusterFinderV2SPD.h"
+#include "AliITSClusterFinderV2SSD.h"
#include "AliITSDetType.h"
#include "AliITSLoader.h"
#include "AliITSRawClusterSPD.h"
#include "AliITSRawClusterSDD.h"
#include "AliITSRawClusterSSD.h"
#include "AliITSRecPoint.h"
+#include "AliITSclusterV2.h"
#include "AliITSdigitSPD.h"
#include "AliITSdigitSDD.h"
#include "AliITSdigitSSD.h"
#include "AliITSDigitizer.h"
#include "AliITSDDLRawData.h"
#include "AliRun.h"
+#include "AliRawReader.h"
ClassImp(AliITS)
fNctype(0),
fRecPoints(0),
fNRecPoints(0),
+ fClustersV2(0),
+ fNClustersV2(0),
fSelectedVertexer(0){
// Default initializer for ITS
// The default constructor of the AliITS class. In addition to
fNctype(0),
fRecPoints(0),
fNRecPoints(0),
+ fClustersV2(0),
+ fNClustersV2(0),
fSelectedVertexer(0){
// The standard Constructor for the ITS class. In addition to
// creating the AliITS class, it allocates memory for the TClonesArrays
fRecPoints = new TClonesArray("AliITSRecPoint",1000);
fNRecPoints = 0;
+ fClustersV2 = new TClonesArray("AliITSclusterV2",3000);
+ fNClustersV2=0;
Int_t i;
for(i=0;i<fNDetTypes;i++) {
delete fRecPoints;
fRecPoints=0;
}
+ if(fClustersV2){
+ fClustersV2->Delete();
+ delete fClustersV2;
+ fClustersV2=0;
+ }
delete[] fIdName; // Array of TStrings
delete[] fIdSens;
if(fITSmodules) {
} // end if
} // end if iDetType
}
+
+
+//______________________________________________________________________
+void AliITS::SetDefaultClusterFindersV2(){
+ // Sets the default cluster finders. Used in finding RecPoints.
+ // Inputs:
+ // none.
+ // Outputs:
+ // none.
+ // Return:
+ // none.
+ AliITSDetType *iDetType;
+ AliITSsegmentation *seg;
+ AliITSClusterFinder *clf;
+
+ MakeTreeC();
+
+ // SPD
+ iDetType=DetType(kSPD);
+ if(iDetType){
+ if (!iDetType->GetReconstructionModel()) {
+ seg =(AliITSsegmentation*)iDetType->GetSegmentationModel();
+ clf = new AliITSClusterFinderV2SPD();
+ clf->SetSegmentation(seg);
+ clf->SetDigits(DigitsAddress(0));
+ SetReconstructionModel(kSPD,clf);
+ } // end if
+ } // end if iDetType
+
+ // SDD
+ iDetType=DetType(kSDD);
+ if(iDetType){
+ if (!iDetType->GetReconstructionModel()) {
+ seg = (AliITSsegmentation*)iDetType->GetSegmentationModel();
+ clf = new AliITSClusterFinderV2SDD();
+ clf->SetSegmentation(seg);
+ clf->SetDigits(DigitsAddress(1));
+ SetReconstructionModel(kSDD,clf);
+ } // end if
+ } // end if iDetType
+
+ // SSD
+ iDetType=DetType(kSSD);
+ if(iDetType){
+ if (!iDetType->GetReconstructionModel()) {
+ seg = (AliITSsegmentation*)iDetType->GetSegmentationModel();
+ clf = new AliITSClusterFinderV2SSD();
+ clf->SetSegmentation(seg);
+ clf->SetDigits(DigitsAddress(2));
+ SetReconstructionModel(kSSD,clf);
+ } // end if
+ } // end if iDetType
+}
+
+
//______________________________________________________________________
void AliITS::MakeBranch(Option_t* option){
// Creates Tree branches for the ITS.
Bool_t cD = (strstr(option,"D")!=0);
Bool_t cR = (strstr(option,"R")!=0);
Bool_t cRF = (strstr(option,"RF")!=0);
+ Bool_t v2 = (strstr(option,"v2")!=0);
+
if(cRF)cR = kFALSE;
if(cH && (fHits == 0x0)) fHits = new TClonesArray("AliITShit", 1560);
if(cD) MakeBranchD(0);
if(cR) MakeBranchR(0);
if(cRF) MakeBranchRF(0);
+ if(v2) MakeBranchR(0,"v2");
}
//______________________________________________________________________
void AliITS::SetTreeAddress(){
// Return:
// none.
+
if(!GetITSgeom()) return; // need transformations to do digitization.
AliITSgeom *geom = GetITSgeom();
// only one branch for rec points for all detector types
Bool_t oFast= (strstr(opt,"Fast")!=0);
+ Bool_t v2 = (strstr(opt,"v2")!=0);
+
+
if(oFast){
sprintf(branchname,"%sRecPointsF",GetName());
+ } else if(v2){
+ sprintf(branchname,"Clusters");
} else {
sprintf(branchname,"%sRecPoints",GetName());
}
- if(!fRecPoints)fRecPoints = new TClonesArray("AliITSRecPoint",1000);
- if (fLoader->TreeR()) {
+ if(v2){
+
+ if(!fClustersV2)fClustersV2 = new TClonesArray("AliITSclusterV2",3000);
+ if(fLoader->TreeR()){
+ if(fClustersV2==0x0) fClustersV2 = new TClonesArray("AliITSclusterV2",3000);
+ MakeBranchInTree(fLoader->TreeR(),branchname,&fClustersV2,buffsz,file);
+
+ }
+ }else{
+ if(!fRecPoints)fRecPoints = new TClonesArray("AliITSRecPoint",1000);
+ if (fLoader->TreeR()) {
if(fRecPoints==0x0) fRecPoints = new TClonesArray("AliITSRecPoint",
1000);
MakeBranchInTree(fLoader->TreeR(),branchname,&fRecPoints,buffsz,file);
- } // end if
+ } // end if
+ }
}
//______________________________________________________________________
void AliITS::SetTreeAddressR(TTree *treeR){
if(!treeR) return;
if(fRecPoints==0x0) fRecPoints = new TClonesArray("AliITSRecPoint",1000);
- TBranch *branch;
- sprintf(branchname,"%sRecPoints",GetName());
- branch = treeR->GetBranch(branchname);
- if (branch) {
+ TBranch *branch1;
+ sprintf(branchname,"Clusters");
+ branch1 = treeR->GetBranch(branchname);
+ if(branch1){
+ if(fClustersV2==0x0) fClustersV2 = new TClonesArray("AliITSclusterV2",3000);
+ branch1->SetAddress(&fClustersV2);
+ }
+ else{
+ TBranch *branch;
+ sprintf(branchname,"%sRecPoints",GetName());
+ branch = treeR->GetBranch(branchname);
+ if (branch) {
branch->SetAddress(&fRecPoints);
- }else {
+ }else {
sprintf(branchname,"%sRecPointsF",GetName());
branch = treeR->GetBranch(branchname);
if (branch) {
- branch->SetAddress(&fRecPoints);
+ branch->SetAddress(&fRecPoints);
}
+ }
}
}
//______________________________________________________________________
new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
}
//______________________________________________________________________
+void AliITS::AddClusterV2(const AliITSclusterV2 &r){
+ // Add a reconstructed space point to the list
+ // Inputs:
+ // const AliITSClusterV2 &r class to be added to the tree
+ // of reconstructed points TreeR.
+ // Outputs:
+ // none.
+ // Return:
+ // none.
+
+ TClonesArray &lrecp = *fClustersV2;
+ new(lrecp[fNClustersV2++]) AliITSclusterV2(r);
+}
+//______________________________________________________________________
void AliITS::HitsToFastRecPoints(Int_t evNumber,Int_t bgrev,Int_t size,
Option_t *opt0,Option_t *opt1,
const char *flnm){
rec->SetClusters(ClustersAddress(id));
rec->FindRawClusters(module);
} // end if
- pITSloader->TreeR()->Fill();
+ pITSloader->TreeR()->Fill();
ResetRecPoints();
+ ResetClustersV2();
treeC->Fill();
ResetClusters();
} // end for module
pITSloader->WriteRecPoints("OVERWRITE");
pITSloader->WriteRawClusters("OVERWRITE");
}
+//______________________________________________________________________
+void AliITS::DigitsToRecPoints(AliRawReader* rawReader){
+ // cluster finding and reconstruction of space points
+ // the condition below will disappear when the geom class will be
+ // initialized for all versions - for the moment it is only for v5 !
+ // 7 is the SDD beam test version
+ // Inputs:
+ // Int_t evNumber Event number to be processed.
+ // Int_t lastentry Offset for module when not all of the modules
+ // are processed.
+ // Option_t *opt String indicating which ITS sub-detectors should
+ // be processed. If ="All" then all of the ITS
+ // sub detectors are processed.
+ // Outputs:
+ // none.
+ // Return:
+ // none.
+ if(!GetITSgeom()) return;
+ AliITSgeom *geom = GetITSgeom();
+
+
+ SetDefaultClusterFindersV2();
+
+ AliITSLoader *pITSloader = (AliITSLoader*)fLoader;
+ AliITSClusterFinderV2 *rec = 0;
+ AliITSDetType *iDetType = 0;
+ Int_t id=0;
+
+ if(!pITSloader->TreeR()) pITSloader->MakeTree("R");
+ TTree* cTree = pITSloader->TreeR();
+ TClonesArray *array=new TClonesArray("AliITSclusterV2",1000);
+ cTree->Branch("Clusters",&array);
+ delete array;
+
+ TClonesArray** clusters = new TClonesArray*[geom->GetIndexMax()];
+ for (Int_t iModule = 0; iModule < geom->GetIndexMax(); iModule++) {
+ clusters[iModule] = NULL;
+ }
+ for(id=0;id<3;id++){
+ iDetType = DetType(id);
+ rec = (AliITSClusterFinderV2*)iDetType->GetReconstructionModel();
+ if (!rec) {
+ Error("DigitsToRecPoints",
+ "The reconstruction class was not instanciated");
+ exit(1);
+ }
+ rec->RawdataToClusters(rawReader,clusters);
+ }
+ Int_t nClusters =0;
+ for(Int_t iModule=0;iModule<geom->GetIndexMax();iModule++){
+ array = clusters[iModule];
+ if(!array){
+ Error("DigitsToRecPoints","data for module %d missing!",iModule);
+ array = new TClonesArray("AliITSclusterV2");
+ }
+ cTree->SetBranchAddress("Clusters",&array);
+ cTree->Fill();
+ nClusters+=array->GetEntriesFast();
+ delete array;
+ }
+ pITSloader->WriteRecPoints("OVERWRITE");
+
+ delete[] clusters;
+ Info("DigitsToRecPoints", "total number of found clustersV2 in ITS: %d\n",
+ nClusters);
+
+}
+
//______________________________________________________________________
AliLoader* AliITS::MakeLoader(const char* topfoldername){
//builds ITSgetter (AliLoader type)
class AliITSsimulation;
class AliITSClusterFinder;
+class AliITSclusterV2;
class AliITSLoader;
class AliITSsegmentation;
class AliITSresponse;
class AliVertexer;
class AliDigitizer;
class AliRunDigitizer;
+class AliRawReader;
const Int_t kNTYPES=3;
virtual void SetDefaults();
virtual void SetDefaultSimulation();
virtual void SetDefaultClusterFinders();
+ virtual void SetDefaultClusterFindersV2();
virtual void MakeBranch(Option_t *opt=" ");
virtual void SetTreeAddress();
#ifndef NEWVERSION
void ResetRecPoints(){if(fRecPoints) fRecPoints->Clear();fNRecPoints = 0;};
// Return pointer to rec points
TClonesArray *RecPoints() {return fRecPoints;}
+
+ void AddClusterV2(const AliITSclusterV2 &cl);
+ void ResetClustersV2(){if(fClustersV2) fClustersV2->Clear();fNClustersV2=0;}
+ Int_t GetNClustersV2()const {return fNClustersV2;}
+// Return pointer to clustersV2
+TClonesArray *ClustersV2() {return fClustersV2;}
#endif
#ifdef NEWVERSION
void MakeBranchR(const char *file, Option_t *opt=" ");
void Digits2Reco(){
DigitsToRecPoints(fLoader->GetRunLoader()->GetEventNumber(),0,fOpt);};
void DigitsToRecPoints(Int_t evNumber,Int_t lastEntry,Option_t *det);
-
+void DigitsToRecPoints(AliRawReader* rawReader);
protected:
//================== Data Members ==================================
#ifdef NEWVERSION
TClonesArray *fRecPoints; //! List of reconstructed points
Int_t fNRecPoints; // Number of rec points
+
+ TClonesArray *fClustersV2; //!List of reconstructed clusters v2
+ Int_t fNClustersV2; //Number of clusters v2
+
#endif
TString fSelectedVertexer; // Vertexer selected in CreateVertexer
#ifndef NEWVERSION
- ClassDef(AliITS,4) // Base class for ITS
+ ClassDef(AliITS,5) // Base class for ITS
#endif
#ifdef NEWVERSION
ClassDef(AliITS,5) // Base class for ITS
* about the suitability of this software for any purpose. It is *
* provided "as is" without express or implied warranty. *
**************************************************************************/
-//
-// Base Class used to find
-// the reconstructed points for ITS
-// See also AliITSClusterFinderSPD, AliITSClusterFinderSDD,
-// AliITSClusterFinderSDD
-//
+////////////////////////////////////////////////////////////////////////////
+// //
+// Base Class used to find //
+// the reconstructed points for ITS //
+// See also AliITSClusterFinderSPD, AliITSClusterFinderSDD, //
+// AliITSClusterFinderSDD AliITSClusterFinderV2 //
+////////////////////////////////////////////////////////////////////////////
#include "AliITSClusterFinder.h"
+#include "AliITSclusterV2.h"
+#include "AliITSRecPoint.h"
#include "AliITSdigitSPD.h"
#include "AliITSdigitSDD.h"
#include "AliITSdigitSSD.h"
+#include "AliITSgeom.h"
#include "AliITSMap.h"
#include "AliRun.h"
+#include "AliMC.h"
#include "AliITS.h"
+#include "TParticle.h"
ClassImp(AliITSClusterFinder)
// none.
// Return:
// A default constructed AliITSCulsterFinder
+ Init();
}
//----------------------------------------------------------------------
AliITSClusterFinder::AliITSClusterFinder(AliITSsegmentation *seg,
SetNperMax();
SetClusterSize();
SetDeclusterFlag();
+ Init();
}
//----------------------------------------------------------------------
AliITSClusterFinder::AliITSClusterFinder(AliITSsegmentation *seg,
SetNperMax();
SetClusterSize();
SetDeclusterFlag();
+ Init();
}
+//______________________________________________________________________
+AliITSClusterFinder::AliITSClusterFinder(const AliITSClusterFinder &source) : TObject(source) {
+ // Copy constructor
+ // Copies are not allowed. The method is protected to avoid misuse.
+ Fatal("AliITSClusterFinder","Copy constructor not allowed\n");
+}
+
+//______________________________________________________________________
+AliITSClusterFinder& AliITSClusterFinder::operator=(const AliITSClusterFinder& /* source */){
+ // Assignment operator
+ // Assignment is not allowed. The method is protected to avoid misuse.
+ Fatal("= operator","Assignment operator not allowed\n");
+ return *this;
+}
+
//----------------------------------------------------------------------
AliITSClusterFinder::~AliITSClusterFinder(){
// destructor cluster finder
fDeclusterFlag= 0;
fClusterSize = 0;
fNPeaks = 0;
+ fITS = 0;
}
+
//__________________________________________________________________________
-AliITSClusterFinder::AliITSClusterFinder(const AliITSClusterFinder &source) :
- TObject(source){
- // Copy Constructor
- // Input:
- // AliITSClusterFinder &source The class which will become a copy of
- // this class
- // Output:
- // none.
- // Return:
- // A copy of this class
+void AliITSClusterFinder::Init(){
+
+ //Initialisation of ITS geometry
+
+ fITS = (AliITS*)gAlice->GetModule("ITS");
+ AliITSgeom *geom=(AliITSgeom*)fITS->GetITSgeom();
+
+ Int_t mmax=geom->GetIndexMax();
+ if (mmax>2200) {
+ Fatal("AliITSClusterFinder","Too many ITS subdetectors !");
+ }
+ Int_t m;
+ for (m=0; m<mmax; m++) {
+ Int_t lay,lad,det; geom->GetModuleId(m,lay,lad,det);
+ Float_t x,y,z; geom->GetTrans(lay,lad,det,x,y,z);
+ Double_t rot[9]; geom->GetRotMatrix(lay,lad,det,rot);
+ Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
+ Double_t ca=TMath::Cos(alpha), sa=TMath::Sin(alpha);
+ fYshift[m] = x*ca + y*sa;
+ fZshift[m] = (Double_t)z;
+ fNdet[m] = (lad-1)*geom->GetNdetectors(lay) + (det-1);
+ fNlayer[m] = lay-1;
+ }
- if(&source == this) return;
- *this = source;
- return;
-}
-//______________________________________________________________________
-AliITSClusterFinder& AliITSClusterFinder::operator=(const AliITSClusterFinder
- &source) {
- // Assignment operator
- // Input:
- // AliITSClusterFinder &source The class which will become a copy of
- // this class
- // Output:
- // none.
- // Return:
- // A copy of this class
-
- if(&source == this) return *this;
- this->fDigits = source.fDigits;
- this->fNdigits = source.fNdigits;
- this->fResponse = source.fResponse;
- this->fSegmentation = source.fSegmentation;
- this->fNRawClusters = source.fNRawClusters;
- this->fMap = source.fMap;
- this->fNperMax = source.fNperMax;
- this->fDeclusterFlag = source.fDeclusterFlag;
- this->fClusterSize = source.fClusterSize;
- this->fNPeaks = source.fNPeaks;
- return *this;
}
+
+
//----------------------------------------------------------------------
void AliITSClusterFinder::AddCluster(Int_t branch, AliITSRawCluster *c){
// Add a raw cluster copy to the list
// Return:
// none.
- AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
- iTS->AddCluster(branch,c);
+ fITS->AddCluster(branch,c);
fNRawClusters++;
}
//----------------------------------------------------------------------
// Return:
// none.
- AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
- iTS->AddCluster(branch,c);
+ fITS->AddCluster(branch,c);
fNRawClusters++;
- iTS->AddRecPoint(rp);
+ fITS->AddRecPoint(rp);
}
+
+//______________________________________________________________________
+void AliITSClusterFinder::CheckLabels(Int_t lab[3]) {
+ //------------------------------------------------------------
+ // Tries to find mother's labels
+ //------------------------------------------------------------
+
+ if(lab[0]<0 && lab[1]<0 && lab[2]<0) return; // In case of no labels just exit
+ // Check if simulation
+ AliMC* mc = gAlice->GetMCApp();
+ if(!mc)return;
+
+ Int_t ntracks = mc->GetNtrack();
+ for (Int_t i=0;i<3;i++){
+ Int_t label = lab[i];
+ if (label>=0 && label<ntracks) {
+ TParticle *part=(TParticle*)mc->Particle(label);
+ if (part->P() < 0.005) {
+ Int_t m=part->GetFirstMother();
+ if (m<0) {
+ continue;
+ }
+ if (part->GetStatusCode()>0) {
+ continue;
+ }
+ lab[i]=m;
+ }
+ }
+ }
+
+}
+
//______________________________________________________________________
void AliITSClusterFinder::FindRawClusters(Int_t module){
// Default Cluster finder.
*os << fClusterSize<<",";
*os << fNPeaks<<endl;
}
+//______________________________________________________________________
+void AliITSClusterFinder::RecPoints2Clusters
+(const TClonesArray *points, Int_t idx, TClonesArray *clusters) {
+ //------------------------------------------------------------
+ // Conversion AliITSRecPoint -> AliITSclusterV2 for the ITS
+ // subdetector indexed by idx
+ //------------------------------------------------------------
+ AliITSgeom* geom = (AliITSgeom*)fITS->GetITSgeom();
+ Int_t lastSPD1=geom->GetModuleIndex(2,1,1)-1;
+ TClonesArray &cl=*clusters;
+ Int_t ncl=points->GetEntriesFast();
+ for (Int_t i=0; i<ncl; i++) {
+ AliITSRecPoint *p = (AliITSRecPoint *)points->UncheckedAt(i);
+ Float_t lp[5];
+ lp[0]=-(-p->GetX()+fYshift[idx]); if (idx<=lastSPD1) lp[0]*=-1; //SPD1
+ lp[1]= -p->GetZ()+fZshift[idx];
+ lp[2]=p->GetSigmaX2();
+ lp[3]=p->GetSigmaZ2();
+ lp[4]=p->GetQ()*36./23333.; //electrons -> ADC
+ Int_t lab[4];
+ lab[0]=p->GetLabel(0); lab[1]=p->GetLabel(1); lab[2]=p->GetLabel(2);
+ lab[3]=fNdet[idx];
+ CheckLabels(lab);
+ Int_t dummy[3]={0,0,0};
+ new (cl[i]) AliITSclusterV2(lab,lp, dummy);
+ }
+}
+
//______________________________________________________________________
void AliITSClusterFinder::Read(istream *is){
//Standard input for this class
source.Read(&is);
return is;
}
+/*
+void AliITSClusterFinder::RecPoints2Clusters
+(const TClonesArray *points, Int_t idx, TClonesArray *clusters) {
+ //------------------------------------------------------------
+ // Conversion AliITSRecPoint -> AliITSclusterV2 for the ITS
+ // subdetector indexed by idx
+ //------------------------------------------------------------
+ TClonesArray &cl=*clusters;
+ Int_t ncl=points->GetEntriesFast();
+ for (Int_t i=0; i<ncl; i++) {
+ AliITSRecPoint *p = (AliITSRecPoint *)points->UncheckedAt(i);
+ Float_t lp[5];
+ lp[0]=-(-p->GetX()+fYshift[idx]); if (idx<=fLastSPD1) lp[0]*=-1; //SPD1
+ lp[1]= -p->GetZ()+fZshift[idx];
+ lp[2]=p->GetSigmaX2();
+ lp[3]=p->GetSigmaZ2();
+ lp[4]=p->GetQ()*36./23333.; //electrons -> ADC
+ Int_t lab[4];
+ lab[0]=p->GetLabel(0); lab[1]=p->GetLabel(1); lab[2]=p->GetLabel(2);
+ lab[3]=fNdet[idx];
+ CheckLabels(lab);
+ Int_t dummy[3]={0,0,0};
+ new (cl[i]) AliITSclusterV2(lab,lp, dummy);
+ }
+}
+*/
////////////////////////////////////////////////
// ITS Cluster Finder Class //
+// //
+// //
////////////////////////////////////////////////
-#include <Riostream.h>
#include <TObject.h>
#include <TClonesArray.h>
AliITSClusterFinder(AliITSsegmentation *seg, AliITSresponse *resp,
TClonesArray *digits);// Standard+ Constructor
virtual ~AliITSClusterFinder(); // Destructor
- AliITSClusterFinder(const AliITSClusterFinder &source); // copy constructor
- // assignment operator
- AliITSClusterFinder& operator=(const AliITSClusterFinder &source);
//
// Do the Reconstruction.
virtual void FindRawClusters(Int_t mod=0); // Finds cluster of digits.
- //
+ //
// Sets the debug flag for debugging output
void SetDebug(Int_t level=1){fDebug=level;}
// Clears the debug flag so no debugging output will be generated
void Print(ostream *os); // Class ascii print function
void Read(istream *os); // Class ascii read function
- // data members
+ // Conversion from RecPoints to V2Clusters
+ void RecPoints2Clusters(const TClonesArray *points, Int_t idx, TClonesArray *clusters);
+
+
protected:
+ // methods
+ static void CheckLabels(Int_t lab[3]);
+ void Init();
+ AliITSClusterFinder(const AliITSClusterFinder &source); // copy constructor
+ // assignment operator
+ AliITSClusterFinder& operator=(const AliITSClusterFinder &source);
+ // data members
Int_t fDebug; //! Debug flag/level
Int_t fModule; //! Module number to be reconstuctted
TClonesArray *fDigits; //! digits
// and want to keep track of
// the cluster which was splitted
AliITSMap *fMap; //! map
+ AliITS *fITS; //! pointer to the ITS
Int_t fNperMax; //! NperMax
Int_t fDeclusterFlag; //! DeclusterFlag
Int_t fClusterSize; //! ClusterSize
Int_t fNPeaks; //! NPeaks
+ // Data members needed to fill AliCluster objects
+ Float_t fYshift[2200]; // y-shifts of detector local coor. systems
+ Float_t fZshift[2200]; // z-shifts of detector local coor. systems
+ Int_t fNdet[2200]; // detector index
+ Int_t fNlayer[2200]; // detector layer
- ClassDef(AliITSClusterFinder,3) //Class for clustering and reconstruction of space points
+ ClassDef(AliITSClusterFinder,4) //Class for clustering and reconstruction of space points
};
// Input and output functions for standard C++ input/output.
ostream &operator<<(ostream &os,AliITSClusterFinder &source);
/*
$Id$
$Log$
+ Revision 1.37 2004/06/10 21:00:24 nilsen
+ Modifications associated with remerging the Ba/Sa and Dubna pixel simulations,
+ some cleaning of general code (including coding convensions), and adding some
+ protections associated with SetDefaults/SetDefaultSimulations which should help
+ with the Test beam simulations. Details below. The default SPD simulation for
+ the general ITS runs/geometry is still the Ba/Sa, but for the Test beam
+ geometries this has been changed to the merged versions.
+ File: AliITS.cxx Modified
+ File: AliITS.h Modified
+ In lined many one-two line functions. Added some protection to
+ SetDefaults(), SetDefaultSimulation(), and SetDefaultClusterFinders(),
+ such that they should now even work when only one detector type has
+ been defined (as it should be for the test beams...). Some mostly
+ cosmetic issues associated with getting branch names for digits. And
+ Generally some cleaning up of the code.
+ File: AliITSClusterFinder.cxx Modified
+ File: AliITSClusterFinder.h Modified
+ Did some additional consolidation of data into the base class, added
+ TClonesArray *fClusters, a fDebug, and fModule variables. Otherwise
+ some cosmetic and coding conversion changes.
+ File: AliITSClusterFinderSDD.cxx Modified
+ File: AliITSClusterFinderSDD.h Modified
+ Changes to be consistent with the modified base class, and cosmetic
+ and coding conversion changes.
+ File: AliITSClusterFinderSPD.cxx Modified
+ File: AliITSClusterFinderSPD.h Modified
+ Changes to be consistent with the modified base class, and cosmetic
+ and coding conversion changes.
+ File: AliITSClusterFinderSPDdubna.h Removed
+ File: AliITSClusterFinderSPDdubna.cxx Removed
+ Since we have ClusterFinderSPD and V2 and this version isn't being
+ maintained, it is being retired.
+ File: AliITSClusterFinderSSD.cxx Modified
+ File: AliITSClusterFinderSSD.h Modified
+ Changes to be consistent with the modified base class, and cosmetic
+ and coding conversion changes.
+ File: AliITSDetType.cxx Modified
+ File: AliITSDetType.h Modified
+ Added a new class variable to indicate what the detector type is
+ AliITSDetector fDetType; values of kSPD, kSDD, kSSD, .... Otherwise
+ cosmetic and Coding convention changes.
+ File: AliITSLoader.cxx Modified
+ File: AliITSLoader.h Modified
+ Some changes which are not complete. The idea is to be able to get,
+ simply via one call, a specific hit, Sdigit, digit, RecPoint,...
+ without all of the usual over head of initializing TClonesArrays setting
+ branch addresses and the like. Work is far form ready.
+ File: AliITSdcsSSD.cxx Modified
+ Some nearly cosmetic changes necessary due to changes to response and
+ segmentation class'.
+ File: AliITSgeom.h Modified
+ In the definition of AliITSDetector type, added kND=-1, no detector
+ defined. Expect to use it later(?).
+ File: AliITSresponse.h Modified
+ Basically cosmetic. Mostly changing Float_t to Double_t.
+ File: AliITSresponseSDD.cxx Modified
+ File: AliITSresponseSDD.h Modified
+ Basically the cosmetic and Float_t to Double_t
+ File: AliITSresponseSPD.cxx Modified
+ File: AliITSresponseSPD.h Modified
+ Mostly Float_t to Double_t and added in the IsPixelDead function for
+ the dubna version (otherwise the merging had been done).
+ File: AliITSresponseSPDdubna.h Removed
+ File: AliITSresponseSPDdubna.cxx Removed
+ We should be able to remove this class now. AliITSresponseSPD is now
+ used for both the Bari-Salerno and the dubna models.
+ File: AliITSresponseSSD.cxx Modified
+ File: AliITSresponseSSD.h Modified
+ Float_t to Double_t changes.
+ File: AliITSsegmentation.h Modified
+ Made LocaltoDet return a Bool_t. Now if the x,z location is outside
+ of the volume, it returns kFALSE. see below.
+ File: AliITSsegmentationSDD.cxx Modified
+ File: AliITSsegmentationSDD.h Modified
+ Made LocaltoDet return a Bool_t. Now if the x,z location is outside
+ of the volume, it returns kFALSE.
+ File: AliITSsegmentationSPD.cxx Modified
+ File: AliITSsegmentationSPD.h Modified
+ Made LocaltoDet return a Bool_t. Now if the x,z location is outside
+ of the volume, it returns kFALSE.
+ File: AliITSsegmentationSSD.cxx Modified
+ File: AliITSsegmentationSSD.h Modified
+ Made LocaltoDet return a Bool_t. Now if the x,z location is outside
+ of the volume, it returns kFALSE. see below.
+ File: AliITSsimulation.cxx Modified
+ File: AliITSsimulation.h Modified
+ Added fDebug variable, new Constructor for use below. Cosmetic and
+ coding convention changes
+ File: AliITSsimulationSDD.cxx Modified
+ File: AliITSsimulationSDD.h Modified
+ Added new Constructor, removed redundant variables and Cosmetic and
+ coding convention changes.
+ File: AliITSsimulationSPD.cxx Modified
+ File: AliITSsimulationSPD.h Modified
+ Removed some dead code, made changes as needed by the changes above
+ (response and segmentation classes...). a few cosmetic and coding
+ convention changes.
+ File: AliITSsimulationSPDdubna.cxx Modified
+ File: AliITSsimulationSPDdubna.h Modified
+ New merged version, implemented new and old coupling with switch,
+ coding convention and similar changes. (found 1 bugs, missing
+ ! in front of if(mod-LineSegmentL(....,).
+ File: AliITSsimulationSSD.cxx Modified
+ File: AliITSsimulationSSD.h Modified
+ removed redundant variables with base class. Fixed for coding
+ convention and other cosmetic changes.
+ File: AliITSvSDD03.cxx Modified
+ File: AliITSvSPD02.cxx Modified
+ File: AliITSvSSD03.cxx Modified
+ These two have their private versions of SetDefaults and
+ SetDefaultSimulation which have been similarly protected as in AliITS.cxx
+ File: ITSLinkDef.h Modified
+ File: libITS.pkg Modified
+ Versions which include v11 geometry and other private changes
+
Revision 1.36 2004/01/27 16:12:03 masera
Coding conventions for AliITSdigitXXX classes and AliITSTrackerV1
noise values).
*/
-//
-// Cluster finder
-// for Silicon
-// Drift Detector
-//
+///////////////////////////////////////////////////////////////////////////
+// Cluster finder //
+// for Silicon //
+// Drift Detector //
+//////////////////////////////////////////////////////////////////////////
#include <Riostream.h>
#include <TMath.h>
//______________________________________________________________________
void AliITSClusterFinderSDD::Find1DClusters(){
// find 1D clusters
- static AliITS *iTS = (AliITS*)gAlice->GetModule("ITS");
// retrieve the parameters
Int_t fNofMaps = GetSeg()->Npz();
clusteranodePath, //f
clusterMult, //i
0,0,0,0,0,0,0);//7*i
- iTS->AddCluster(1,&clust);
+ fITS->AddCluster(1,&clust);
it = tstop;
} // ilcl
it++;
//______________________________________________________________________
void AliITSClusterFinderSDD::Find1DClustersE(){
// find 1D clusters
- static AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
// retrieve the parameters
Int_t fNofMaps = GetSeg()->Npz();
Int_t fMaxNofSamples = GetSeg()->Npx();
driftPath,anodePath,
nTsteps,start,stop,
start, stop, 1, k, k );
- iTS->AddCluster( 1, &clust );
+ fITS->AddCluster( 1, &clust );
if(GetDebug(5)) clust.PrintInfo();
nClu++;
} // end if nTsteps
}
//______________________________________________________________________
-void AliITSClusterFinderSDD::ResolveClustersE(){
+void AliITSClusterFinderSDD::ResolveClusters(){
// The function to resolve clusters if the clusters overlapping exists
Int_t i;
- static AliITS *iTS = (AliITS*)gAlice->GetModule( "ITS" );
// get number of clusters for this module
Int_t nofClusters = NClusters();
nofClusters -= fNclusters;
Int_t xdim = tstop-tstart+3;
Int_t zdim = astop-astart+3;
if( xdim > 50 || zdim > 30 ) {
- Warning("ResolveClustersE","xdim: %d , zdim: %d ",xdim,zdim);
+ Warning("ResolveClusters","xdim: %d , zdim: %d ",xdim,zdim);
continue;
}
Double_t *sp = new Double_t[ xdim*zdim+1 ];
}
if( peakpos < 0 ) {
- //Warning("ResolveClustersE",
+ //Warning("ResolveClusters",
// "Digit not found for cluster");
//if(GetDebug(3)) clusterI.PrintInfo();
continue;
clusterI.SetTsigma( tau[i]*fTimeStep );
clusterI.SetQ( integral[i] );
- iTS->AddCluster( 1, &clusterI );
+ fITS->AddCluster( 1, &clusterI );
} // end for i
Clusters()->RemoveAt( j );
delete [] par;
} else { // something odd
- Warning( "ResolveClustersE",
+ Warning( "ResolveClusters",
"--- Peak not found!!!! minpeak=%d ,cluster peak= %f"
" , module= %d",
fMinPeak, clusterJ->PeakAmpl(),GetModule());
clusterJ->PrintInfo();
- Warning( "ResolveClustersE"," xdim= %d zdim= %d", xdim-2, zdim-2 );
+ Warning( "ResolveClusters"," xdim= %d zdim= %d", xdim-2, zdim-2 );
}
delete [] sp;
} // cluster loop
Clusters()->Compress();
return;
}
-//__________________________________________________________________________
-void AliITSClusterFinderSDD::ResolveClusters(){
- // The function to resolve clusters if the clusters overlapping exists
-/* AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
- // get number of clusters for this module
- Int_t nofClusters = NClusters();
- nofClusters -= fNclusters;
- //cout<<"Resolve Cl: nofClusters, fNclusters ="<<nofClusters<<","
- // <<fNclusters<<endl;
- Int_t fNofMaps = GetSeg()->Npz();
- Int_t fNofAnodes = fNofMaps/2;
- Int_t dummy=0;
- Double_t fTimeStep = GetSeg()->Dpx(dummy);
- Double_t fSddLength = GetSeg()->Dx();
- Double_t fDriftSpeed = GetResp()->DriftSpeed();
- Double_t anodePitch = GetSeg()->Dpz(dummy);
- Double_t n, baseline;
- GetResp()->GetNoiseParam(n,baseline);
- Double_t dzz_1A = anodePitch * anodePitch / 12;
- // fill Map of signals
- Map()->FillMap();
- Int_t j,i,ii,ianode,anode,itime;
- Int_t wing,astart,astop,tstart,tstop,nanode;
- Double_t fadc,ClusterTime;
- Double_t q[400],x[400],z[400]; // digit charges and coordinates
- for(j=0; j<nofClusters; j++) {
- AliITSRawClusterSDD *clusterJ=(AliITSRawClusterSDD*) Cluster(j);
- Int_t ndigits = 0;
- astart=clusterJ->Astart();
- astop=clusterJ->Astop();
- tstart=clusterJ->Tstartf();
- tstop=clusterJ->Tstopf();
- nanode=clusterJ->Anodes(); // <- Ernesto
- wing=(Int_t)clusterJ->W();
- if(wing == 2) {
- astart += fNofAnodes;
- astop += fNofAnodes;
- } // end if
- // cout<<"astart,astop,tstart,tstop ="<<astart<<","<<astop<<","
- // <<tstart<<","<<tstop<<endl;
- // clear the digit arrays
- for(ii=0; ii<400; ii++) {
- q[ii] = 0.;
- x[ii] = 0.;
- z[ii] = 0.;
- } // end for ii
- for(ianode=astart; ianode<=astop; ianode++) {
- for(itime=tstart; itime<=tstop; itime++) {
- fadc=Map()->GetSignal(ianode,itime);
- if(fadc>baseline) {
- fadc-=(Double_t)baseline;
- q[ndigits] = fadc*(fTimeStep/160); // KeV
- anode = ianode;
- if(wing == 2) anode -= fNofAnodes;
- z[ndigits] = (anode + 0.5 - fNofAnodes/2)*anodePitch;
- ClusterTime = itime*fTimeStep;
- if(ClusterTime > fTimeCorr) ClusterTime -= fTimeCorr;// ns
- x[ndigits] = fSddLength - ClusterTime*fDriftSpeed;
- if(wing == 1) x[ndigits] *= (-1);
- // cout<<"ianode,itime,fadc ="<<ianode<<","<<itime<<","
- // <<fadc<<endl;
- // cout<<"wing,anode,ndigits,charge ="<<wing<<","
- // <<anode<<","<<ndigits<<","<<q[ndigits]<<endl;
- ndigits++;
- continue;
- } // end if
- fadc=0;
- // cout<<"fadc=0, ndigits ="<<ndigits<<endl;
- } // time loop
- } // anode loop
- // cout<<"for new cluster ndigits ="<<ndigits<<endl;
- // Fit cluster to resolve for two separate ones --------------------
- Double_t qq=0., xm=0., zm=0., xx=0., zz=0., xz=0.;
- Double_t dxx=0., dzz=0., dxz=0.;
- Double_t scl = 0., tmp, tga, elps = -1.;
- Double_t xfit[2], zfit[2], qfit[2];
- Double_t pitchz = anodePitch*1.e-4; // cm
- Double_t pitchx = fTimeStep*fDriftSpeed*1.e-4; // cm
- Double_t sigma2;
- Int_t nfhits;
- Int_t nbins = ndigits;
- Int_t separate = 0;
- // now, all lengths are in microns
- for (ii=0; ii<nbins; ii++) {
- qq += q[ii];
- xm += x[ii]*q[ii];
- zm += z[ii]*q[ii];
- xx += x[ii]*x[ii]*q[ii];
- zz += z[ii]*z[ii]*q[ii];
- xz += x[ii]*z[ii]*q[ii];
- } // end for ii
- xm /= qq;
- zm /= qq;
- xx /= qq;
- zz /= qq;
- xz /= qq;
- dxx = xx - xm*xm;
- dzz = zz - zm*zm;
- dxz = xz - xm*zm;
-
- // shrink the cluster in the time direction proportionaly to the
- // dxx/dzz, which lineary depends from the drift path
- // new Ernesto........
- if( nanode == 1 ){
- dzz = dzz_1A; // for one anode cluster dzz = anode**2/12
- scl = TMath::Sqrt( 7.2/(-0.57*xm*1.e-3+71.8) );
- } // end if
- if( nanode == 2 ){
- scl = TMath::Sqrt( (-0.18*xm*1.e-3+21.3)/(-0.57*xm*1.e-3+71.8) );
- } // end if
- if( nanode == 3 ){
- scl = TMath::Sqrt( (-0.5*xm*1.e-3+34.5)/(-0.57*xm*1.e-3+71.8) );
- } // end if
- if( nanode > 3 ){
- scl = TMath::Sqrt( (1.3*xm*1.e-3+49.)/(-0.57*xm*1.e-3+71.8) );
- } // end if
- // cout<<"1 microns: zm,dzz,xm,dxx,dxz,qq ="<<zm<<","<<dzz<<","
- // <<xm<<","<<dxx<<","<<dxz<<","<<qq<<endl;
- // old Boris.........
- // tmp=29730. - 585.*fabs(xm/1000.);
- // scl=TMath::Sqrt(tmp/130000.);
-
- xm *= scl;
- xx *= scl*scl;
- xz *= scl;
-
- dxx = xx - xm*xm;
- // dzz = zz - zm*zm;
- dxz = xz - xm*zm;
- // cout<<"microns: zm,dzz,xm,dxx,xz,dxz,qq ="<<zm<<","<<dzz<<","
- // <<xm<<","<<dxx<<","<<xz<<","<<dxz<<","<<qq<<endl;
- // if(dzz < 7200.) dzz=7200.;//for one anode cluster dzz = anode**2/12
-
- if (dxx < 0.) dxx=0.;
- // the data if no cluster overlapping (the coordunates are in cm)
- nfhits = 1;
- xfit[0] = xm*1.e-4;
- zfit[0] = zm*1.e-4;
- qfit[0] = qq;
- // if(nbins < 7) cout<<"**** nbins ="<<nbins<<endl;
-
- if (nbins >= 7) {
- if (dxz==0.) tga=0.;
- else {
- tmp=0.5*(dzz-dxx)/dxz;
- tga = (dxz<0.) ? tmp-TMath::Sqrt(tmp*tmp+1) :
- tmp+TMath::Sqrt(tmp*tmp+1);
- } // end if dxz
- elps=(tga*tga*dxx-2*tga*dxz+dzz)/(dxx+2*tga*dxz+tga*tga*dzz);
- // change from microns to cm
- xm *= 1.e-4;
- zm *= 1.e-4;
- zz *= 1.e-8;
- xx *= 1.e-8;
- xz *= 1.e-8;
- dxz *= 1.e-8;
- dxx *= 1.e-8;
- dzz *= 1.e-8;
- // cout<<"cm: zm,dzz,xm,dxx,xz,dxz,qq ="<<zm<<","<<dzz<<","
- // <<xm<<","<<dxx<<","<<xz<<","<<dxz<<","<<qq<<endl;
- for (i=0; i<nbins; i++) {
- x[i] = x[i] *= scl;
- x[i] = x[i] *= 1.e-4;
- z[i] = z[i] *= 1.e-4;
- } // end for i
- // cout<<"!!! elps ="<<elps<<endl;
- if (elps < 0.3) { // try to separate hits
- separate = 1;
- tmp=atan(tga);
- Double_t cosa=cos(tmp),sina=sin(tmp);
- Double_t a1=0., x1=0., xxx=0.;
- for (i=0; i<nbins; i++) {
- tmp=x[i]*cosa + z[i]*sina;
- if (q[i] > a1) {
- a1=q[i];
- x1=tmp;
- } // end if
- xxx += tmp*tmp*tmp*q[i];
- } // end for i
- xxx /= qq;
- Double_t z12=-sina*xm + cosa*zm;
- sigma2=(sina*sina*xx-2*cosa*sina*xz+cosa*cosa*zz) - z12*z12;
- xm=cosa*xm + sina*zm;
- xx=cosa*cosa*xx + 2*cosa*sina*xz + sina*sina*zz;
- Double_t x2=(xx - xm*x1 - sigma2)/(xm - x1);
- Double_t r=a1*2*TMath::ACos(-1.)*sigma2/(qq*pitchx*pitchz);
- for (i=0; i<33; i++) { // solve a system of equations
- Double_t x1_old=x1, x2_old=x2, r_old=r;
- Double_t c11=x1-x2;
- Double_t c12=r;
- Double_t c13=1-r;
- Double_t c21=x1*x1 - x2*x2;
- Double_t c22=2*r*x1;
- Double_t c23=2*(1-r)*x2;
- Double_t c31=3*sigma2*(x1-x2) + x1*x1*x1 - x2*x2*x2;
- Double_t c32=3*r*(sigma2 + x1*x1);
- Double_t c33=3*(1-r)*(sigma2 + x2*x2);
- Double_t f1=-(r*x1 + (1-r)*x2 - xm);
- Double_t f2=-(r*(sigma2+x1*x1)+(1-r)*(sigma2+x2*x2)- xx);
- Double_t f3=-(r*x1*(3*sigma2+x1*x1)+(1-r)*x2*
- (3*sigma2+x2*x2)-xxx);
- Double_t d=c11*c22*c33+c21*c32*c13+c12*c23*c31-
- c31*c22*c13 - c21*c12*c33 - c32*c23*c11;
- if (d==0.) {
- cout<<"*********** d=0 ***********\n";
- break;
- } // end if
- Double_t dr=f1*c22*c33 + f2*c32*c13 + c12*c23*f3 -
- f3*c22*c13 - f2*c12*c33 - c32*c23*f1;
- Double_t d1=c11*f2*c33 + c21*f3*c13 + f1*c23*c31 -
- c31*f2*c13 - c21*f1*c33 - f3*c23*c11;
- Double_t d2=c11*c22*f3 + c21*c32*f1 + c12*f2*c31 -
- c31*c22*f1 - c21*c12*f3 - c32*f2*c11;
- r += dr/d;
- x1 += d1/d;
- x2 += d2/d;
- if (fabs(x1-x1_old) > 0.0001) continue;
- if (fabs(x2-x2_old) > 0.0001) continue;
- if (fabs(r-r_old)/5 > 0.001) continue;
- a1=r*qq*pitchx*pitchz/(2*TMath::ACos(-1.)*sigma2);
- Double_t a2=a1*(1-r)/r;
- qfit[0]=a1; xfit[0]=x1*cosa - z12*sina; zfit[0]=x1*sina +
- z12*cosa;
- qfit[1]=a2; xfit[1]=x2*cosa - z12*sina; zfit[1]=x2*sina +
- z12*cosa;
- nfhits=2;
- break; // Ok !
- } // end for i
- if (i==33) cerr<<"No more iterations ! "<<endl;
- } // end of attempt to separate overlapped clusters
- } // end of nbins cut
- if(elps < 0.) cout<<" elps=-1 ="<<elps<<endl;
- if(elps >0. && elps< 0.3 && nfhits == 1) cout<<" small elps, nfh=1 ="
- <<elps<<","<<nfhits<<endl;
- if(nfhits == 2) cout<<" nfhits=2 ="<<nfhits<<endl;
- for (i=0; i<nfhits; i++) {
- xfit[i] *= (1.e+4/scl);
- if(wing == 1) xfit[i] *= (-1);
- zfit[i] *= 1.e+4;
- // cout<<" --------- i,xfiti,zfiti,qfiti ="<<i<<","
- // <<xfit[i]<<","<<zfit[i]<<","<<qfit[i]<<endl;
- } // end for i
- Int_t ncl = nfhits;
- if(nfhits == 1 && separate == 1) {
- cout<<"!!!!! no separate"<<endl;
- ncl = -2;
- } // end if
- if(nfhits == 2) {
- cout << "Split cluster: " << endl;
- clusterJ->PrintInfo();
- cout << " in: " << endl;
- for (i=0; i<nfhits; i++) {
- // AliITSRawClusterSDD *clust = new AliITSRawClusterSDD(wing,
- -1,-1,(Double_t)qfit[i],ncl,0,0,
- (Double_t)xfit[i],
- (Double_t)zfit[i],0,0,0,0,
- tstart,tstop,astart,astop);
- // AliITSRawClusterSDD *clust = new AliITSRawClusterSDD(wing,-1,
- // -1,(Double_t)qfit[i],0,0,0,
- // (Double_t)xfit[i],
- // (Double_t)zfit[i],0,0,0,0,
- // tstart,tstop,astart,astop,ncl);
- // ???????????
- // if(wing == 1) xfit[i] *= (-1);
- Double_t Anode = (zfit[i]/anodePitch+fNofAnodes/2-0.5);
- Double_t Time = (fSddLength - xfit[i])/fDriftSpeed;
- Double_t clusterPeakAmplitude = clusterJ->PeakAmpl();
- Double_t peakpos = clusterJ->PeakPos();
- Double_t clusteranodePath = (Anode - fNofAnodes/2)*anodePitch;
- Double_t clusterDriftPath = Time*fDriftSpeed;
- clusterDriftPath = fSddLength-clusterDriftPath;
- AliITSRawClusterSDD *clust = new AliITSRawClusterSDD(wing,Anode,
- Time,qfit[i],
- clusterPeakAmplitude,peakpos,
- 0.,0.,clusterDriftPath,
- clusteranodePath,clusterJ->Samples()/2
- ,tstart,tstop,0,0,0,astart,astop);
- clust->PrintInfo();
- iTS->AddCluster(1,clust);
- // cout<<"new cluster added: tstart,tstop,astart,astop,x,ncl ="
- // <<tstart<<","<<tstop<<","<<astart<<","<<astop<<","<<xfit[i]
- // <<","<<ncl<<endl;
- delete clust;
- }// nfhits loop
- Clusters()->RemoveAt(j);
- } // if nfhits = 2
-} // cluster loop
-Clusters()->Compress();
-Map()->ClearMap();
-*/
- return;
-}
//______________________________________________________________________
void AliITSClusterFinderSDD::GetRecPoints(){
// get rec points
- static AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
+
// get number of clusters for this module
Int_t nofClusters = NClusters();
nofClusters -= fNclusters;
if(dig) rnew.fTracks[1]=dig->GetTrack(1);
if(dig) rnew.fTracks[2]=dig->GetTrack(2);
- iTS->AddRecPoint(rnew);
+ fITS->AddRecPoint(rnew);
} // I clusters
// Map()->ClearMap();
}
Find1DClustersE();
GroupClusters();
SelectClusters();
- ResolveClustersE();
+ ResolveClusters();
GetRecPoints();
}
//_______________________________________________________________________
void GroupClusters();
void SelectClusters();
void GetRecPoints();
- void ResolveClusters(); // Boris........
- void ResolveClustersE(); // Ernesto
+ void ResolveClusters(); // Ernesto Lopez Torres
Int_t SearchPeak(Double_t *spect,Int_t xdim,Int_t zdim,Int_t *peakX,
Int_t *peakZ,Double_t *peakAmp,Double_t minpeak);//Ernesto
Int_t NoLinearFit( Int_t xdim, Int_t zdim, Double_t *param, Double_t *spe,
* about the suitability of this software for any purpose. It is *
* provided "as is" without express or implied warranty. *
**************************************************************************/
-//
-// Cluster finder
-// for Silicon pixels
-//
-//
+////////////////////////////////////////////////////////////////////////////
+// Cluster finder ///
+// for Silicon pixels //
+// //
+////////////////////////////////////////////////////////////////////////////
#include "AliITSClusterFinderSPD.h"
#include "AliITS.h"
#include "AliITSdigitSPD.h"
SetDx();
SetDz();
}
-//_____________________________________________________________________
-AliITSClusterFinderSPD::AliITSClusterFinderSPD(
- const AliITSClusterFinderSPD &source): AliITSClusterFinder(source){
- // Copy Constructor
- if(&source == this) return;
- this->fDz = source.fDz;
- this->fDx = source.fDx;
- this->fMinNCells = source.fMinNCells;
- return;
-}
//______________________________________________________________________
-AliITSClusterFinderSPD& AliITSClusterFinderSPD::operator=(
- const AliITSClusterFinderSPD &source) {
- // Assignment operator
+AliITSClusterFinderSPD::AliITSClusterFinderSPD(const AliITSClusterFinderSPD &source) : AliITSClusterFinder(source) {
+ // Copy constructor
+ // Copies are not allowed. The method is protected to avoid misuse.
+ Fatal("AliITSClusterFinderSPD","Copy constructor not allowed\n");
+}
- if(&source == this) return *this;
- this->fDz = source.fDz;
- this->fDx = source.fDx;
- this->fMinNCells = source.fMinNCells;
- return *this;
+//______________________________________________________________________
+AliITSClusterFinderSPD& AliITSClusterFinderSPD::operator=(const AliITSClusterFinderSPD& /* source */){
+ // Assignment operator
+ // Assignment is not allowed. The method is protected to avoid misuse.
+ Fatal("= operator","Assignment operator not allowed\n");
+ return *this;
}
//______________________________________________________________________
void AliITSClusterFinderSPD::FindRawClusters(Int_t module){
} // end if iclus[i]
// store the cluster information to the AliITSRawCLusterSPD object
- static AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
+
//put the cluster center in local reference frame of the detector
// and in microns
(Double_t) ndzmin,
(Double_t) ndzmax,
0,GetModule());
- iTS->AddCluster(0,clust);
+ fITS->AddCluster(0,clust);
delete clust;
}//end loop on clusters
delete[] ifpad;
const Double_t kconv = 1.0e-4; // micron -> cm
// get rec points
- static AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
for (Int_t i=0; i<nclus; i++){
l[0] = kconv*xcenter[i];
l[1] = kconv*GetSeg()->Dy()/2.;
rnew.fTracks[0]=tr1clus[i];
rnew.fTracks[1]=tr2clus[i];
rnew.fTracks[2]=tr3clus[i];
- iTS->AddRecPoint(rnew);
+ fITS->AddRecPoint(rnew);
} // end for i
}
AliITSClusterFinderSPD(AliITSsegmentation *segmentation,
TClonesArray *digits,TClonesArray *recpoints);
virtual ~AliITSClusterFinderSPD(){}// destructor
- // copy constructor
- AliITSClusterFinderSPD(const AliITSClusterFinderSPD &source);
- // assignment operator
- AliITSClusterFinderSPD& operator=(const AliITSClusterFinderSPD &source);
virtual AliITSresponseSPD* GetResp()const{
return (AliITSresponseSPD*) GetResponse();}//Return Response
Double_t xcenter[],Double_t zcenter[],
Double_t errxcenter[],Double_t errzcenter[],
Int_t tr1clus[],Int_t tr2clus[], Int_t tr3clus[]);
- private:
+ protected:
+ // copy constructor
+ AliITSClusterFinderSPD(const AliITSClusterFinderSPD &source);
+ // assignment operator
+ AliITSClusterFinderSPD& operator=(const AliITSClusterFinderSPD &source);
+
Double_t fDz; // dz
Double_t fDx; // dx
Int_t fMinNCells; // min num of cells in the cluster
* about the suitability of this software for any purpose. It is *
* provided "as is" without express or implied warranty. *
**************************************************************************/
-
-// **************************************************************************
-// * The package was revised and changed by Boris Batiounia in the time *
-// * period of March - June 2001 *
-// **************************************************************************/
-//
+///////////////////////////////////////////////////////////////////////////////
+// **************************************************************************//
+// * The package was revised and changed by Boris Batiounia in the time //
+// * period of March - June 2001 //
+// **************************************************************************//
+///////////////////////////////////////////////////////////////////////////////
#include <Riostream.h>
#include <TArrayI.h>
#include "AliRun.h"
//______________________________________________________________________
AliITSClusterFinderSSD::AliITSClusterFinderSSD():
AliITSClusterFinder(),
-fITS(0),
fClusterP(0),
fNClusterP(0),
fClusterN(0),
AliITSClusterFinderSSD::AliITSClusterFinderSSD(AliITSsegmentation *seg,
TClonesArray *digits):
AliITSClusterFinder(seg,0),
-fITS(0),
fClusterP(0),
fNClusterP(0),
fClusterN(0),
SetDigits(digits);
SetMap(new AliITSMapA1(GetSeg(),Digits()));
- fITS = (AliITS*)gAlice->GetModule("ITS");
fClusterP = new TClonesArray ("AliITSclusterSSD",200);
fNClusterP = 0;
fClusterN = new TClonesArray ("AliITSclusterSSD",200);
AliITSClusterFinderSSD::AliITSClusterFinderSSD(AliITSsegmentation *seg,
AliITSresponse *res):
AliITSClusterFinder(seg,res),
-fITS(0),
fClusterP(0),
fNClusterP(0),
fClusterN(0),
fSFB(0){
//Standard constructor
- fITS = (AliITS*)gAlice->GetModule("ITS");
fClusterP = new TClonesArray ("AliITSclusterSSD",200);
fNClusterP = 0;
fClusterN = new TClonesArray ("AliITSclusterSSD",200);
fPitch = GetSeg()->Dpx(0);
fPNsignalRatio= 7./8.; // warning: hard-wired number
}
+
+//______________________________________________________________________
+AliITSClusterFinderSSD::AliITSClusterFinderSSD(const AliITSClusterFinderSSD &source) : AliITSClusterFinder(source) {
+ // Copy constructor
+ // Copies are not allowed. The method is protected to avoid misuse.
+ Fatal("AliITSClusterFinderSSD","Copy constructor not allowed\n");
+}
+
+//______________________________________________________________________
+AliITSClusterFinderSSD& AliITSClusterFinderSSD::operator=(const AliITSClusterFinderSSD& /* source */){
+ // Assignment operator
+ // Assignment is not allowed. The method is protected to avoid misuse.
+ Fatal("= operator","Assignment operator not allowed\n");
+ return *this;
+}
+
//______________________________________________________________________
AliITSClusterFinderSSD::~AliITSClusterFinderSSD(){
// Default destructor
- fITS = 0;
delete fClusterP;
delete fClusterN;
delete fPackages;
#define ALIITSCLUSTERFINDERSSD_H
/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
* See cxx source for full Copyright notice */
-
+/////////////////////////////////////////////////////////////////////////////
+// SSD Cluster Finder //
+// //
+/////////////////////////////////////////////////////////////////////////////
//#include <TMath.h>
#include "AliITSClusterFinder.h"
//#include "AliITSsegmentationSSD.h"
void FindRawClusters(Int_t module);
protected:
+ // copy constructor
+ AliITSClusterFinderSSD(const AliITSClusterFinderSSD &source);
+ // assignment operator
+ AliITSClusterFinderSSD& operator=(const AliITSClusterFinderSSD &source);
virtual AliITSresponseSSD* GetResp()const{
return (AliITSresponseSSD*) GetResponse();}//Return Response
//Returns fSegmentation
void GetCrossingError(Double_t& dp, Double_t& dn);
// Data memebers
- AliITS *fITS; //!Pointer to AliITS object
TClonesArray *fClusterP; //!
Int_t fNClusterP; //!Number of P side clusters in the array
TClonesArray *fClusterN; //!Number of N side clusters in the array
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+////////////////////////////////////////////////////////////////////////////
+// Implementation of the ITS clusterer V2 class //
+// //
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
+// //
+///////////////////////////////////////////////////////////////////////////
+
+
+#include "AliRun.h"
+#include "AliITSClusterFinderV2.h"
+#include "AliITSRecPoint.h"
+#include "AliITSclusterV2.h"
+#include "AliITSsimulationFastPoints.h"
+#include "AliITS.h"
+#include "AliITSgeom.h"
+#include <TParticle.h>
+#include "AliMC.h"
+
+ClassImp(AliITSClusterFinderV2)
+
+extern AliRun *gAlice;
+
+AliITSClusterFinderV2::AliITSClusterFinderV2():AliITSClusterFinder(){
+
+ //Default constructor
+ fEvent = 0;
+ fModule = 0;
+
+ AliITSgeom *geom=(AliITSgeom*)fITS->GetITSgeom();
+
+ fNModules = geom->GetIndexMax();
+}
+
+//______________________________________________________________________
+AliITSClusterFinderV2::AliITSClusterFinderV2(const AliITSClusterFinderV2 &source) : AliITSClusterFinder(source) {
+ // Copy constructor
+ // Copies are not allowed. The method is protected to avoid misuse.
+ Fatal("AliITSClusterFinderV2","Copy constructor not allowed\n");
+}
+
+//______________________________________________________________________
+AliITSClusterFinderV2& AliITSClusterFinderV2::operator=(const AliITSClusterFinderV2& /* source */){
+ // Assignment operator
+ // Assignment is not allowed. The method is protected to avoid misuse.
+ Fatal("= operator","Assignment operator not allowed\n");
+ return *this;
+}
+
+
+//______________________________________________________________________
+void AliITSClusterFinderV2::CheckLabels2(Int_t lab[10]) {
+ //------------------------------------------------------------
+ // Tries to find mother's labels
+ //------------------------------------------------------------
+ Int_t nlabels =0;
+ for (Int_t i=0;i<10;i++) if (lab[i]>=0) nlabels++;
+ if(nlabels == 0) return; // In case of no labels just exit
+
+
+ Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
+
+ for (Int_t i=0;i<10;i++){
+ Int_t label = lab[i];
+ if (label>=0 && label<ntracks) {
+ TParticle *part=(TParticle*)gAlice->GetMCApp()->Particle(label);
+
+ if (part->P() < 0.02) {
+ Int_t m=part->GetFirstMother();
+ if (m<0) {
+ continue;
+ }
+ if (part->GetStatusCode()>0) {
+ continue;
+ }
+ lab[i]=m;
+ }
+ else
+ if (part->P() < 0.12 && nlabels>3) {
+ lab[i]=-2;
+ nlabels--;
+ }
+ }
+ else{
+ if ( (label>ntracks||label <0) && nlabels>3) {
+ lab[i]=-2;
+ nlabels--;
+ }
+ }
+ }
+ if (nlabels>3){
+ for (Int_t i=0;i<10;i++){
+ if (nlabels>3){
+ Int_t label = lab[i];
+ if (label>=0 && label<ntracks) {
+ TParticle *part=(TParticle*)gAlice->GetMCApp()->Particle(label);
+ if (part->P() < 0.1) {
+ lab[i]=-2;
+ nlabels--;
+ }
+ }
+ }
+ }
+ }
+
+ //compress labels -- if multi-times the same
+ Int_t lab2[10];
+ for (Int_t i=0;i<10;i++) lab2[i]=-2;
+ for (Int_t i=0;i<10 ;i++){
+ if (lab[i]<0) continue;
+ for (Int_t j=0;j<10 &&lab2[j]!=lab[i];j++){
+ if (lab2[j]<0) {
+ lab2[j]= lab[i];
+ break;
+ }
+ }
+ }
+ for (Int_t j=0;j<10;j++) lab[j]=lab2[j];
+
+}
+
+//______________________________________________________________________
+void AliITSClusterFinderV2::AddLabel(Int_t lab[10], Int_t label) {
+
+ //add label to the cluster
+
+ if(label<0) return; // In case of no label just exit
+
+ Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
+ if (label>ntracks) return;
+ for (Int_t i=0;i<10;i++){
+ // if (label<0) break;
+ if (lab[i]==label) break;
+ if (lab[i]<0) {
+ lab[i]= label;
+ break;
+ }
+ }
+}
+
+
+//______________________________________________________________________
+void AliITSClusterFinderV2::
+FindCluster(Int_t k,Int_t maxz,AliBin *bins,Int_t &n,Int_t *idx) {
+ //------------------------------------------------------------
+ // returns an array of indices of digits belonging to the cluster
+ // (needed when the segmentation is not regular)
+ //------------------------------------------------------------
+ if (n<200) idx[n++]=bins[k].GetIndex();
+ bins[k].Use();
+
+ if (bins[k-maxz].IsNotUsed()) FindCluster(k-maxz,maxz,bins,n,idx);
+ if (bins[k-1 ].IsNotUsed()) FindCluster(k-1 ,maxz,bins,n,idx);
+ if (bins[k+maxz].IsNotUsed()) FindCluster(k+maxz,maxz,bins,n,idx);
+ if (bins[k+1 ].IsNotUsed()) FindCluster(k+1 ,maxz,bins,n,idx);
+ /*
+ if (bins[k-maxz-1].IsNotUsed()) FindCluster(k-maxz-1,maxz,bins,n,idx);
+ if (bins[k-maxz+1].IsNotUsed()) FindCluster(k-maxz+1,maxz,bins,n,idx);
+ if (bins[k+maxz-1].IsNotUsed()) FindCluster(k+maxz-1,maxz,bins,n,idx);
+ if (bins[k+maxz+1].IsNotUsed()) FindCluster(k+maxz+1,maxz,bins,n,idx);
+ */
+}
+
+//______________________________________________________________________
+Bool_t AliITSClusterFinderV2::IsMaximum(Int_t k,Int_t max,const AliBin *bins) {
+ //------------------------------------------------------------
+ //is this a local maximum ?
+ //------------------------------------------------------------
+ UShort_t q=bins[k].GetQ();
+ if (q==1023) return kFALSE;
+ if (bins[k-max].GetQ() > q) return kFALSE;
+ if (bins[k-1 ].GetQ() > q) return kFALSE;
+ if (bins[k+max].GetQ() > q) return kFALSE;
+ if (bins[k+1 ].GetQ() > q) return kFALSE;
+ if (bins[k-max-1].GetQ() > q) return kFALSE;
+ if (bins[k+max-1].GetQ() > q) return kFALSE;
+ if (bins[k+max+1].GetQ() > q) return kFALSE;
+ if (bins[k-max+1].GetQ() > q) return kFALSE;
+ return kTRUE;
+}
+
+//______________________________________________________________________
+void AliITSClusterFinderV2::
+FindPeaks(Int_t k,Int_t max,AliBin *b,Int_t *idx,UInt_t *msk,Int_t& n) {
+ //------------------------------------------------------------
+ //find local maxima
+ //------------------------------------------------------------
+ if (n<31)
+ if (IsMaximum(k,max,b)) {
+ idx[n]=k; msk[n]=(2<<n);
+ n++;
+ }
+ b[k].SetMask(0);
+ if (b[k-max].GetMask()&1) FindPeaks(k-max,max,b,idx,msk,n);
+ if (b[k-1 ].GetMask()&1) FindPeaks(k-1 ,max,b,idx,msk,n);
+ if (b[k+max].GetMask()&1) FindPeaks(k+max,max,b,idx,msk,n);
+ if (b[k+1 ].GetMask()&1) FindPeaks(k+1 ,max,b,idx,msk,n);
+}
+
+//______________________________________________________________________
+void AliITSClusterFinderV2::
+MarkPeak(Int_t k, Int_t max, AliBin *bins, UInt_t m) {
+ //------------------------------------------------------------
+ //mark this peak
+ //------------------------------------------------------------
+ UShort_t q=bins[k].GetQ();
+
+ bins[k].SetMask(bins[k].GetMask()|m);
+
+ if (bins[k-max].GetQ() <= q)
+ if ((bins[k-max].GetMask()&m) == 0) MarkPeak(k-max,max,bins,m);
+ if (bins[k-1 ].GetQ() <= q)
+ if ((bins[k-1 ].GetMask()&m) == 0) MarkPeak(k-1 ,max,bins,m);
+ if (bins[k+max].GetQ() <= q)
+ if ((bins[k+max].GetMask()&m) == 0) MarkPeak(k+max,max,bins,m);
+ if (bins[k+1 ].GetQ() <= q)
+ if ((bins[k+1 ].GetMask()&m) == 0) MarkPeak(k+1 ,max,bins,m);
+}
+
+//______________________________________________________________________
+void AliITSClusterFinderV2::
+MakeCluster(Int_t k,Int_t max,AliBin *bins,UInt_t m,AliITSclusterV2 &c) {
+ //------------------------------------------------------------
+ //make cluster using digits of this peak
+ //------------------------------------------------------------
+ Float_t q=(Float_t)bins[k].GetQ();
+ Int_t i=k/max, j=k-i*max;
+
+ c.SetQ(c.GetQ()+q);
+ c.SetY(c.GetY()+i*q);
+ c.SetZ(c.GetZ()+j*q);
+ c.SetSigmaY2(c.GetSigmaY2()+i*i*q);
+ c.SetSigmaZ2(c.GetSigmaZ2()+j*j*q);
+
+ bins[k].SetMask(0xFFFFFFFE);
+
+ if (bins[k-max].GetMask() == m) MakeCluster(k-max,max,bins,m,c);
+ if (bins[k-1 ].GetMask() == m) MakeCluster(k-1 ,max,bins,m,c);
+ if (bins[k+max].GetMask() == m) MakeCluster(k+max,max,bins,m,c);
+ if (bins[k+1 ].GetMask() == m) MakeCluster(k+1 ,max,bins,m,c);
+}
+
+//______________________________________________________________________
+Int_t AliITSClusterFinderV2::Hits2Clusters(TTree *hTree, TTree *cTree) {
+ //------------------------------------------------------------
+ // This function creates ITS clusters
+ //------------------------------------------------------------
+
+ AliITSgeom *geom=fITS->GetITSgeom();
+ Int_t mmax=geom->GetIndexMax();
+
+ fITS->InitModules(-1,mmax);
+ fITS->FillModules(hTree,0);
+
+ TClonesArray *clusters=new TClonesArray("AliITSclusterV2",1000);
+ TBranch *branch=cTree->GetBranch("Clusters");
+ if (!branch) cTree->Branch("Clusters",&clusters);
+ else branch->SetAddress(&clusters);
+
+ static TClonesArray *points=fITS->RecPoints();
+ AliITSsimulationFastPoints sim;
+ Int_t ncl=0;
+ for (Int_t m=0; m<mmax; m++) {
+ AliITSmodule *mod=fITS->GetModule(m);
+ sim.CreateFastRecPoints(mod,m,gRandom);
+
+ RecPoints2Clusters(points, m, clusters);
+ fITS->ResetRecPoints();
+
+ ncl+=clusters->GetEntriesFast();
+ cTree->Fill();
+ clusters->Clear();
+ }
+
+ Info("Hits2Clusters","Number of found fast clusters : %d",ncl);
+
+ //cTree->Write();
+
+ delete clusters;
+
+ return 0;
+}
+
--- /dev/null
+#ifndef ALIITSCLUSTERFINDERV2_H
+#define ALIITSCLUSTERFINDERV2_H
+////////////////////////////////////////////////////////////////////
+// ITS clusterer V2 //
+// //
+// //
+// //
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
+////////////////////////////////////////////////////////////////////
+#include "AliITSClusterFinder.h"
+
+
+
+class AliITS;
+class AliITSclusterV2;
+class AliRawReader;
+
+class AliITSClusterFinderV2 : public AliITSClusterFinder {
+public:
+ AliITSClusterFinderV2();
+ virtual ~AliITSClusterFinderV2() {;}
+
+ void SetEvent(Int_t event) { fEvent=event; }
+ virtual void RawdataToClusters(AliRawReader* /*rawReader*/,TClonesArray** /*clusters*/){
+ Warning("RawdataToClusters","Method not implemented in this class ");}
+ virtual Int_t Hits2Clusters(TTree *in, TTree *out);
+
+
+
+protected:
+ class Ali1Dcluster {
+ public:
+ void SetY(Float_t y) {fY=y;}
+ void SetQ(Float_t q) {fQ=q;}
+ void SetNd(Int_t n) {fNd=n;}
+ void SetLabels(Int_t *lab) {fLab[0]=lab[0];fLab[1]=lab[1];fLab[2]=lab[2];}
+ Float_t GetY() const {return fY;}
+ Float_t GetQ() const {return fQ;}
+ Int_t GetNd()const {return fNd;}
+ Int_t GetLabel(Int_t lab) const { return fLab[lab]; }
+ protected:
+ Float_t fY; //cluster position
+ Float_t fQ; //cluster charge
+ Int_t fNd; //number of digits
+ Int_t fLab[3]; //track label
+ };
+ class AliBin {
+ public:
+ AliBin() {fIndex=0; fQ=0; fMask=0xFFFFFFFE;}
+ void SetIndex(UInt_t idx) {fIndex=idx;}
+ void SetQ(UShort_t q) {fQ=q;}
+ void SetMask(UInt_t m) {fMask=m;}
+
+ void Use() {fMask&=0xFFFFFFFE;}
+ Bool_t IsNotUsed() const {return (fMask&1);}
+ Bool_t IsUsed() const {return !(IsNotUsed());}
+
+ UInt_t GetIndex() const {return fIndex;}
+ UShort_t GetQ() const {return fQ;}
+ UInt_t GetMask() const {return fMask;}
+ protected:
+ UInt_t fIndex; //digit index
+ UInt_t fMask; //peak mask
+ UShort_t fQ; //signal
+ };
+ static Bool_t IsMaximum(Int_t k, Int_t max, const AliBin *bins);
+ static void FindPeaks(Int_t k,Int_t m,AliBin*b,Int_t*idx,UInt_t*msk,Int_t&n);
+ static void MarkPeak(Int_t k, Int_t max, AliBin *bins, UInt_t m);
+ static void MakeCluster(Int_t k,Int_t max,AliBin *bins,UInt_t m,
+ AliITSclusterV2 &c);
+
+ static void FindCluster(Int_t k,Int_t maxz,AliBin *bins,Int_t &n,Int_t *idx);
+
+protected:
+ AliITSClusterFinderV2(const AliITSClusterFinderV2 &source); // copy constructor
+ // assignment operator
+ AliITSClusterFinderV2& operator=(const AliITSClusterFinderV2 &source);
+
+ static void CheckLabels2(Int_t lab[10]);
+ static void AddLabel(Int_t lab[10], Int_t label);
+
+ Int_t fNModules; // total number of modules
+ Int_t fEvent; //event number
+
+ ClassDef(AliITSClusterFinderV2,1) // ITS cluster finder V2
+};
+// Input and output functions for standard C++ input/output.
+ostream &operator<<(ostream &os,AliITSClusterFinderV2 &source);
+istream &operator>>(istream &os,AliITSClusterFinderV2 &source);
+#endif
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+////////////////////////////////////////////////////////////////////////////
+// Implementation of the ITS clusterer V2 class //
+// //
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
+// //
+///////////////////////////////////////////////////////////////////////////
+
+
+#include "AliRun.h"
+
+#include "AliITSClusterFinderV2SDD.h"
+#include "AliITSclusterV2.h"
+#include "AliRawReader.h"
+#include "AliITSRawStreamSDD.h"
+
+#include <TClonesArray.h>
+#include "AliITS.h"
+#include "AliITSgeom.h"
+#include "AliITSdigitSDD.h"
+
+ClassImp(AliITSClusterFinderV2SDD)
+
+extern AliRun *gAlice;
+
+AliITSClusterFinderV2SDD::AliITSClusterFinderV2SDD():AliITSClusterFinderV2(){
+
+ //Default constructor
+
+ fNySDD=256; fNzSDD=256;
+ fYpitchSDD=0.01825;
+ fZpitchSDD=0.02940;
+ fHwSDD=3.5085; fHlSDD=3.7632;
+ fYoffSDD=0.0425;
+
+
+
+}
+
+
+void AliITSClusterFinderV2SDD::FindRawClusters(Int_t mod){
+
+ //Find clusters V2
+ SetModule(mod);
+ FindClustersSDD(fDigits);
+
+}
+
+void AliITSClusterFinderV2SDD::FindClustersSDD(TClonesArray *digits) {
+ //------------------------------------------------------------
+ // Actual SDD cluster finder
+ //------------------------------------------------------------
+ Int_t kNzBins = fNzSDD + 2;
+ const Int_t kMAXBIN=kNzBins*(fNySDD+2);
+
+ AliBin *bins[2];
+ bins[0]=new AliBin[kMAXBIN];
+ bins[1]=new AliBin[kMAXBIN];
+
+ AliITSdigitSDD *d=0;
+ Int_t i, ndigits=digits->GetEntriesFast();
+ for (i=0; i<ndigits; i++) {
+ d=(AliITSdigitSDD*)digits->UncheckedAt(i);
+ Int_t y=d->GetCoord2()+1; //y
+ Int_t z=d->GetCoord1()+1; //z
+ Int_t q=d->GetSignal();
+ if (q<3) continue;
+
+ if (z <= fNzSDD) {
+ bins[0][y*kNzBins+z].SetQ(q);
+ bins[0][y*kNzBins+z].SetMask(1);
+ bins[0][y*kNzBins+z].SetIndex(i);
+ } else {
+ z-=fNzSDD;
+ bins[1][y*kNzBins+z].SetQ(q);
+ bins[1][y*kNzBins+z].SetMask(1);
+ bins[1][y*kNzBins+z].SetIndex(i);
+ }
+ }
+
+ FindClustersSDD(bins, kMAXBIN, kNzBins, digits);
+
+ delete[] bins[0];
+ delete[] bins[1];
+
+}
+
+void AliITSClusterFinderV2SDD::
+FindClustersSDD(AliBin* bins[2], Int_t nMaxBin, Int_t nzBins,
+ TClonesArray *digits, TClonesArray *clusters) {
+ //------------------------------------------------------------
+ // Actual SDD cluster finder
+ //------------------------------------------------------------
+ Int_t ncl=0;
+ TClonesArray &cl=*clusters;
+ for (Int_t s=0; s<2; s++)
+ for (Int_t i=0; i<nMaxBin; i++) {
+ if (bins[s][i].IsUsed()) continue;
+ Int_t idx[32]; UInt_t msk[32]; Int_t npeaks=0;
+ FindPeaks(i, nzBins, bins[s], idx, msk, npeaks);
+
+ if (npeaks>30) continue;
+ if (npeaks==0) continue;
+
+ Int_t k,l;
+ for (k=0; k<npeaks-1; k++){//mark adjacent peaks
+ if (idx[k] < 0) continue; //this peak is already removed
+ for (l=k+1; l<npeaks; l++) {
+ if (idx[l] < 0) continue; //this peak is already removed
+ Int_t ki=idx[k]/nzBins, kj=idx[k] - ki*nzBins;
+ Int_t li=idx[l]/nzBins, lj=idx[l] - li*nzBins;
+ Int_t di=TMath::Abs(ki - li);
+ Int_t dj=TMath::Abs(kj - lj);
+ if (di>1 || dj>1) continue;
+ if (bins[s][idx[k]].GetQ() > bins[s][idx[l]].GetQ()) {
+ msk[l]=msk[k];
+ idx[l]*=-1;
+ } else {
+ msk[k]=msk[l];
+ idx[k]*=-1;
+ break;
+ }
+ }
+ }
+
+ for (k=0; k<npeaks; k++) {
+ MarkPeak(TMath::Abs(idx[k]), nzBins, bins[s], msk[k]);
+ }
+
+ for (k=0; k<npeaks; k++) {
+ if (idx[k] < 0) continue; //removed peak
+ AliITSclusterV2 c;
+ MakeCluster(idx[k], nzBins, bins[s], msk[k], c);
+ //mi change
+ Int_t milab[10];
+ for (Int_t ilab=0;ilab<10;ilab++){
+ milab[ilab]=-2;
+ }
+ Int_t maxi=0,mini=0,maxj=0,minj=0;
+ //AliBin *bmax=&bins[s][idx[k]];
+ //Float_t max = TMath::Max(TMath::Abs(bmax->GetQ())/5.,3.);
+ Float_t max=3;
+ for (Int_t di=-2; di<=2;di++)
+ for (Int_t dj=-3;dj<=3;dj++){
+ Int_t index = idx[k]+di+dj*nzBins;
+ if (index<0) continue;
+ if (index>=nMaxBin) continue;
+ AliBin *b=&bins[s][index];
+ if (TMath::Abs(b->GetQ())>max){
+ if (di>maxi) maxi=di;
+ if (di<mini) mini=di;
+ if (dj>maxj) maxj=dj;
+ if (dj<minj) minj=dj;
+ //
+ if(digits) {
+ if (TMath::Abs(di)<2&&TMath::Abs(dj)<2){
+ AliITSdigitSDD* d=(AliITSdigitSDD*)digits->UncheckedAt(b->GetIndex());
+ for (Int_t itrack=0;itrack<10;itrack++){
+ Int_t track = (d->GetTracks())[itrack];
+ if (track>=0) {
+ AddLabel(milab, track);
+ }
+ }
+ }
+ }
+ }
+ }
+
+ /*
+ Float_t s2 = c.GetSigmaY2()/c.GetQ() - c.GetY()*c.GetY();
+ Float_t w=par->GetPadPitchWidth(sec);
+ c.SetSigmaY2(s2);
+ if (s2 != 0.) {
+ c.SetSigmaY2(c.GetSigmaY2()*0.108);
+ if (sec<par->GetNInnerSector()) c.SetSigmaY2(c.GetSigmaY2()*2.07);
+ }
+ s2 = c.GetSigmaZ2()/c.GetQ() - c.GetZ()*c.GetZ();
+ w=par->GetZWidth();
+ c.SetSigmaZ2(s2);
+
+ if (s2 != 0.) {
+ c.SetSigmaZ2(c.GetSigmaZ2()*0.169);
+ if (sec<par->GetNInnerSector()) c.SetSigmaZ2(c.GetSigmaZ2()*1.77);
+ }
+ */
+
+ c.SetSigmaY2(0.0030*0.0030);
+ c.SetSigmaZ2(0.0020*0.0020);
+ c.SetDetectorIndex(fNdet[fModule]);
+
+ Float_t y=c.GetY(),z=c.GetZ(), q=c.GetQ();
+ y/=q; z/=q;
+ //
+ //Float_t s2 = c.GetSigmaY2()/c.GetQ() - y*y;
+ // c.SetSigmaY2(s2);
+ //s2 = c.GetSigmaZ2()/c.GetQ() - z*z;
+ //c.SetSigmaZ2(s2);
+ //
+ y=(y-0.5)*fYpitchSDD;
+ y-=fHwSDD;
+ y-=fYoffSDD; //delay ?
+ if (s) y=-y;
+
+ z=(z-0.5)*fZpitchSDD;
+ z-=fHlSDD;
+
+ y=-(-y+fYshift[fModule]);
+ z= -z+fZshift[fModule];
+ c.SetY(y);
+ c.SetZ(z);
+ c.SetNy(maxj-minj+1);
+ c.SetNz(maxi-mini+1);
+ c.SetType(npeaks);
+ c.SetQ(q/12.7); //to be consistent with the SSD charges
+
+ if (c.GetQ() < 20.) continue; //noise cluster
+
+ if (digits) {
+ // AliBin *b=&bins[s][idx[k]];
+ // AliITSdigitSDD* d=(AliITSdigitSDD*)digits->UncheckedAt(b->GetIndex());
+ {
+ //Int_t lab[3];
+ //lab[0]=(d->GetTracks())[0];
+ //lab[1]=(d->GetTracks())[1];
+ //lab[2]=(d->GetTracks())[2];
+ //CheckLabels(lab);
+ CheckLabels2(milab);
+ c.SetLabel(milab[0],0);
+ c.SetLabel(milab[1],1);
+ c.SetLabel(milab[2],2);
+ c.SetLayer(fNlayer[fModule]);
+ }
+ }
+ if(clusters) new (cl[ncl]) AliITSclusterV2(c);
+ else {
+ fITS->AddClusterV2(c);
+ }
+ ncl++;
+ }
+ }
+}
+
+
+
+void AliITSClusterFinderV2SDD::RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters){
+ //------------------------------------------------------------
+ // This function creates ITS clusters from raw data
+ //------------------------------------------------------------
+ rawReader->Reset();
+ AliITSRawStreamSDD inputSDD(rawReader);
+ FindClustersSDD(&inputSDD,clusters);
+
+}
+
+void AliITSClusterFinderV2SDD::FindClustersSDD(AliITSRawStream* input,
+ TClonesArray** clusters)
+{
+ //------------------------------------------------------------
+ // Actual SDD cluster finder for raw data
+ //------------------------------------------------------------
+ Int_t nClustersSDD = 0;
+ Int_t kNzBins = fNzSDD + 2;
+ Int_t kMaxBin = kNzBins * (fNySDD+2);
+ AliBin* bins[2] = {NULL, NULL};
+
+ // read raw data input stream
+ while (kTRUE) {
+ Bool_t next = input->Next();
+ if (!next || input->IsNewModule()) {
+ Int_t iModule = input->GetPrevModuleID();
+
+ // when all data from a module was read, search for clusters
+ if (bins[0]) {
+ clusters[iModule] = new TClonesArray("AliITSclusterV2");
+ fModule = iModule;
+ FindClustersSDD(bins, kMaxBin, kNzBins, NULL, clusters[iModule]);
+ Int_t nClusters = clusters[iModule]->GetEntriesFast();
+ nClustersSDD += nClusters;
+ delete[] bins[0];
+ delete[] bins[1];
+ }
+
+ if (!next) break;
+ bins[0] = new AliBin[kMaxBin];
+ bins[1] = new AliBin[kMaxBin];
+ }
+
+ // fill the current digit into the bins array
+ if(input->GetSignal()>=3) {
+ Int_t iz = input->GetCoord1()+1;
+ Int_t side = ((iz <= fNzSDD) ? 0 : 1);
+ iz -= side*fNzSDD;
+ Int_t index = (input->GetCoord2()+1) * kNzBins + iz;
+ bins[side][index].SetQ(input->GetSignal());
+ bins[side][index].SetMask(1);
+ bins[side][index].SetIndex(index);
+ }
+ }
+
+ Info("FindClustersSDD", "found clusters in ITS SDD: %d", nClustersSDD);
+
+}
+
+
--- /dev/null
+#ifndef ALIITSCLUSTERFINDERV2SDD_H
+#define ALIITSCLUSTERFINDERV2SDD_H
+//--------------------------------------------------------------
+// ITS clusterer V2 for SDD
+//
+// This can be a "wrapping" for the V1 cluster finding classes
+// if compiled with uncommented "#define V1" line
+// in the AliITSclustererV2.cxx file.
+//
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
+//--------------------------------------------------------------
+#include "AliITSClusterFinderV2.h"
+
+class TClonesArray;
+class AliRawReader;
+class AliITSRawStream;
+
+class AliITSClusterFinderV2SDD : public AliITSClusterFinderV2 {
+public:
+ AliITSClusterFinderV2SDD();
+ virtual ~AliITSClusterFinderV2SDD(){;}
+ virtual void FindRawClusters(Int_t mod);
+ virtual void RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters);
+ protected:
+
+ void FindClustersSDD(TClonesArray *digits);
+ void FindClustersSDD(AliBin* bins[2], Int_t nMaxBin, Int_t nMaxZ,
+ TClonesArray *dig, TClonesArray *clusters=0x0);
+
+ void FindClustersSDD(AliITSRawStream* input,TClonesArray** clusters);
+
+ Int_t fNySDD; //number of "pixels" in Y
+ Int_t fNzSDD; //number of "pixels" in Z
+ Float_t fYpitchSDD; //"pixel size" in Y (drift direction)
+ Float_t fZpitchSDD; //"pixel sizes" in Z
+ Float_t fHwSDD; //half width of the SDD detector
+ Float_t fHlSDD; //half length of the SDD detector
+ Float_t fYoffSDD; //some delay in the drift channel
+
+ ClassDef(AliITSClusterFinderV2SDD,1) // ITS cluster finder V2 for SDD
+};
+
+#endif
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+////////////////////////////////////////////////////////////////////////////
+// Implementation of the ITS clusterer V2 class //
+// //
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
+// //
+///////////////////////////////////////////////////////////////////////////
+
+#include "AliRun.h"
+
+#include "AliITSClusterFinderV2SPD.h"
+#include "AliITSclusterV2.h"
+#include "AliRawReader.h"
+#include "AliITSRawStreamSPD.h"
+
+#include <TClonesArray.h>
+#include "AliITS.h"
+#include "AliITSgeom.h"
+#include "AliITSdigitSPD.h"
+
+ClassImp(AliITSClusterFinderV2SPD)
+
+extern AliRun *gAlice;
+
+AliITSClusterFinderV2SPD::AliITSClusterFinderV2SPD():AliITSClusterFinderV2(){
+
+ //Default constructor
+ AliITSgeom *geom=(AliITSgeom*)fITS->GetITSgeom();
+
+
+ fLastSPD1=geom->GetModuleIndex(2,1,1)-1;
+ fNySPD=256; fNzSPD=160;
+ fYpitchSPD=0.0050;
+ fZ1pitchSPD=0.0425; fZ2pitchSPD=0.0625;
+ fHwSPD=0.64; fHlSPD=3.48;
+ fYSPD[0]=0.5*fYpitchSPD;
+ for (Int_t m=1; m<fNySPD; m++) fYSPD[m]=fYSPD[m-1]+fYpitchSPD;
+ fZSPD[0]=fZ1pitchSPD;
+ for (Int_t m=1; m<fNzSPD; m++) {
+ Double_t dz=fZ1pitchSPD;
+ if (m==31 || m==32 || m==63 || m==64 || m==95 || m==96 ||
+ m==127 || m==128) dz=fZ2pitchSPD;
+ fZSPD[m]=fZSPD[m-1]+dz;
+ }
+ for (Int_t m=0; m<fNzSPD; m++) {
+ Double_t dz=0.5*fZ1pitchSPD;
+ if (m==31 || m==32 || m==63 || m==64 || m==95 || m==96 ||
+ m==127 || m==128) dz=0.5*fZ2pitchSPD;
+ fZSPD[m]-=dz;
+ }
+
+}
+
+
+void AliITSClusterFinderV2SPD::FindRawClusters(Int_t mod){
+
+ //Find clusters V2
+ SetModule(mod);
+ FindClustersSPD(fDigits);
+
+}
+
+void AliITSClusterFinderV2SPD::RawdataToClusters(AliRawReader* rawReader, TClonesArray** clusters){
+ //------------------------------------------------------------
+ // This function creates ITS clusters from raw data
+ //------------------------------------------------------------
+ rawReader->Reset();
+ AliITSRawStreamSPD inputSPD(rawReader);
+ FindClustersSPD(&inputSPD, clusters);
+
+}
+
+Int_t AliITSClusterFinderV2SPD::ClustersSPD(AliBin* bins, TClonesArray* digits,TClonesArray* clusters,Int_t maxBins,Int_t nzbins,Int_t iModule,Bool_t rawdata){
+
+ //Cluster finder for SPD (from digits and from rawdata)
+
+ Int_t nclu=0;
+ for(Int_t iBin =0; iBin < maxBins;iBin++){
+ if(bins[iBin].IsUsed()) continue;
+ Int_t nBins = 0;
+ Int_t idxBins[200];
+ FindCluster(iBin,nzbins,bins,nBins,idxBins);
+ if (nBins == 200){
+ Error("ClustersSPD","SPD Too big cluster !\n");
+ continue;
+ }
+ Int_t milab[10];
+ for(Int_t ilab=0;ilab<10;ilab++){
+ milab[ilab]=-2;
+ }
+ if(rawdata){
+ milab[3]=fNdet[iModule];
+ }
+ Int_t ymin,ymax,zmin,zmax;
+ if(rawdata){
+ ymin = (idxBins[0] / nzbins) - 1;
+ ymax = ymin;
+ zmin = (idxBins[0] % nzbins) - 1;
+ zmax = zmin;
+ }
+ else{
+ AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[0]);
+ ymin=dig->GetCoord2();
+ ymax=ymin;
+ zmin=dig->GetCoord1();
+ zmax=zmin;
+ }
+ for (Int_t idx = 0; idx < nBins; idx++) {
+ Int_t iy;
+ Int_t iz;
+ if(rawdata){
+ iy = (idxBins[idx] / nzbins) - 1;
+ iz = (idxBins[idx] % nzbins) - 1;
+ }
+ else{
+ AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[idx]);
+ iy = dig->GetCoord2();
+ iz = dig->GetCoord1();
+ }
+ if (ymin > iy) ymin = iy;
+ if (ymax < iy) ymax = iy;
+ if (zmin > iz) zmin = iz;
+ if (zmax < iz) zmax = iz;
+
+ }
+ if(!rawdata){
+ for(Int_t l=0;l<nBins;l++){
+ AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[l]);
+ for(Int_t dlab=0;dlab<10;dlab++){
+ Int_t digitlab = (dig->GetTracks())[dlab];
+ if(digitlab<0) continue;
+ AddLabel(milab,digitlab);
+ }
+ if (milab[9]>0) CheckLabels2(milab);
+ }
+ CheckLabels2(milab);
+ }
+
+ Int_t idy =0; //max 2 clusters
+ if((iModule <= fLastSPD1) &&idy<3) idy=3;
+ if((iModule > fLastSPD1) &&idy<4) idy=4;
+ Int_t idz=3;
+ for(Int_t iiz=zmin; iiz<=zmax;iiz+=idz){
+ for(Int_t iiy=ymin;iiy<=ymax;iiy+=idy){
+
+ Int_t ndigits=0;
+ Float_t y=0.,z=0.,q=0.;
+ for(Int_t idx=0;idx<nBins;idx++){
+ Int_t iy;
+ Int_t iz;
+ if(rawdata){
+ iy = (idxBins[idx] / nzbins)-1;
+ iz = (idxBins[idx] % nzbins)-1;
+ }
+ else{
+ AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[idx]);
+ iy = dig->GetCoord2();
+ iz = dig->GetCoord1();
+ }
+ if(zmax-zmin>=idz || ymax-ymin>=idy){
+ if(TMath::Abs(iy-iiy)>0.75*idy) continue;
+ if(TMath::Abs(iz-iiz)>0.75*idz) continue;
+ }
+ ndigits++;
+ Float_t qBin;
+ if(rawdata) qBin = bins[idxBins[idx]].GetQ();
+ if(!rawdata){
+ AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[idx]);
+ qBin = (Float_t)dig->GetSignal();
+ }
+ y+= qBin * fYSPD[iy];
+ z+= qBin * fZSPD[iz];
+ q+= qBin;
+ }// for idx
+ if(ndigits==0) continue;
+ y /= q;
+ z /= q;
+ y -= fHwSPD;
+ z -= fHlSPD;
+ Float_t hit[5]; //y,z,sigma(y)^2, sigma(z)^2, charge
+ hit[0] = -(-y+fYshift[iModule]);
+ if(iModule <= fLastSPD1) hit[0] = -hit[0];
+ hit[1] = -z+fZshift[iModule];
+ hit[2] = fYpitchSPD*fYpitchSPD/12.;
+ hit[3] = fZ1pitchSPD*fZ1pitchSPD/12.;
+ hit[4] = (zmax-zmin+1)*100 + (ymax-ymin+1);
+ if(!rawdata) milab[3]=fNdet[iModule];
+ Int_t info[3] = {ymax-ymin+1,zmax-zmin+1,fNlayer[iModule]};
+ if(!rawdata){
+ AliITSclusterV2 cl(milab,hit,info);
+ fITS->AddClusterV2(cl);
+ }
+ else{
+ Int_t label[4]={milab[0],milab[1],milab[2],milab[3]};
+ new (clusters->AddrAt(nclu))
+ AliITSclusterV2(label, hit,info);
+ }
+ nclu++;
+ }// for iiy
+ }// for iiz
+ }//end for iBin
+ return nclu;
+
+}
+
+
+
+
+void AliITSClusterFinderV2SPD::FindClustersSPD(AliITSRawStream* input,
+ TClonesArray** clusters)
+{
+ //------------------------------------------------------------
+ // Actual SPD cluster finder for raw data
+ //------------------------------------------------------------
+ Int_t nClustersSPD = 0;
+ Int_t kNzBins = fNzSPD + 2;
+ Int_t kNyBins = fNySPD + 2;
+ Int_t kMaxBin = kNzBins * kNyBins;
+ AliBin* bins = NULL;
+
+ // read raw data input stream
+ while (kTRUE) {
+ Bool_t next = input->Next();
+ if (!next || input->IsNewModule()) {
+ Int_t iModule = input->GetPrevModuleID();
+
+ // when all data from a module was read, search for clusters
+ if (bins) {
+ clusters[iModule] = new TClonesArray("AliITSclusterV2");
+ Int_t nClusters = ClustersSPD(bins,0,clusters[iModule],kMaxBin,kNzBins,iModule,kTRUE);
+ nClustersSPD += nClusters;
+ delete bins;
+ }
+
+ if (!next) break;
+ bins = new AliBin[kMaxBin];
+ }
+
+ // fill the current digit into the bins array
+ Int_t index = (input->GetCoord2()+1) * kNzBins + (input->GetCoord1()+1);
+ bins[index].SetIndex(index);
+ bins[index].SetMask(1);
+ bins[index].SetQ(1);
+ }
+
+ Info("FindClustersSPD", "found clusters in ITS SPD: %d", nClustersSPD);
+}
+
+
+
+void AliITSClusterFinderV2SPD::FindClustersSPD(TClonesArray *digits) {
+ //------------------------------------------------------------
+ // Actual SPD cluster finder
+ //------------------------------------------------------------
+
+
+ Int_t kNzBins = fNzSPD + 2;
+ const Int_t kMAXBIN=kNzBins*(fNySPD+2);
+
+ Int_t ndigits=digits->GetEntriesFast();
+ AliBin *bins=new AliBin[kMAXBIN];
+
+ Int_t k;
+ AliITSdigitSPD *d=0;
+ for (k=0; k<ndigits; k++) {
+ d=(AliITSdigitSPD*)digits->UncheckedAt(k);
+ Int_t i=d->GetCoord2()+1; //y
+ Int_t j=d->GetCoord1()+1;
+ Int_t index=i*kNzBins+j;
+ bins[index].SetIndex(k);
+ bins[index].SetMask(1);
+ }
+
+ ClustersSPD(bins,digits,0,kMAXBIN,kNzBins,fModule,kFALSE);
+ delete [] bins;
+}
+
+
--- /dev/null
+#ifndef ALIITSCLUSTERFINDERV2SPD_H
+#define ALIITSCLUSTERFINDERV2SPD_H
+//--------------------------------------------------------------
+// ITS clusterer V2 for SPD
+//
+// This can be a "wrapping" for the V1 cluster finding classes
+// if compiled with uncommented "#define V1" line
+// in the AliITSclustererV2.cxx file.
+//
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
+//--------------------------------------------------------------
+#include "AliITSClusterFinderV2.h"
+
+class TClonesArray;
+class AliRawReader;
+class AliITSRawStream;
+
+class AliITSClusterFinderV2SPD : public AliITSClusterFinderV2 {
+public:
+ AliITSClusterFinderV2SPD();
+ virtual ~AliITSClusterFinderV2SPD(){;}
+ virtual void FindRawClusters(Int_t mod);
+ virtual void RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters);
+
+
+ protected:
+
+ void FindClustersSPD(TClonesArray *digits);
+ void FindClustersSPD(AliITSRawStream* input,TClonesArray** clusters);
+ Int_t ClustersSPD(AliBin* bins, TClonesArray* digits,TClonesArray* clusters,Int_t maxBins, Int_t nzbins,Int_t iModule,Bool_t rawdata=kFALSE);
+
+ Int_t fLastSPD1; //index of the last SPD1 detector
+ Int_t fNySPD; //number of pixels in Y
+ Int_t fNzSPD; //number of pixels in Z
+ Float_t fYpitchSPD; //pixel size in Y
+ Float_t fZ1pitchSPD,fZ2pitchSPD; //pixel sizes in Z
+ Float_t fHwSPD; //half width of the SPD detector
+ Float_t fHlSPD; //half length of the SPD detector
+ Float_t fYSPD[260]; //Y-coordinates of pixel centers
+ Float_t fZSPD[170]; //Z-coordinates of pixel centers
+
+ ClassDef(AliITSClusterFinderV2SPD,1) // ITS cluster finder V2 for SPD
+};
+
+#endif
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+////////////////////////////////////////////////////////////////////////////
+// Implementation of the ITS clusterer V2 class //
+// //
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
+// //
+///////////////////////////////////////////////////////////////////////////
+
+#include "AliRun.h"
+
+#include "AliITSClusterFinderV2SSD.h"
+#include "AliITSclusterV2.h"
+#include "AliRawReader.h"
+#include "AliITSRawStreamSSD.h"
+
+#include <TClonesArray.h>
+#include "AliITS.h"
+#include "AliITSgeom.h"
+#include "AliITSdigitSSD.h"
+
+ClassImp(AliITSClusterFinderV2SSD)
+
+
+AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD():AliITSClusterFinderV2(){
+
+ //Default constructor
+
+ AliITSgeom* geom = (AliITSgeom*)fITS->GetITSgeom();
+
+ fLastSSD1=geom->GetModuleIndex(6,1,1)-1;
+ fYpitchSSD=0.0095;
+ fHwSSD=3.65;
+ fHlSSD=2.00;
+ fTanP=0.0275;
+ fTanN=0.0075;
+
+
+
+}
+
+
+void AliITSClusterFinderV2SSD::FindRawClusters(Int_t mod){
+
+ //Find clusters V2
+ SetModule(mod);
+ FindClustersSSD(fDigits);
+
+}
+
+void AliITSClusterFinderV2SSD::FindClustersSSD(TClonesArray *alldigits) {
+ //------------------------------------------------------------
+ // Actual SSD cluster finder
+ //------------------------------------------------------------
+ Int_t smaxall=alldigits->GetEntriesFast();
+ if (smaxall==0) return;
+ TObjArray *digits = new TObjArray;
+ for (Int_t i=0;i<smaxall; i++){
+ AliITSdigitSSD *d=(AliITSdigitSSD*)alldigits->UncheckedAt(i);
+ if (d->GetSignal()<3) continue;
+ digits->AddLast(d);
+ }
+ Int_t smax = digits->GetEntriesFast();
+ if (smax==0) return;
+
+ const Int_t kMax=1000;
+ Int_t np=0, nn=0;
+ Ali1Dcluster pos[kMax], neg[kMax];
+ Float_t y=0., q=0., qmax=0.;
+ Int_t lab[4]={-2,-2,-2,-2};
+
+ AliITSdigitSSD *d=(AliITSdigitSSD*)digits->UncheckedAt(0);
+ q += d->GetSignal();
+ y += d->GetCoord2()*d->GetSignal();
+ qmax=d->GetSignal();
+ lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
+ Int_t curr=d->GetCoord2();
+ Int_t flag=d->GetCoord1();
+ Int_t *n=&nn;
+ Ali1Dcluster *c=neg;
+ Int_t nd=1;
+ Int_t milab[10];
+ for (Int_t ilab=0;ilab<10;ilab++){
+ milab[ilab]=-2;
+ }
+ milab[0]=d->GetTrack(0); milab[1]=d->GetTrack(1); milab[2]=d->GetTrack(2);
+
+ for (Int_t s=1; s<smax; s++) {
+ d=(AliITSdigitSSD*)digits->UncheckedAt(s);
+ Int_t strip=d->GetCoord2();
+ if ((strip-curr) > 1 || flag!=d->GetCoord1()) {
+ c[*n].SetY(y/q);
+ c[*n].SetQ(q);
+ c[*n].SetNd(nd);
+ CheckLabels2(milab);
+ c[*n].SetLabels(milab);
+ //Split suspiciously big cluster
+ /*
+ if (nd>10&&nd<16){
+ c[*n].SetY(y/q-0.3*nd);
+ c[*n].SetQ(0.5*q);
+ (*n)++;
+ if (*n==MAX) {
+ Error("FindClustersSSD","Too many 1D clusters !");
+ return;
+ }
+ c[*n].SetY(y/q-0.0*nd);
+ c[*n].SetQ(0.5*q);
+ c[*n].SetNd(nd);
+ (*n)++;
+ if (*n==MAX) {
+ Error("FindClustersSSD","Too many 1D clusters !");
+ return;
+ }
+ //
+ c[*n].SetY(y/q+0.3*nd);
+ c[*n].SetQ(0.5*q);
+ c[*n].SetNd(nd);
+ c[*n].SetLabels(milab);
+ }
+ else{
+ */
+ if (nd>4&&nd<25) {
+ c[*n].SetY(y/q-0.25*nd);
+ c[*n].SetQ(0.5*q);
+ (*n)++;
+ if (*n==kMax) {
+ Error("FindClustersSSD","Too many 1D clusters !");
+ return;
+ }
+ c[*n].SetY(y/q+0.25*nd);
+ c[*n].SetQ(0.5*q);
+ c[*n].SetNd(nd);
+ c[*n].SetLabels(milab);
+ }
+ (*n)++;
+ if (*n==kMax) {
+ Error("FindClustersSSD","Too many 1D clusters !");
+ return;
+ }
+ y=q=qmax=0.;
+ nd=0;
+ lab[0]=lab[1]=lab[2]=-2;
+ //
+ for (Int_t ilab=0;ilab<10;ilab++){
+ milab[ilab]=-2;
+ }
+ //
+ if (flag!=d->GetCoord1()) { n=&np; c=pos; }
+ }
+ flag=d->GetCoord1();
+ q += d->GetSignal();
+ y += d->GetCoord2()*d->GetSignal();
+ nd++;
+ if (d->GetSignal()>qmax) {
+ qmax=d->GetSignal();
+ lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
+ }
+ for (Int_t ilab=0;ilab<10;ilab++) {
+ if (d->GetTrack(ilab)>=0) AddLabel(milab, (d->GetTrack(ilab)));
+ }
+ curr=strip;
+ }
+ c[*n].SetY(y/q);
+ c[*n].SetQ(q);
+ c[*n].SetNd(nd);
+ c[*n].SetLabels(lab);
+ //Split suspiciously big cluster
+ if (nd>4 && nd<25) {
+ c[*n].SetY(y/q-0.25*nd);
+ c[*n].SetQ(0.5*q);
+ (*n)++;
+ if (*n==kMax) {
+ Error("FindClustersSSD","Too many 1D clusters !");
+ return;
+ }
+ c[*n].SetY(y/q+0.25*nd);
+ c[*n].SetQ(0.5*q);
+ c[*n].SetNd(nd);
+ c[*n].SetLabels(lab);
+ }
+ (*n)++;
+ if (*n==kMax) {
+ Error("FindClustersSSD","Too many 1D clusters !");
+ return;
+ }
+
+ FindClustersSSD(neg, nn, pos, np);
+}
+
+
+void AliITSClusterFinderV2SSD::RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters){
+
+ //------------------------------------------------------------
+ // This function creates ITS clusters from raw data
+ //------------------------------------------------------------
+ rawReader->Reset();
+ AliITSRawStreamSSD inputSSD(rawReader);
+ FindClustersSSD(&inputSSD,clusters);
+
+}
+
+void AliITSClusterFinderV2SSD::FindClustersSSD(AliITSRawStream* input,
+ TClonesArray** clusters)
+{
+ //------------------------------------------------------------
+ // Actual SSD cluster finder for raw data
+ //------------------------------------------------------------
+ Int_t nClustersSSD = 0;
+ const Int_t kMax = 1000;
+ Ali1Dcluster clusters1D[2][kMax];
+ Int_t nClusters[2] = {0, 0};
+ Int_t lab[3]={-2,-2,-2};
+ Float_t q = 0.;
+ Float_t y = 0.;
+ Int_t nDigits = 0;
+ Int_t prevStrip = -1;
+ Int_t prevFlag = -1;
+ Int_t prevModule = -1;
+
+ // read raw data input stream
+ while (kTRUE) {
+ Bool_t next = input->Next();
+
+ if(input->GetSignal()<3 && next) continue;
+ // check if a new cluster starts
+ Int_t strip = input->GetCoord2();
+ Int_t flag = input->GetCoord1();
+ if ((!next || (input->GetModuleID() != prevModule)||
+ (strip-prevStrip > 1) || (flag != prevFlag)) &&
+ (nDigits > 0)) {
+ if (nClusters[prevFlag] == kMax) {
+ Error("FindClustersSSD", "Too many 1D clusters !");
+ return;
+ }
+ Ali1Dcluster& cluster = clusters1D[prevFlag][nClusters[prevFlag]++];
+ cluster.SetY(y/q);
+ cluster.SetQ(q);
+ cluster.SetNd(nDigits);
+ cluster.SetLabels(lab);
+
+ //Split suspiciously big cluster
+ if (nDigits > 4&&nDigits < 25) {
+ cluster.SetY(y/q - 0.25*nDigits);
+ cluster.SetQ(0.5*q);
+ if (nClusters[prevFlag] == kMax) {
+ Error("FindClustersSSD", "Too many 1D clusters !");
+ return;
+ }
+ Ali1Dcluster& cluster2 = clusters1D[prevFlag][nClusters[prevFlag]++];
+ cluster2.SetY(y/q + 0.25*nDigits);
+ cluster2.SetQ(0.5*q);
+ cluster2.SetNd(nDigits);
+ cluster2.SetLabels(lab);
+ }
+ y = q = 0.;
+ nDigits = 0;
+ }
+
+ if (!next || (input->GetModuleID() != prevModule)) {
+ Int_t iModule = prevModule;
+
+ // when all data from a module was read, search for clusters
+ if (prevFlag >= 0) {
+ clusters[iModule] = new TClonesArray("AliITSclusterV2");
+ fModule = iModule;
+ FindClustersSSD(&clusters1D[0][0], nClusters[0],
+ &clusters1D[1][0], nClusters[1], clusters[iModule]);
+ Int_t nClusters = clusters[iModule]->GetEntriesFast();
+ nClustersSSD += nClusters;
+ }
+
+ if (!next) break;
+ nClusters[0] = nClusters[1] = 0;
+ y = q = 0.;
+ nDigits = 0;
+ }
+
+ // add digit to current cluster
+ q += input->GetSignal();
+ y += strip * input->GetSignal();
+ nDigits++;
+ prevStrip = strip;
+ prevFlag = flag;
+ prevModule = input->GetModuleID();
+
+ }
+
+ Info("FindClustersSSD", "found clusters in ITS SSD: %d", nClustersSSD);
+}
+
+void AliITSClusterFinderV2SSD::
+FindClustersSSD(Ali1Dcluster* neg, Int_t nn,
+ Ali1Dcluster* pos, Int_t np,
+ TClonesArray *clusters) {
+ //------------------------------------------------------------
+ // Actual SSD cluster finder
+ //------------------------------------------------------------
+ TClonesArray &cl=*clusters;
+ //
+ Float_t tanp=fTanP, tann=fTanN;
+ if (fModule>fLastSSD1) {tann=fTanP; tanp=fTanN;}
+ Int_t idet=fNdet[fModule];
+ Int_t ncl=0;
+ //
+ Int_t negativepair[30000];
+ Int_t cnegative[3000];
+ Int_t cused1[3000];
+ Int_t positivepair[30000];
+ Int_t cpositive[3000];
+ Int_t cused2[3000];
+ for (Int_t i=0;i<3000;i++) {cnegative[i]=0; cused1[i]=0;}
+ for (Int_t i=0;i<3000;i++) {cpositive[i]=0; cused2[i]=0;}
+ static Short_t pairs[1000][1000];
+ memset(pairs,0,sizeof(Short_t)*1000000);
+// Short_t ** pairs = new Short_t*[1000];
+// for (Int_t i=0; i<1000; i++) {
+// pairs[i] = new Short_t[1000];
+// memset(pairs[i],0,sizeof(Short_t)*1000);
+// }
+ //
+ // find available pairs
+ //
+ for (Int_t i=0; i<np; i++) {
+ Float_t yp=pos[i].GetY()*fYpitchSSD;
+ if (pos[i].GetQ()<3) continue;
+ for (Int_t j=0; j<nn; j++) {
+ if (neg[j].GetQ()<3) continue;
+ Float_t yn=neg[j].GetY()*fYpitchSSD;
+ Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Float_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ if (TMath::Abs(yt)<fHwSSD+0.01)
+ if (TMath::Abs(zt)<fHlSSD+0.01*(neg[j].GetNd()+pos[i].GetNd())) {
+ negativepair[i*10+cnegative[i]] =j; //index
+ positivepair[j*10+cpositive[j]] =i;
+ cnegative[i]++; //counters
+ cpositive[j]++;
+ pairs[i][j]=100;
+ }
+ }
+ }
+ //
+ for (Int_t i=0; i<np; i++) {
+ Float_t yp=pos[i].GetY()*fYpitchSSD;
+ if (pos[i].GetQ()<3) continue;
+ for (Int_t j=0; j<nn; j++) {
+ if (neg[j].GetQ()<3) continue;
+ if (cpositive[j]&&cnegative[i]) continue;
+ Float_t yn=neg[j].GetY()*fYpitchSSD;
+ Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Float_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ if (TMath::Abs(yt)<fHwSSD+0.1)
+ if (TMath::Abs(zt)<fHlSSD+0.15) {
+ if (cnegative[i]==0) pos[i].SetNd(100); // not available pair
+ if (cpositive[j]==0) neg[j].SetNd(100); // not available pair
+ negativepair[i*10+cnegative[i]] =j; //index
+ positivepair[j*10+cpositive[j]] =i;
+ cnegative[i]++; //counters
+ cpositive[j]++;
+ pairs[i][j]=100;
+ }
+ }
+ }
+ //
+ Float_t lp[5];
+ Int_t milab[10];
+ Double_t ratio;
+
+ //
+ // sign gold tracks
+ //
+ for (Int_t ip=0;ip<np;ip++){
+ Float_t ybest=1000,zbest=1000,qbest=0;
+ //
+ // select gold clusters
+ if ( (cnegative[ip]==1) && cpositive[negativepair[10*ip]]==1){
+ Float_t yp=pos[ip].GetY()*fYpitchSSD;
+ Int_t j = negativepair[10*ip];
+ ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
+ //
+ Float_t yn=neg[j].GetY()*fYpitchSSD;
+ Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Float_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ ybest=yt; zbest=zt;
+ qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
+ lp[0]=-(-ybest+fYshift[fModule]);
+ lp[1]= -zbest+fZshift[fModule];
+ lp[2]=0.0025*0.0025; //SigmaY2
+ lp[3]=0.110*0.110; //SigmaZ2
+
+ lp[4]=qbest; //Q
+ for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
+ for (Int_t ilab=0;ilab<3;ilab++){
+ milab[ilab] = pos[ip].GetLabel(ilab);
+ milab[ilab+3] = neg[j].GetLabel(ilab);
+ }
+ //
+ CheckLabels2(milab);
+ milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
+ Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
+ AliITSclusterV2 * cl2;
+ if(clusters) cl2 = new (cl[ncl]) AliITSclusterV2(milab,lp,info);
+ else{
+ cl2 = new AliITSclusterV2(milab,lp,info);
+ fITS->AddClusterV2(*cl2);
+ }
+ ncl++;
+ cl2->SetChargeRatio(ratio);
+ cl2->SetType(1);
+ pairs[ip][j]=1;
+ if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
+ cl2->SetType(2);
+ pairs[ip][j]=2;
+ }
+ cused1[ip]++;
+ cused2[j]++;
+ }
+ }
+
+ for (Int_t ip=0;ip<np;ip++){
+ Float_t ybest=1000,zbest=1000,qbest=0;
+ //
+ //
+ // select "silber" cluster
+ if ( cnegative[ip]==1 && cpositive[negativepair[10*ip]]==2){
+ Int_t in = negativepair[10*ip];
+ Int_t ip2 = positivepair[10*in];
+ if (ip2==ip) ip2 = positivepair[10*in+1];
+ Float_t pcharge = pos[ip].GetQ()+pos[ip2].GetQ();
+ if (TMath::Abs(pcharge-neg[in].GetQ())<10){
+ //
+ // add first pair
+ if (pairs[ip][in]==100){ //
+ Float_t yp=pos[ip].GetY()*fYpitchSSD;
+ Float_t yn=neg[in].GetY()*fYpitchSSD;
+ Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Float_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ ybest =yt; zbest=zt;
+ qbest =pos[ip].GetQ();
+ lp[0]=-(-ybest+fYshift[fModule]);
+ lp[1]= -zbest+fZshift[fModule];
+ lp[2]=0.0025*0.0025; //SigmaY2
+ lp[3]=0.110*0.110; //SigmaZ2
+
+ lp[4]=qbest; //Q
+ for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
+ for (Int_t ilab=0;ilab<3;ilab++){
+ milab[ilab] = pos[ip].GetLabel(ilab);
+ milab[ilab+3] = neg[in].GetLabel(ilab);
+ }
+ //
+ CheckLabels2(milab);
+ ratio = (pos[ip].GetQ()-neg[in].GetQ())/(pos[ip].GetQ()+neg[in].GetQ());
+ milab[3]=(((ip<<10) + in)<<10) + idet; // pos|neg|det
+ Int_t info[3] = {pos[ip].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
+
+ AliITSclusterV2 * cl2;
+ if(clusters) cl2 = new (cl[ncl]) AliITSclusterV2(milab,lp,info);
+ else{
+ cl2 = new AliITSclusterV2(milab,lp,info);
+ fITS->AddClusterV2(*cl2);
+ }
+ ncl++;
+ cl2->SetChargeRatio(ratio);
+ cl2->SetType(5);
+ pairs[ip][in] = 5;
+ if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster
+ cl2->SetType(6);
+ pairs[ip][in] = 6;
+ }
+ }
+ //
+ // add second pair
+
+ // if (!(cused1[ip2] || cused2[in])){ //
+ if (pairs[ip2][in]==100){
+ Float_t yp=pos[ip2].GetY()*fYpitchSSD;
+ Float_t yn=neg[in].GetY()*fYpitchSSD;
+ Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Float_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ ybest =yt; zbest=zt;
+ qbest =pos[ip2].GetQ();
+ lp[0]=-(-ybest+fYshift[fModule]);
+ lp[1]= -zbest+fZshift[fModule];
+ lp[2]=0.0025*0.0025; //SigmaY2
+ lp[3]=0.110*0.110; //SigmaZ2
+
+ lp[4]=qbest; //Q
+ for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
+ for (Int_t ilab=0;ilab<3;ilab++){
+ milab[ilab] = pos[ip2].GetLabel(ilab);
+ milab[ilab+3] = neg[in].GetLabel(ilab);
+ }
+ //
+ CheckLabels2(milab);
+ ratio = (pos[ip2].GetQ()-neg[in].GetQ())/(pos[ip2].GetQ()+neg[in].GetQ());
+ milab[3]=(((ip2<<10) + in)<<10) + idet; // pos|neg|det
+ Int_t info[3] = {pos[ip2].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
+
+ AliITSclusterV2 * cl2;
+ if(clusters) cl2 = new (cl[ncl]) AliITSclusterV2(milab,lp,info);
+ else{
+ cl2 = new AliITSclusterV2(milab,lp,info);
+ fITS->AddClusterV2(*cl2);
+ }
+ ncl++;
+ cl2->SetChargeRatio(ratio);
+ cl2->SetType(5);
+ pairs[ip2][in] =5;
+ if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster
+ cl2->SetType(6);
+ pairs[ip2][in] =6;
+ }
+ }
+ cused1[ip]++;
+ cused1[ip2]++;
+ cused2[in]++;
+ }
+ }
+ }
+
+ //
+ for (Int_t jn=0;jn<nn;jn++){
+ if (cused2[jn]) continue;
+ Float_t ybest=1000,zbest=1000,qbest=0;
+ // select "silber" cluster
+ if ( cpositive[jn]==1 && cnegative[positivepair[10*jn]]==2){
+ Int_t ip = positivepair[10*jn];
+ Int_t jn2 = negativepair[10*ip];
+ if (jn2==jn) jn2 = negativepair[10*ip+1];
+ Float_t pcharge = neg[jn].GetQ()+neg[jn2].GetQ();
+ //
+ if (TMath::Abs(pcharge-pos[ip].GetQ())<10){
+ //
+ // add first pair
+ // if (!(cused1[ip]||cused2[jn])){
+ if (pairs[ip][jn]==100){
+ Float_t yn=neg[jn].GetY()*fYpitchSSD;
+ Float_t yp=pos[ip].GetY()*fYpitchSSD;
+ Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Float_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ ybest =yt; zbest=zt;
+ qbest =neg[jn].GetQ();
+ lp[0]=-(-ybest+fYshift[fModule]);
+ lp[1]= -zbest+fZshift[fModule];
+ lp[2]=0.0025*0.0025; //SigmaY2
+ lp[3]=0.110*0.110; //SigmaZ2
+
+ lp[4]=qbest; //Q
+ for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
+ for (Int_t ilab=0;ilab<3;ilab++){
+ milab[ilab] = pos[ip].GetLabel(ilab);
+ milab[ilab+3] = neg[jn].GetLabel(ilab);
+ }
+ //
+ CheckLabels2(milab);
+ ratio = (pos[ip].GetQ()-neg[jn].GetQ())/(pos[ip].GetQ()+neg[jn].GetQ());
+ milab[3]=(((ip<<10) + jn)<<10) + idet; // pos|neg|det
+ Int_t info[3] = {pos[ip].GetNd(),neg[jn].GetNd(),fNlayer[fModule]};
+
+ AliITSclusterV2 * cl2;
+ if(clusters) cl2 = new (cl[ncl]) AliITSclusterV2(milab,lp,info);
+ else{
+ cl2 = new AliITSclusterV2(milab,lp,info);
+ fITS->AddClusterV2(*cl2);
+ }
+ ncl++;
+ cl2->SetChargeRatio(ratio);
+ cl2->SetType(7);
+ pairs[ip][jn] =7;
+ if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster
+ cl2->SetType(8);
+ pairs[ip][jn]=8;
+ }
+ }
+ //
+ // add second pair
+ // if (!(cused1[ip]||cused2[jn2])){
+ if (pairs[ip][jn2]==100){
+ Float_t yn=neg[jn2].GetY()*fYpitchSSD;
+ Double_t yp=pos[ip].GetY()*fYpitchSSD;
+ Double_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Double_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ ybest =yt; zbest=zt;
+ qbest =neg[jn2].GetQ();
+ lp[0]=-(-ybest+fYshift[fModule]);
+ lp[1]= -zbest+fZshift[fModule];
+ lp[2]=0.0025*0.0025; //SigmaY2
+ lp[3]=0.110*0.110; //SigmaZ2
+
+ lp[4]=qbest; //Q
+ for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
+ for (Int_t ilab=0;ilab<3;ilab++){
+ milab[ilab] = pos[ip].GetLabel(ilab);
+ milab[ilab+3] = neg[jn2].GetLabel(ilab);
+ }
+ //
+ CheckLabels2(milab);
+ ratio = (pos[ip].GetQ()-neg[jn2].GetQ())/(pos[ip].GetQ()+neg[jn2].GetQ());
+ milab[3]=(((ip<<10) + jn2)<<10) + idet; // pos|neg|det
+ Int_t info[3] = {pos[ip].GetNd(),neg[jn2].GetNd(),fNlayer[fModule]};
+ AliITSclusterV2 * cl2;
+ if(clusters) cl2 = new (cl[ncl]) AliITSclusterV2(milab,lp,info);
+ else{
+ cl2 = new AliITSclusterV2(milab,lp,info);
+ fITS->AddClusterV2(*cl2);
+ }
+
+ ncl++;
+ cl2->SetChargeRatio(ratio);
+ pairs[ip][jn2]=7;
+ cl2->SetType(7);
+ if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster
+ cl2->SetType(8);
+ pairs[ip][jn2]=8;
+ }
+ }
+ cused1[ip]++;
+ cused2[jn]++;
+ cused2[jn2]++;
+ }
+ }
+ }
+
+ for (Int_t ip=0;ip<np;ip++){
+ Float_t ybest=1000,zbest=1000,qbest=0;
+ //
+ // 2x2 clusters
+ //
+ if ( (cnegative[ip]<5) && cpositive[negativepair[10*ip]]<5){
+ Float_t minchargediff =4.;
+ Int_t j=-1;
+ for (Int_t di=0;di<cnegative[ip];di++){
+ Int_t jc = negativepair[ip*10+di];
+ Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
+ if (TMath::Abs(chargedif)<minchargediff){
+ j =jc;
+ minchargediff = TMath::Abs(chargedif);
+ }
+ }
+ if (j<0) continue; // not proper cluster
+ Int_t count =0;
+ for (Int_t di=0;di<cnegative[ip];di++){
+ Int_t jc = negativepair[ip*10+di];
+ Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
+ if (TMath::Abs(chargedif)<minchargediff+3.) count++;
+ }
+ if (count>1) continue; // more than one "proper" cluster for positive
+ //
+ count =0;
+ for (Int_t dj=0;dj<cpositive[j];dj++){
+ Int_t ic = positivepair[j*10+dj];
+ Float_t chargedif = pos[ic].GetQ()-neg[j].GetQ();
+ if (TMath::Abs(chargedif)<minchargediff+3.) count++;
+ }
+ if (count>1) continue; // more than one "proper" cluster for negative
+
+ Int_t jp = 0;
+
+ count =0;
+ for (Int_t dj=0;dj<cnegative[jp];dj++){
+ Int_t ic = positivepair[jp*10+dj];
+ Float_t chargedif = pos[ic].GetQ()-neg[jp].GetQ();
+ if (TMath::Abs(chargedif)<minchargediff+4.) count++;
+ }
+ if (count>1) continue;
+ if (pairs[ip][j]<100) continue;
+ //
+ //almost gold clusters
+ Float_t yp=pos[ip].GetY()*fYpitchSSD;
+ Float_t yn=neg[j].GetY()*fYpitchSSD;
+ Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Float_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ ybest=yt; zbest=zt;
+ qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
+ lp[0]=-(-ybest+fYshift[fModule]);
+ lp[1]= -zbest+fZshift[fModule];
+ lp[2]=0.0025*0.0025; //SigmaY2
+ lp[3]=0.110*0.110; //SigmaZ2
+ lp[4]=qbest; //Q
+ for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
+ for (Int_t ilab=0;ilab<3;ilab++){
+ milab[ilab] = pos[ip].GetLabel(ilab);
+ milab[ilab+3] = neg[j].GetLabel(ilab);
+ }
+ //
+ CheckLabels2(milab);
+ ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
+ milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
+ Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
+ AliITSclusterV2 * cl2;
+ if(clusters) cl2 = new (cl[ncl]) AliITSclusterV2(milab,lp,info);
+ else{
+ cl2 = new AliITSclusterV2(milab,lp,info);
+ fITS->AddClusterV2(*cl2);
+ }
+
+
+ ncl++;
+ cl2->SetChargeRatio(ratio);
+ cl2->SetType(10);
+ pairs[ip][j]=10;
+ if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
+ cl2->SetType(11);
+ pairs[ip][j]=11;
+ }
+ cused1[ip]++;
+ cused2[j]++;
+ }
+
+ }
+
+ //
+ for (Int_t i=0; i<np; i++) {
+ Float_t ybest=1000,zbest=1000,qbest=0;
+ Float_t yp=pos[i].GetY()*fYpitchSSD;
+ if (pos[i].GetQ()<3) continue;
+ for (Int_t j=0; j<nn; j++) {
+ // for (Int_t di = 0;di<cpositive[i];di++){
+ // Int_t j = negativepair[10*i+di];
+ if (neg[j].GetQ()<3) continue;
+ if (cused2[j]||cused1[i]) continue;
+ if (pairs[i][j]>0 &&pairs[i][j]<100) continue;
+ ratio = (pos[i].GetQ()-neg[j].GetQ())/(pos[i].GetQ()+neg[j].GetQ());
+ Float_t yn=neg[j].GetY()*fYpitchSSD;
+ Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
+ Float_t yt=yn + tann*zt;
+ zt-=fHlSSD; yt-=fHwSSD;
+ if (TMath::Abs(yt)<fHwSSD+0.01)
+ if (TMath::Abs(zt)<fHlSSD+0.01*(neg[j].GetNd()+pos[i].GetNd())) {
+ ybest=yt; zbest=zt;
+ qbest=0.5*(pos[i].GetQ()+neg[j].GetQ());
+ lp[0]=-(-ybest+fYshift[fModule]);
+ lp[1]= -zbest+fZshift[fModule];
+ lp[2]=0.0025*0.0025; //SigmaY2
+ lp[3]=0.110*0.110; //SigmaZ2
+
+ lp[4]=qbest; //Q
+ for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
+ for (Int_t ilab=0;ilab<3;ilab++){
+ milab[ilab] = pos[i].GetLabel(ilab);
+ milab[ilab+3] = neg[j].GetLabel(ilab);
+ }
+ //
+ CheckLabels2(milab);
+ milab[3]=(((i<<10) + j)<<10) + idet; // pos|neg|det
+ Int_t info[3] = {pos[i].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
+ AliITSclusterV2 * cl2;
+ if(clusters) cl2 = new (cl[ncl]) AliITSclusterV2(milab,lp,info);
+ else{
+ cl2 = new AliITSclusterV2(milab,lp,info);
+ fITS->AddClusterV2(*cl2);
+ }
+ ncl++;
+ cl2->SetChargeRatio(ratio);
+ cl2->SetType(100+cpositive[j]+cnegative[i]);
+ //cl2->SetType(0);
+ /*
+ if (pairs[i][j]<100){
+ printf("problem:- %d\n", pairs[i][j]);
+ }
+ if (cnegative[i]<2&&cpositive[j]<2){
+ printf("problem:- %d\n", pairs[i][j]);
+ }
+ */
+ }
+ }
+ }
+
+// for (Int_t i=0; i<1000; i++) delete [] pairs[i];
+// delete [] pairs;
+
+}
+
+
+
--- /dev/null
+#ifndef ALIITSCLUSTERFINDERV2SSD_H
+#define ALIITSCLUSTERFINDERV2SSD_H
+//--------------------------------------------------------------
+// ITS clusterer V2 for SSD
+//
+// This can be a "wrapping" for the V1 cluster finding classes
+// if compiled with uncommented "#define V1" line
+// in the AliITSclustererV2.cxx file.
+//
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
+//--------------------------------------------------------------
+#include "AliITSClusterFinderV2.h"
+
+class TClonesArray;
+class AliRawReader;
+class AliITSRawStream;
+
+class AliITSClusterFinderV2SSD : public AliITSClusterFinderV2 {
+public:
+ AliITSClusterFinderV2SSD();
+ virtual ~AliITSClusterFinderV2SSD(){;}
+ virtual void FindRawClusters(Int_t mod);
+ virtual void RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters);
+ protected:
+
+ void FindClustersSSD(TClonesArray *digits);
+ void FindClustersSSD(Ali1Dcluster* neg, Int_t nn,
+ Ali1Dcluster* pos, Int_t np,
+ TClonesArray *clusters=0x0);
+
+ void FindClustersSSD(AliITSRawStream* input,TClonesArray** clusters);
+
+ Int_t fLastSSD1; //index of the last SSD1 detector
+ Float_t fYpitchSSD; //strip pitch (cm)
+ Float_t fHwSSD; //half-width of an SSD detector (cm)
+ Float_t fHlSSD; //half-length of an SSD detector (cm)
+ Float_t fTanP; //tangent of the stereo angle on the P side
+ Float_t fTanN; //tangent of the stereo angle on the N side
+
+
+ ClassDef(AliITSClusterFinderV2SSD,1) // ITS cluster finder V2 for SDD
+};
+
+#endif
#pragma link C++ class AliITSEventHeader+;
#pragma link C++ class AliITSReconstructor+;
+#pragma link C++ class AliITSClusterFinderV2+;
+#pragma link C++ class AliITSClusterFinderV2SDD+;
+#pragma link C++ class AliITSClusterFinderV2SPD+;
+#pragma link C++ class AliITSClusterFinderV2SSD+;
+
//
//#pragma link C++ class AliACORDEFunction+;
AliITSBeamTestDigitizer.cxx \
AliITSEventHeader.cxx \
AliITSRawStreamSSDv1.cxx \
+ AliITSClusterFinderV2.cxx \
+ AliITSClusterFinderV2SDD.cxx \
+ AliITSClusterFinderV2SPD.cxx \
+ AliITSClusterFinderV2SSD.cxx \
+
# AliITSv11.cxx \
# AliITSv11Geometry.cxx \
# AliITSv11GeometrySupport.cxx \
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff