1 /***************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Conributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
20 ////////////////////////////////////////////////////////////////////////
21 // This class defines the "Standard" reconstruction for the ITS //
24 ////////////////////////////////////////////////////////////////////////
25 #include "TObjArray.h"
28 #include "AliCDBManager.h"
29 #include "AliCDBEntry.h"
30 #include "AliITSClusterFinder.h"
31 #include "AliITSClusterFinderV2SPD.h"
32 #include "AliITSClusterFinderV2SDD.h"
33 #include "AliITSClusterFinderV2SSD.h"
34 #include "AliITSDetTypeRec.h"
35 #include "AliITSDDLModuleMapSDD.h"
36 #include "AliITSRecPoint.h"
37 #include "AliITSCalibrationSDD.h"
38 #include "AliITSMapSDD.h"
39 #include "AliITSCalibrationSSD.h"
40 #include "AliITSNoiseSSDv2.h"
41 #include "AliITSGainSSDv2.h"
42 #include "AliITSBadChannelsSSDv2.h"
43 #include "AliITSNoiseSSD.h"
44 #include "AliITSGainSSD.h"
45 #include "AliITSBadChannelsSSD.h"
46 #include "AliITSresponseSDD.h"
47 #include "AliITSsegmentationSPD.h"
48 #include "AliITSsegmentationSDD.h"
49 #include "AliITSsegmentationSSD.h"
51 #include "AliITSRawStreamSPD.h"
52 #include "AliITSTriggerConditions.h"
53 #include "AliITSFOSignalsSPD.h"
54 #include "AliRunLoader.h"
55 #include "AliDataLoader.h"
56 #include "AliITSLoader.h"
58 class AliITSDriftSpeedArraySDD;
59 class AliITSCorrMapSDD;
60 class AliITSRecoParam;
62 const Int_t AliITSDetTypeRec::fgkNdettypes = 3;
63 const Int_t AliITSDetTypeRec::fgkDefaultNModulesSPD = 240;
64 const Int_t AliITSDetTypeRec::fgkDefaultNModulesSDD = 260;
65 const Int_t AliITSDetTypeRec::fgkDefaultNModulesSSD = 1698;
67 ClassImp(AliITSDetTypeRec)
69 //________________________________________________________________
70 AliITSDetTypeRec::AliITSDetTypeRec(): TObject(),
78 fTriggerConditions(0),
88 fFastOrFiredMap(1200){
89 // Standard Constructor
97 fReconstruction = new TObjArray(fgkNdettypes);
98 fDigits = new TObjArray(fgkNdettypes);
99 for(Int_t i=0; i<3; i++){
102 fSSDCalibration=new AliITSCalibrationSSD();
103 fNMod = new Int_t [fgkNdettypes];
104 fNMod[0] = fgkDefaultNModulesSPD;
105 fNMod[1] = fgkDefaultNModulesSDD;
106 fNMod[2] = fgkDefaultNModulesSSD;
107 fRecPoints = new TClonesArray("AliITSRecPoint",3000);
113 //______________________________________________________________________
114 AliITSDetTypeRec::AliITSDetTypeRec(const AliITSDetTypeRec & rec):TObject(rec),
116 fITSgeom(rec.fITSgeom),
117 fReconstruction(rec.fReconstruction),
118 fSegmentation(rec.fSegmentation),
119 fCalibration(rec.fCalibration),
120 fSSDCalibration(rec.fSSDCalibration),
121 fSPDDead(rec.fSPDDead),
122 fTriggerConditions(rec.fTriggerConditions),
123 fDigits(rec.fDigits),
124 fFOSignals(rec.fFOSignals),
125 fDDLMapSDD(rec.fDDLMapSDD),
126 fRespSDD(rec.fRespSDD),
127 fAveGainSDD(rec.fAveGainSDD),
128 fRecPoints(rec.fRecPoints),
129 fNRecPoints(rec.fNRecPoints),
130 fFirstcall(rec.fFirstcall),
131 fLoadOnlySPDCalib(rec.fLoadOnlySPDCalib),
132 fFastOrFiredMap(rec.fFastOrFiredMap){
137 //______________________________________________________________________
138 AliITSDetTypeRec& AliITSDetTypeRec::operator=(const AliITSDetTypeRec& source){
139 // Assignment operator.
140 this->~AliITSDetTypeRec();
141 new(this) AliITSDetTypeRec(source);
146 //_____________________________________________________________________
147 AliITSDetTypeRec::~AliITSDetTypeRec(){
151 fReconstruction->Delete();
152 delete fReconstruction;
156 fSegmentation->Delete();
157 delete fSegmentation;
161 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
162 fCalibration->Delete();
165 if(fRespSDD) delete fRespSDD;
166 if(fDDLMapSDD) delete fDDLMapSDD;
169 if(fSSDCalibration) delete fSSDCalibration;
171 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
177 if(fTriggerConditions){
178 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
179 fTriggerConditions->Delete();
180 delete fTriggerConditions;
181 fTriggerConditions = 0;
190 fRecPoints->Delete();
196 if (fITSgeom) delete fITSgeom;
200 //___________________________________________________________________
201 void AliITSDetTypeRec::SetReconstructionModel(Int_t dettype,AliITSClusterFinder *clf){
203 //Set reconstruction model for detector type
205 if(fReconstruction==0) fReconstruction = new TObjArray(fgkNdettypes);
206 if(fReconstruction->At(dettype)!=0) delete fReconstruction->At(dettype);
207 fReconstruction->AddAt(clf,dettype);
209 //______________________________________________________________________
210 AliITSClusterFinder* AliITSDetTypeRec::GetReconstructionModel(Int_t dettype) const{
212 //Get reconstruction model for detector type
213 if(fReconstruction==0) {
214 Warning("GetReconstructionModel","fReconstruction is 0!");
217 return (AliITSClusterFinder*)fReconstruction->At(dettype);
220 //______________________________________________________________________
221 void AliITSDetTypeRec::SetSegmentationModel(Int_t dettype,AliITSsegmentation *seg){
223 //Set segmentation model for detector type
225 if(fSegmentation==0) fSegmentation = new TObjArray(fgkNdettypes);
226 if(fSegmentation->At(dettype)!=0) delete fSegmentation->At(dettype);
227 fSegmentation->AddAt(seg,dettype);
230 //______________________________________________________________________
231 AliITSsegmentation* AliITSDetTypeRec::GetSegmentationModel(Int_t dettype) const {
233 //Get segmentation model for detector type
235 if(fSegmentation==0) {
236 Warning("GetSegmentationModel","fSegmentation is 0!");
239 return (AliITSsegmentation*)fSegmentation->At(dettype);
242 //_______________________________________________________________________
243 void AliITSDetTypeRec::SetCalibrationModel(Int_t iMod, AliITSCalibration *cal){
245 //Set calibration (response) for the module iMod of type dettype
246 if (fCalibration==0) {
247 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
248 fCalibration->SetOwner(kTRUE);
249 fCalibration->Clear();
252 if (fCalibration->At(iMod) != 0)
253 delete (AliITSCalibration*) fCalibration->At(iMod);
254 fCalibration->AddAt(cal,iMod);
257 //_______________________________________________________________________
258 void AliITSDetTypeRec::SetSPDDeadModel(Int_t iMod, AliITSCalibration *cal){
260 //Set dead pixel info for the SPD module iMod
262 fSPDDead = new TObjArray(fgkDefaultNModulesSPD);
263 fSPDDead->SetOwner(kTRUE);
267 if (fSPDDead->At(iMod) != 0)
268 delete (AliITSCalibration*) fSPDDead->At(iMod);
269 fSPDDead->AddAt(cal,iMod);
271 //_______________________________________________________________________
272 AliITSCalibration* AliITSDetTypeRec::GetCalibrationModel(Int_t iMod) const {
274 //Get calibration model for module type
276 if(fCalibration==0) {
277 Warning("GetalibrationModel","fCalibration is 0!");
281 if(iMod<fgkDefaultNModulesSPD+fgkDefaultNModulesSDD){
282 return (AliITSCalibration*)fCalibration->At(iMod);
284 Int_t i=iMod-(fgkDefaultNModulesSPD+fgkDefaultNModulesSDD);
285 fSSDCalibration->SetModule(i);
286 return (AliITSCalibration*)fSSDCalibration;
290 //_______________________________________________________________________
291 AliITSCalibration* AliITSDetTypeRec::GetSPDDeadModel(Int_t iMod) const {
293 //Get SPD dead for module iMod
296 AliWarning("fSPDDead is 0!");
299 return (AliITSCalibration*)fSPDDead->At(iMod);
301 //_______________________________________________________________________
302 AliITSTriggerConditions* AliITSDetTypeRec::GetTriggerConditions() const {
303 //Get Pixel Trigger Conditions
304 if (fTriggerConditions==0) {
305 AliWarning("fTriggerConditions is 0!");
307 return fTriggerConditions;
309 //______________________________________________________________________
310 void AliITSDetTypeRec::SetTreeAddressD(TTree* const treeD){
311 // Set branch address for the tree of digits.
313 const char *det[4] = {"SPD","SDD","SSD","ITS"};
315 const Char_t* digclass;
320 if (fDigits == 0x0) fDigits = new TObjArray(fgkNdettypes);
321 for (i=0; i<fgkNdettypes; i++) {
322 digclass = GetDigitClassName(i);
323 if(!(fDigits->At(i))) {
324 fDigits->AddAt(new TClonesArray(digclass,1000),i);
329 if (fgkNdettypes==3) sprintf(branchname,"%sDigits%s",det[3],det[i]);
330 else sprintf(branchname,"%sDigits%d",det[3],i+1);
332 branch = treeD->GetBranch(branchname);
333 if (branch) branch->SetAddress(&((*fDigits)[i]));
339 //_______________________________________________________________________
340 TBranch* AliITSDetTypeRec::MakeBranchInTree(TTree* const tree,
341 const char* name, const char *classname,
342 void* address,Int_t size,Int_t splitlevel)
345 // Makes branch in given tree and diverts them to a separate file
351 Error("MakeBranchInTree","Making Branch %s Tree is NULL",name);
354 TBranch *branch = tree->GetBranch(name);
359 branch = tree->Branch(name,classname,address,size,splitlevel);
362 branch = tree->Bronch(name, "TClonesArray", address, size, splitlevel);
368 //____________________________________________________________________
369 void AliITSDetTypeRec::SetDefaults(){
371 //Set defaults for segmentation and response
374 Warning("SetDefaults","null pointer to AliITSgeomGeom !");
378 AliITSsegmentation* seg;
379 if(!GetCalibration()) {AliFatal("Exit");exit(0);}
381 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
383 seg = new AliITSsegmentationSPD();
384 SetSegmentationModel(dettype,seg);
385 SetDigitClassName(dettype,"AliITSdigitSPD");
387 if(fLoadOnlySPDCalib==kFALSE){
389 seg = new AliITSsegmentationSDD();
390 AliITSCalibrationSDD* cal=(AliITSCalibrationSDD*)GetCalibrationModel(fgkDefaultNModulesSPD+1);
391 if(cal->IsAMAt20MHz()){
392 seg->SetPadSize(seg->Dpz(0),20.);
393 seg->SetNPads(seg->Npz()/2,128);
395 SetSegmentationModel(dettype,seg);
396 SetDigitClassName(dettype,"AliITSdigitSDD");
400 AliITSsegmentationSSD* seg2 = new AliITSsegmentationSSD();
401 SetSegmentationModel(dettype,seg2);
402 SetDigitClassName(dettype,"AliITSdigitSSD");
406 //______________________________________________________________________
407 Bool_t AliITSDetTypeRec::GetCalibration() {
408 // Get Default calibration if a storage is not defined.
411 AliITSCalibration* cal = GetCalibrationModel(0);
417 // SetRunNumber((Int_t)AliCDBManager::Instance()->GetRun());
418 // Int_t run=GetRunNumber();
420 Bool_t cacheStatus = AliCDBManager::Instance()->GetCacheFlag();
421 if (fCalibration==0) {
422 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
423 fCalibration->SetOwner(!cacheStatus);
424 fCalibration->Clear();
427 Bool_t retCode=GetCalibrationSPD(cacheStatus);
428 if(retCode==kFALSE) return kFALSE;
430 if(fLoadOnlySPDCalib==kFALSE){
431 retCode=GetCalibrationSDD(cacheStatus);
432 if(retCode==kFALSE) return kFALSE;
433 retCode=GetCalibrationSSD(cacheStatus);
434 if(retCode==kFALSE) return kFALSE;
437 AliInfo(Form("%i SPD, %i SDD and %i SSD in calibration database",
438 fNMod[0], fNMod[1], fNMod[2]));
441 //______________________________________________________________________
442 Bool_t AliITSDetTypeRec::GetCalibrationSPD(Bool_t cacheStatus) {
443 // Get SPD calibration objects from OCDB
444 // dead pixel are not used for local reconstruction
447 AliCDBEntry *noisySPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDNoisy");
448 AliCDBEntry *deadSPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDDead");
449 AliCDBEntry *pitCond = AliCDBManager::Instance()->Get("TRIGGER/SPD/PITConditions");
450 if(!noisySPD || !deadSPD || !pitCond ){
451 AliFatal("SPD Calibration object retrieval failed! ");
455 TObjArray *calNoisySPD = (TObjArray*) noisySPD->GetObject();
456 if (!cacheStatus) noisySPD->SetObject(NULL);
457 noisySPD->SetOwner(kTRUE);
459 TObjArray *calDeadSPD = (TObjArray*) deadSPD->GetObject();
460 if (!cacheStatus) deadSPD->SetObject(NULL);
461 deadSPD->SetOwner(kTRUE);
463 AliITSTriggerConditions *calPitCond = (AliITSTriggerConditions*) pitCond->GetObject();
464 if (!cacheStatus) pitCond->SetObject(NULL);
465 pitCond->SetOwner(kTRUE);
472 if ((!calNoisySPD) || (!calDeadSPD) || (!calPitCond)){
473 AliWarning("Can not get SPD calibration from calibration database !");
477 fNMod[0] = calNoisySPD->GetEntries();
479 AliITSCalibration* cal;
480 for (Int_t i=0; i<fNMod[0]; i++) {
481 cal = (AliITSCalibration*) calNoisySPD->At(i);
482 SetCalibrationModel(i, cal);
483 cal = (AliITSCalibration*) calDeadSPD->At(i);
484 SetSPDDeadModel(i, cal);
486 fTriggerConditions = calPitCond;
491 //______________________________________________________________________
492 Bool_t AliITSDetTypeRec::GetCalibrationSDD(Bool_t cacheStatus) {
493 // Get SDD calibration objects from OCDB
495 AliCDBEntry *entrySDD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSDD");
496 AliCDBEntry *entry2SDD = AliCDBManager::Instance()->Get("ITS/Calib/RespSDD");
497 AliCDBEntry *drSpSDD = AliCDBManager::Instance()->Get("ITS/Calib/DriftSpeedSDD");
498 AliCDBEntry *ddlMapSDD = AliCDBManager::Instance()->Get("ITS/Calib/DDLMapSDD");
499 // AliCDBEntry *mapASDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsAnodeSDD");
500 AliCDBEntry *mapTSDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsTimeSDD");
502 if(!entrySDD || !entry2SDD || !drSpSDD || !ddlMapSDD || !mapTSDD ){
503 AliFatal("SDD Calibration object retrieval failed! ");
509 TObjArray *calSDD = (TObjArray *)entrySDD->GetObject();
510 if(!cacheStatus)entrySDD->SetObject(NULL);
511 entrySDD->SetOwner(kTRUE);
513 AliITSresponseSDD *pSDD = (AliITSresponseSDD*)entry2SDD->GetObject();
514 if(!cacheStatus)entry2SDD->SetObject(NULL);
515 entry2SDD->SetOwner(kTRUE);
517 TObjArray *drSp = (TObjArray *)drSpSDD->GetObject();
518 if(!cacheStatus)drSpSDD->SetObject(NULL);
519 drSpSDD->SetOwner(kTRUE);
521 AliITSDDLModuleMapSDD *ddlsdd=(AliITSDDLModuleMapSDD*)ddlMapSDD->GetObject();
522 if(!cacheStatus)ddlMapSDD->SetObject(NULL);
523 ddlMapSDD->SetOwner(kTRUE);
525 // TObjArray *mapAn = (TObjArray *)mapASDD->GetObject();
526 // if(!cacheStatus)mapASDD->SetObject(NULL);
527 // mapASDD->SetOwner(kTRUE);
529 TObjArray *mapT = (TObjArray *)mapTSDD->GetObject();
530 if(!cacheStatus)mapTSDD->SetObject(NULL);
531 mapTSDD->SetOwner(kTRUE);
534 // DB entries are deleted. In this way metadeta objects are deleted as well
544 if ((!pSDD)||(!calSDD) || (!drSp) || (!ddlsdd) || (!mapT) ){
545 AliWarning("Can not get SDD calibration from calibration database !");
549 fNMod[1] = calSDD->GetEntries();
553 AliITSCalibration* cal;
556 Bool_t oldMapFormat=kFALSE;
557 TObject* objmap=(TObject*)mapT->At(0);
558 TString cname(objmap->ClassName());
559 if(cname.CompareTo("AliITSMapSDD")==0){
561 AliInfo("SDD Maps converted to new format");
563 for(Int_t iddl=0; iddl<AliITSDDLModuleMapSDD::GetNDDLs(); iddl++){
564 for(Int_t icar=0; icar<AliITSDDLModuleMapSDD::GetNModPerDDL();icar++){
565 Int_t iMod=fDDLMapSDD->GetModuleNumber(iddl,icar);
566 if(iMod==-1) continue;
567 Int_t i=iMod - fgkDefaultNModulesSPD;
568 cal = (AliITSCalibration*) calSDD->At(i);
571 for(Int_t iAnode=0;iAnode< ((AliITSCalibrationSDD*)cal)->NOfAnodes(); iAnode++){
572 if(((AliITSCalibrationSDD*)cal)->IsBadChannel(iAnode)) continue;
573 avegain+= ((AliITSCalibrationSDD*)cal)->GetChannelGain(iAnode);
576 AliITSDriftSpeedArraySDD* arr0 = (AliITSDriftSpeedArraySDD*) drSp->At(i0);
577 AliITSDriftSpeedArraySDD* arr1 = (AliITSDriftSpeedArraySDD*) drSp->At(i1);
579 AliITSCorrMapSDD* mt0 = 0;
580 AliITSCorrMapSDD* mt1 = 0;
582 AliITSMapSDD* oldmap0=(AliITSMapSDD*)mapT->At(i0);
583 AliITSMapSDD* oldmap1=(AliITSMapSDD*)mapT->At(i1);
584 mt0=oldmap0->ConvertToNewFormat();
585 mt1=oldmap1->ConvertToNewFormat();
587 mt0=(AliITSCorrMapSDD*)mapT->At(i0);
588 mt1=(AliITSCorrMapSDD*)mapT->At(i1);
590 cal->SetDriftSpeed(0,arr0);
591 cal->SetDriftSpeed(1,arr1);
594 SetCalibrationModel(iMod, cal);
597 if(nGdAnodes) fAveGainSDD=avegain/nGdAnodes;
602 //______________________________________________________________________
603 Bool_t AliITSDetTypeRec::GetCalibrationSSD(Bool_t cacheStatus) {
604 // Get SSD calibration objects from OCDB
605 // AliCDBEntry *entrySSD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSSD");
607 AliCDBEntry *entryNoiseSSD = AliCDBManager::Instance()->Get("ITS/Calib/NoiseSSD");
608 AliCDBEntry *entryGainSSD = AliCDBManager::Instance()->Get("ITS/Calib/GainSSD");
609 AliCDBEntry *entryBadChannelsSSD = AliCDBManager::Instance()->Get("ITS/Calib/BadChannelsSSD");
611 if(!entryNoiseSSD || !entryGainSSD || !entryBadChannelsSSD){
612 AliFatal("SSD Calibration object retrieval failed! ");
616 TObject *emptyssd = 0; TString ssdobjectname;
617 AliITSNoiseSSDv2 *noiseSSD = new AliITSNoiseSSDv2();
618 emptyssd = (TObject *)entryNoiseSSD->GetObject();
619 ssdobjectname = emptyssd->GetName();
620 if(ssdobjectname=="TObjArray") {
621 TObjArray *noiseSSDOld = (TObjArray *)entryNoiseSSD->GetObject();
622 ReadOldSSDNoise(noiseSSDOld, noiseSSD);
624 else if(ssdobjectname=="AliITSNoiseSSDv2")
625 noiseSSD = (AliITSNoiseSSDv2 *)entryNoiseSSD->GetObject();
626 if(!cacheStatus)entryNoiseSSD->SetObject(NULL);
627 entryNoiseSSD->SetOwner(kTRUE);
629 AliITSGainSSDv2 *gainSSD = new AliITSGainSSDv2();
630 emptyssd = (TObject *)entryGainSSD->GetObject();
631 ssdobjectname = emptyssd->GetName();
632 if(ssdobjectname=="Gain") {
633 TObjArray *gainSSDOld = (TObjArray *)entryGainSSD->GetObject();
634 ReadOldSSDGain(gainSSDOld, gainSSD);
636 else if(ssdobjectname=="AliITSGainSSDv2")
637 gainSSD = (AliITSGainSSDv2 *)entryGainSSD->GetObject();
638 if(!cacheStatus)entryGainSSD->SetObject(NULL);
639 entryGainSSD->SetOwner(kTRUE);
641 AliITSBadChannelsSSDv2 *badChannelsSSD = new AliITSBadChannelsSSDv2();
642 emptyssd = (TObject *)entryBadChannelsSSD->GetObject();
643 ssdobjectname = emptyssd->GetName();
644 if(ssdobjectname=="TObjArray") {
645 TObjArray *badChannelsSSDOld = (TObjArray *)entryBadChannelsSSD->GetObject();
646 ReadOldSSDBadChannels(badChannelsSSDOld, badChannelsSSD);
648 else if(ssdobjectname=="AliITSBadChannelsSSDv2")
649 badChannelsSSD = (AliITSBadChannelsSSDv2*)entryBadChannelsSSD->GetObject();
650 if(!cacheStatus)entryBadChannelsSSD->SetObject(NULL);
651 entryBadChannelsSSD->SetOwner(kTRUE);
653 // DB entries are deleted. In this way metadeta objects are deleted as well
655 delete entryNoiseSSD;
657 delete entryBadChannelsSSD;
660 if ((!noiseSSD)|| (!gainSSD)|| (!badChannelsSSD)) {
661 AliWarning("Can not get SSD calibration from calibration database !");
665 fSSDCalibration->SetNoise(noiseSSD);
666 fSSDCalibration->SetGain(gainSSD);
667 fSSDCalibration->SetBadChannels(badChannelsSSD);
668 //fSSDCalibration->FillBadChipMap();
673 //________________________________________________________________
674 void AliITSDetTypeRec::SetDefaultClusterFindersV2(Bool_t rawdata){
676 //Set defaults for cluster finder V2
679 Warning("SetDefaults","Null pointer to AliITSgeom !");
683 AliITSClusterFinder *clf;
685 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
688 if(!GetReconstructionModel(dettype)){
689 clf = new AliITSClusterFinderV2SPD(this);
691 if(!rawdata) clf->SetDigits(DigitsAddress(0));
692 SetReconstructionModel(dettype,clf);
698 if(!GetReconstructionModel(dettype)){
699 clf = new AliITSClusterFinderV2SDD(this);
701 if(!rawdata) clf->SetDigits(DigitsAddress(1));
702 SetReconstructionModel(dettype,clf);
709 if(!GetReconstructionModel(dettype)){
710 clf = new AliITSClusterFinderV2SSD(this);
712 if(!rawdata) clf->SetDigits(DigitsAddress(2));
713 SetReconstructionModel(dettype,clf);
720 //______________________________________________________________________
721 void AliITSDetTypeRec::MakeBranch(TTree* tree, Option_t* option){
723 //Creates branches for clusters and recpoints
724 Bool_t cR = (strstr(option,"R")!=0);
725 Bool_t cRF = (strstr(option,"RF")!=0);
729 if(cR) MakeBranchR(tree);
730 if(cRF) MakeBranchRF(tree);
734 //___________________________________________________________________
735 void AliITSDetTypeRec::ResetDigits(){
736 // Reset number of digits and the digits array for the ITS detector.
739 for(Int_t i=0;i<fgkNdettypes;i++){
743 //___________________________________________________________________
744 void AliITSDetTypeRec::ResetDigits(Int_t branch){
745 // Reset number of digits and the digits array for this branch.
747 if(fDigits->At(branch)) ((TClonesArray*)fDigits->At(branch))->Clear();
750 //__________________________________________________________________
751 void AliITSDetTypeRec::MakeBranchR(TTree *treeR, Option_t *opt){
753 //Creates tree branches for recpoints
755 // cont char *file File name where RecPoints branch is to be written
756 // to. If blank it write the SDigits to the same
757 // file in which the Hits were found.
762 // only one branch for rec points for all detector types
763 Bool_t oFast= (strstr(opt,"Fast")!=0);
765 Char_t detname[10] = "ITS";
769 sprintf(branchname,"%sRecPointsF",detname);
771 sprintf(branchname,"%sRecPoints",detname);
774 if(!fRecPoints)fRecPoints = new TClonesArray("AliITSRecPoint",1000);
776 MakeBranchInTree(treeR,branchname,0,&fRecPoints,buffsz,99);
778 //______________________________________________________________________
779 void AliITSDetTypeRec::SetTreeAddressR(TTree* const treeR){
780 // Set branch address for the Reconstructed points Trees.
782 // TTree *treeR Tree containing the RecPoints.
788 Char_t namedet[10]="ITS";
791 if(fRecPoints==0x0) fRecPoints = new TClonesArray("AliITSRecPoint",1000);
793 sprintf(branchname,"%sRecPoints",namedet);
794 branch = treeR->GetBranch(branchname);
796 branch->SetAddress(&fRecPoints);
799 sprintf(branchname,"%sRecPointsF",namedet);
800 branch = treeR->GetBranch(branchname);
802 branch->SetAddress(&fRecPoints);
806 //____________________________________________________________________
807 void AliITSDetTypeRec::AddRecPoint(const AliITSRecPoint &r){
808 // Add a reconstructed space point to the list
810 // const AliITSRecPoint &r RecPoint class to be added to the tree
811 // of reconstructed points TreeR.
817 TClonesArray &lrecp = *fRecPoints;
818 new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
821 //______________________________________________________________________
822 void AliITSDetTypeRec::DigitsToRecPoints(TTree *treeD,TTree *treeR,Int_t lastentry,Option_t *opt, Int_t optCluFind){
823 // cluster finding and reconstruction of space points
824 // the condition below will disappear when the geom class will be
825 // initialized for all versions - for the moment it is only for v5 !
826 // 7 is the SDD beam test version
828 // TTree *treeD Digits tree
829 // TTree *treeR Clusters tree
830 // Int_t lastentry Offset for module when not all of the modules
832 // Option_t *opt String indicating which ITS sub-detectors should
833 // be processed. If ="All" then all of the ITS
834 // sub detectors are processed.
836 const char *all = strstr(opt,"All");
837 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
840 SetDefaultClusterFindersV2();
841 AliDebug(1,"V2 cluster finder has been selected \n");
843 SetDefaultClusterFindersV2();
844 AliInfo("Cluster Finder Option not implemented, V2 cluster finder will be used \n");
848 // Reset Fast-OR fired map
849 ResetFastOrFiredMap();
851 if (all || det[0]) { // SPD present
852 // Get the FO signals for this event
853 AliRunLoader* runLoader = AliRunLoader::Instance();
854 AliITSLoader* itsLoader = (AliITSLoader*) runLoader->GetLoader("ITSLoader");
856 AliError("ITS loader is NULL.");
859 fFOSignals = (AliITSFOSignalsSPD*)itsLoader->TreeD()->GetUserInfo()->FindObject("AliITSFOSignalsSPD");
860 if(!fFOSignals) AliError("FO signals not retrieved");
866 AliITSClusterFinder *rec = 0;
867 Int_t id,module,first=0;
868 for(module=0;module<GetITSgeom()->GetIndexMax();module++){
869 id = GetITSgeom()->GetModuleType(module);
870 if (!all && !det[id]) continue;
871 if(det[id]) first = GetITSgeom()->GetStartDet(id);
872 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
873 TClonesArray *itsDigits = DigitsAddress(id);
875 AliFatal("The reconstruction class was not instanciated!");
876 ResetDigits(); // MvL: Not sure we neeed this when rereading anyways
878 treeD->GetEvent(lastentry+module);
881 treeD->GetEvent(lastentry+(module-first));
883 Int_t ndigits = itsDigits->GetEntriesFast();
884 if (ndigits>0 || id==0) { // for SPD we always want to call FindRawClusters (to process FO signals)
885 rec->SetDetTypeRec(this);
886 rec->SetDigits(DigitsAddress(id));
887 // rec->SetClusters(ClustersAddress(id));
888 rec->FindRawClusters(module);
894 // Remove PIT in-active chips from Fast-OR fired map
895 if (all || det[0]) { // SPD present
896 RemoveFastOrFiredInActive();
899 //______________________________________________________________________
900 void AliITSDetTypeRec::DigitsToRecPoints(AliRawReader* rawReader,TTree *treeR,Option_t *opt){
901 // cluster finding and reconstruction of space points
902 // the condition below will disappear when the geom class will be
903 // initialized for all versions - for the moment it is only for v5 !
904 // 7 is the SDD beam test version
906 // AliRawReader *rawReader Pointer to the raw-data reader
907 // TTree *treeR Clusters tree
912 const char *all = strstr(opt,"All");
913 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
918 TClonesArray *array=new TClonesArray("AliITSRecPoint",1000);
919 TBranch *branch = treeR->Branch("ITSRecPoints",&array);
922 TClonesArray** clusters = new TClonesArray*[GetITSgeom()->GetIndexMax()];
923 for (Int_t iModule = 0; iModule < GetITSgeom()->GetIndexMax(); iModule++) {
924 clusters[iModule] = NULL;
927 DigitsToRecPoints(rawReader,clusters,opt);
930 TClonesArray *emptyArray=new TClonesArray("AliITSRecPoint");
931 for(Int_t iModule=0;iModule<GetITSgeom()->GetIndexMax();iModule++){
932 id = GetITSgeom()->GetModuleType(iModule);
933 if (!all && !det[id]) continue;
934 array = clusters[iModule];
936 AliDebug(1,Form("data for module %d missing!",iModule));
939 branch->SetAddress(&array);
941 nClusters+=array->GetEntriesFast();
943 if (array != emptyArray) {
951 Info("DigitsToRecPoints", "total number of found recpoints in ITS: %d\n",
955 //______________________________________________________________________
956 void AliITSDetTypeRec::DigitsToRecPoints(AliRawReader* rawReader,TClonesArray** clusters,Option_t *opt){
957 // cluster finding and reconstruction of space points
958 // the condition below will disappear when the geom class will be
959 // initialized for all versions - for the moment it is only for v5 !
960 // 7 is the SDD beam test version
962 // AliRawReader *rawReader Pointer to the raw-data reader
963 // TClonesArray **clusters Clusters Array
968 const char *all = strstr(opt,"All");
969 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
972 // Reset Fast-OR fired map
973 ResetFastOrFiredMap();
975 AliITSClusterFinder *rec = 0;
979 if (!all && !det[id]) continue;
980 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
982 AliFatal("The reconstruction class was not instantiated");
983 rec->SetDetTypeRec(this);
984 rec->RawdataToClusters(rawReader,clusters);
987 // Remove PIT in-active chips from Fast-OR fired map
988 if (all || det[0]) { // SPD present
989 RemoveFastOrFiredInActive();
992 //______________________________________________________________________
993 void AliITSDetTypeRec::ReadOldSSDNoise(const TObjArray *array,
994 AliITSNoiseSSDv2 *noiseSSD) {
995 //Reads the old SSD calibration object and converts it to the new format
996 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
997 const Int_t fgkSSDPSIDESTRIPSPERMODULE = 768;
999 Int_t gNMod = array->GetEntries();
1000 cout<<"Converting old calibration object for noise..."<<endl;
1003 Double_t noise = 0.0;
1004 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1005 AliITSNoiseSSD *noiseModule = (AliITSNoiseSSD*) (array->At(iModule));
1006 for(Int_t iStrip = 0; iStrip < fgkSSDSTRIPSPERMODULE; iStrip++) {
1007 noise = (iStrip < fgkSSDPSIDESTRIPSPERMODULE) ? noiseModule->GetNoiseP(iStrip) : noiseModule->GetNoiseN(1535 - iStrip);
1008 if(iStrip < fgkSSDPSIDESTRIPSPERMODULE)
1009 noiseSSD->AddNoiseP(iModule,iStrip,noise);
1010 if(iStrip >= fgkSSDPSIDESTRIPSPERMODULE)
1011 noiseSSD->AddNoiseN(iModule,1535 - iStrip,noise);
1013 }//loop over modules
1016 //______________________________________________________________________
1017 void AliITSDetTypeRec::ReadOldSSDBadChannels(const TObjArray *array,
1018 AliITSBadChannelsSSDv2 *badChannelsSSD) {
1019 //Reads the old SSD calibration object and converts it to the new format
1020 Int_t gNMod = array->GetEntries();
1021 cout<<"Converting old calibration object for bad channels..."<<endl;
1022 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1023 //for (Int_t iModule = 0; iModule < 1; iModule++) {
1024 AliITSBadChannelsSSD *bad = (AliITSBadChannelsSSD*) (array->At(iModule));
1025 TArrayI arrayPSide = bad->GetBadPChannelsList();
1026 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1027 badChannelsSSD->AddBadChannelP(iModule,
1029 (Char_t)arrayPSide.At(iPCounter));
1031 TArrayI arrayNSide = bad->GetBadNChannelsList();
1032 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1033 badChannelsSSD->AddBadChannelN(iModule,
1035 (Char_t)arrayNSide.At(iNCounter));
1037 }//loop over modules
1040 //______________________________________________________________________
1041 void AliITSDetTypeRec::ReadOldSSDGain(const TObjArray *array,
1042 AliITSGainSSDv2 *gainSSD) {
1043 //Reads the old SSD calibration object and converts it to the new format
1045 Int_t gNMod = array->GetEntries();
1046 cout<<"Converting old calibration object for gain..."<<endl;
1049 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1050 AliITSGainSSD *gainModule = (AliITSGainSSD*) (array->At(iModule));
1051 TArrayF arrayPSide = gainModule->GetGainP();
1052 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1053 gainSSD->AddGainP(iModule,
1055 arrayPSide.At(iPCounter));
1056 TArrayF arrayNSide = gainModule->GetGainN();
1057 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1058 gainSSD->AddGainN(iModule,
1060 arrayNSide.At(iNCounter));
1061 }//loop over modules
1063 //______________________________________________________________________
1064 void AliITSDetTypeRec::RemoveFastOrFiredInActive() {
1065 // Removes the chips that were in-active in the pixel trigger (from fast-or fired map)
1066 if (fTriggerConditions==NULL) {
1067 AliError("Pixel trigger conditions are missing.");
1073 while (fTriggerConditions->GetNextInActiveChip(eq,hs,chip)) {
1074 UInt_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1075 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1078 //______________________________________________________________________
1079 void AliITSDetTypeRec::SetFastOrFiredMapOnline(UInt_t eq, UInt_t hs, UInt_t chip) {
1080 // Set fast-or fired map for this chip
1081 Int_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1082 return SetFastOrFiredMap(chipKey);