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 "AliITSRecPointContainer.h"
38 #include "AliITSCalibrationSDD.h"
39 #include "AliITSMapSDD.h"
40 #include "AliITSCalibrationSSD.h"
41 #include "AliITSNoiseSSDv2.h"
42 #include "AliITSGainSSDv2.h"
43 #include "AliITSBadChannelsSSDv2.h"
44 #include "AliITSNoiseSSD.h"
45 #include "AliITSGainSSD.h"
46 #include "AliITSBadChannelsSSD.h"
47 #include "AliITSresponseSDD.h"
48 #include "AliITSsegmentationSPD.h"
49 #include "AliITSsegmentationSDD.h"
50 #include "AliITSsegmentationSSD.h"
52 #include "AliITSRawStreamSPD.h"
53 #include "AliITSTriggerConditions.h"
54 #include "AliITSFOSignalsSPD.h"
55 #include "AliRunLoader.h"
56 #include "AliDataLoader.h"
57 #include "AliITSLoader.h"
59 class AliITSDriftSpeedArraySDD;
60 class AliITSCorrMapSDD;
61 class AliITSRecoParam;
63 const Int_t AliITSDetTypeRec::fgkNdettypes = 3;
64 const Int_t AliITSDetTypeRec::fgkDefaultNModulesSPD = 240;
65 const Int_t AliITSDetTypeRec::fgkDefaultNModulesSDD = 260;
66 const Int_t AliITSDetTypeRec::fgkDefaultNModulesSSD = 1698;
68 ClassImp(AliITSDetTypeRec)
70 //________________________________________________________________
71 AliITSDetTypeRec::AliITSDetTypeRec(): TObject(),
79 fTriggerConditions(0),
89 fFastOrFiredMap(1200){
90 // Standard Constructor
98 fReconstruction = new TObjArray(fgkNdettypes);
99 fDigits = new TObjArray(fgkNdettypes);
100 for(Int_t i=0; i<3; i++){
103 fSSDCalibration=new AliITSCalibrationSSD();
104 fNMod = new Int_t [fgkNdettypes];
105 fNMod[0] = fgkDefaultNModulesSPD;
106 fNMod[1] = fgkDefaultNModulesSDD;
107 fNMod[2] = fgkDefaultNModulesSSD;
108 fRecPoints = new TClonesArray("AliITSRecPoint",3000);
114 //______________________________________________________________________
115 AliITSDetTypeRec::AliITSDetTypeRec(const AliITSDetTypeRec & rec):TObject(rec),
117 fITSgeom(rec.fITSgeom),
118 fReconstruction(rec.fReconstruction),
119 fSegmentation(rec.fSegmentation),
120 fCalibration(rec.fCalibration),
121 fSSDCalibration(rec.fSSDCalibration),
122 fSPDDead(rec.fSPDDead),
123 fTriggerConditions(rec.fTriggerConditions),
124 fDigits(rec.fDigits),
125 fFOSignals(rec.fFOSignals),
126 fDDLMapSDD(rec.fDDLMapSDD),
127 fRespSDD(rec.fRespSDD),
128 fAveGainSDD(rec.fAveGainSDD),
129 fRecPoints(rec.fRecPoints),
130 fNRecPoints(rec.fNRecPoints),
131 fFirstcall(rec.fFirstcall),
132 fLoadOnlySPDCalib(rec.fLoadOnlySPDCalib),
133 fFastOrFiredMap(rec.fFastOrFiredMap){
136 for(Int_t i=0; i<3; i++){
137 fkDigClassName[i]=rec.fkDigClassName[i]; // NB only copies Char_t*, so not so safe, but this code should never be reached anyways
140 //______________________________________________________________________
141 AliITSDetTypeRec& AliITSDetTypeRec::operator=(const AliITSDetTypeRec& source){
142 // Assignment operator.
143 this->~AliITSDetTypeRec();
144 new(this) AliITSDetTypeRec(source);
149 //_____________________________________________________________________
150 AliITSDetTypeRec::~AliITSDetTypeRec(){
154 fReconstruction->Delete();
155 delete fReconstruction;
159 fSegmentation->Delete();
160 delete fSegmentation;
164 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
165 fCalibration->Delete();
168 if(fRespSDD) delete fRespSDD;
169 if(fDDLMapSDD) delete fDDLMapSDD;
173 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
174 delete fSSDCalibration;
175 fSSDCalibration = NULL;
179 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
185 if(fTriggerConditions){
186 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
187 fTriggerConditions->Delete();
188 delete fTriggerConditions;
189 fTriggerConditions = 0;
198 fRecPoints->Delete();
204 if (fITSgeom) delete fITSgeom;
208 //___________________________________________________________________
209 void AliITSDetTypeRec::SetReconstructionModel(Int_t dettype,AliITSClusterFinder *clf){
211 //Set reconstruction model for detector type
213 if(fReconstruction==0) fReconstruction = new TObjArray(fgkNdettypes);
214 if(fReconstruction->At(dettype)!=0) delete fReconstruction->At(dettype);
215 fReconstruction->AddAt(clf,dettype);
217 //______________________________________________________________________
218 AliITSClusterFinder* AliITSDetTypeRec::GetReconstructionModel(Int_t dettype) const{
220 //Get reconstruction model for detector type
221 if(fReconstruction==0) {
222 Warning("GetReconstructionModel","fReconstruction is 0!");
225 return (AliITSClusterFinder*)fReconstruction->At(dettype);
228 //______________________________________________________________________
229 void AliITSDetTypeRec::SetSegmentationModel(Int_t dettype,AliITSsegmentation *seg){
231 //Set segmentation model for detector type
233 if(fSegmentation==0) fSegmentation = new TObjArray(fgkNdettypes);
234 if(fSegmentation->At(dettype)!=0) delete fSegmentation->At(dettype);
235 fSegmentation->AddAt(seg,dettype);
238 //______________________________________________________________________
239 AliITSsegmentation* AliITSDetTypeRec::GetSegmentationModel(Int_t dettype) const {
241 //Get segmentation model for detector type
243 if(fSegmentation==0) {
244 Warning("GetSegmentationModel","fSegmentation is 0!");
247 return (AliITSsegmentation*)fSegmentation->At(dettype);
250 //_______________________________________________________________________
251 void AliITSDetTypeRec::SetCalibrationModel(Int_t iMod, AliITSCalibration *cal){
253 //Set calibration (response) for the module iMod of type dettype
254 if (fCalibration==0) {
255 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
256 fCalibration->SetOwner(kTRUE);
257 fCalibration->Clear();
260 if (fCalibration->At(iMod) != 0)
261 delete (AliITSCalibration*) fCalibration->At(iMod);
262 fCalibration->AddAt(cal,iMod);
265 //_______________________________________________________________________
266 void AliITSDetTypeRec::SetSPDDeadModel(Int_t iMod, AliITSCalibration *cal){
268 //Set dead pixel info for the SPD module iMod
270 fSPDDead = new TObjArray(fgkDefaultNModulesSPD);
271 fSPDDead->SetOwner(kTRUE);
275 if (fSPDDead->At(iMod) != 0)
276 delete (AliITSCalibration*) fSPDDead->At(iMod);
277 fSPDDead->AddAt(cal,iMod);
279 //_______________________________________________________________________
280 AliITSCalibration* AliITSDetTypeRec::GetCalibrationModel(Int_t iMod) const {
282 //Get calibration model for module type
284 if(fCalibration==0) {
285 Warning("GetalibrationModel","fCalibration is 0!");
289 if(iMod<fgkDefaultNModulesSPD+fgkDefaultNModulesSDD){
290 return (AliITSCalibration*)fCalibration->At(iMod);
292 Int_t i=iMod-(fgkDefaultNModulesSPD+fgkDefaultNModulesSDD);
293 fSSDCalibration->SetModule(i);
294 return (AliITSCalibration*)fSSDCalibration;
298 //_______________________________________________________________________
299 AliITSCalibration* AliITSDetTypeRec::GetSPDDeadModel(Int_t iMod) const {
301 //Get SPD dead for module iMod
304 AliWarning("fSPDDead is 0!");
307 return (AliITSCalibration*)fSPDDead->At(iMod);
309 //_______________________________________________________________________
310 AliITSTriggerConditions* AliITSDetTypeRec::GetTriggerConditions() const {
311 //Get Pixel Trigger Conditions
312 if (fTriggerConditions==0) {
313 AliWarning("fTriggerConditions is 0!");
315 return fTriggerConditions;
317 //______________________________________________________________________
318 void AliITSDetTypeRec::SetTreeAddressD(TTree* const treeD){
319 // Set branch address for the tree of digits.
321 const char *det[4] = {"SPD","SDD","SSD","ITS"};
323 const Char_t* digclass;
328 if (fDigits == 0x0) fDigits = new TObjArray(fgkNdettypes);
329 for (i=0; i<fgkNdettypes; i++) {
330 digclass = GetDigitClassName(i);
331 if(!(fDigits->At(i))) {
332 fDigits->AddAt(new TClonesArray(digclass,1000),i);
337 if (fgkNdettypes==3) sprintf(branchname,"%sDigits%s",det[3],det[i]);
338 else sprintf(branchname,"%sDigits%d",det[3],i+1);
340 branch = treeD->GetBranch(branchname);
341 if (branch) branch->SetAddress(&((*fDigits)[i]));
347 //_______________________________________________________________________
348 TBranch* AliITSDetTypeRec::MakeBranchInTree(TTree* const tree,
349 const char* name, const char *classname,
350 void* address,Int_t size,Int_t splitlevel)
353 // Makes branch in given tree and diverts them to a separate file
359 Error("MakeBranchInTree","Making Branch %s Tree is NULL",name);
362 TBranch *branch = tree->GetBranch(name);
367 branch = tree->Branch(name,classname,address,size,splitlevel);
370 branch = tree->Bronch(name, "TClonesArray", address, size, splitlevel);
376 //____________________________________________________________________
377 void AliITSDetTypeRec::SetDefaults(){
379 //Set defaults for segmentation and response
382 Warning("SetDefaults","null pointer to AliITSgeomGeom !");
386 AliITSsegmentation* seg;
387 if(!GetCalibration()) {AliFatal("Exit");exit(0);}
389 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
391 seg = new AliITSsegmentationSPD();
392 SetSegmentationModel(dettype,seg);
393 SetDigitClassName(dettype,"AliITSdigitSPD");
396 seg = new AliITSsegmentationSDD();
397 if(fLoadOnlySPDCalib==kFALSE){
398 AliITSCalibrationSDD* cal=(AliITSCalibrationSDD*)GetCalibrationModel(fgkDefaultNModulesSPD+1);
399 if(cal->IsAMAt20MHz()){
400 seg->SetPadSize(seg->Dpz(0),20.);
401 seg->SetNPads(seg->Npz()/2,128);
404 SetSegmentationModel(dettype,seg);
405 SetDigitClassName(dettype,"AliITSdigitSDD");
408 AliITSsegmentationSSD* seg2 = new AliITSsegmentationSSD();
409 SetSegmentationModel(dettype,seg2);
410 SetDigitClassName(dettype,"AliITSdigitSSD");
414 //______________________________________________________________________
415 Bool_t AliITSDetTypeRec::GetCalibration() {
416 // Get Default calibration if a storage is not defined.
419 AliITSCalibration* cal = GetCalibrationModel(0);
425 // SetRunNumber((Int_t)AliCDBManager::Instance()->GetRun());
426 // Int_t run=GetRunNumber();
428 Bool_t cacheStatus = AliCDBManager::Instance()->GetCacheFlag();
429 if (fCalibration==0) {
430 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
431 fCalibration->SetOwner(!cacheStatus);
432 fCalibration->Clear();
435 Bool_t retCode=GetCalibrationSPD(cacheStatus);
436 if(retCode==kFALSE) return kFALSE;
438 if(fLoadOnlySPDCalib==kFALSE){
439 retCode=GetCalibrationSDD(cacheStatus);
440 if(retCode==kFALSE) return kFALSE;
441 retCode=GetCalibrationSSD(cacheStatus);
442 if(retCode==kFALSE) return kFALSE;
445 AliInfo(Form("%i SPD, %i SDD and %i SSD in calibration database",
446 fNMod[0], fNMod[1], fNMod[2]));
449 //______________________________________________________________________
450 Bool_t AliITSDetTypeRec::GetCalibrationSPD(Bool_t cacheStatus) {
451 // Get SPD calibration objects from OCDB
452 // dead pixel are not used for local reconstruction
455 AliCDBEntry *noisySPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDNoisy");
456 AliCDBEntry *deadSPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDDead");
457 AliCDBEntry *pitCond = AliCDBManager::Instance()->Get("TRIGGER/SPD/PITConditions");
458 if(!noisySPD || !deadSPD || !pitCond ){
459 AliFatal("SPD Calibration object retrieval failed! ");
463 TObjArray *calNoisySPD = (TObjArray*) noisySPD->GetObject();
464 if (!cacheStatus) noisySPD->SetObject(NULL);
465 noisySPD->SetOwner(kTRUE);
467 TObjArray *calDeadSPD = (TObjArray*) deadSPD->GetObject();
468 if (!cacheStatus) deadSPD->SetObject(NULL);
469 deadSPD->SetOwner(kTRUE);
471 AliITSTriggerConditions *calPitCond = (AliITSTriggerConditions*) pitCond->GetObject();
472 if (!cacheStatus) pitCond->SetObject(NULL);
473 pitCond->SetOwner(kTRUE);
480 if ((!calNoisySPD) || (!calDeadSPD) || (!calPitCond)){
481 AliWarning("Can not get SPD calibration from calibration database !");
485 fNMod[0] = calNoisySPD->GetEntries();
487 AliITSCalibration* cal;
488 for (Int_t i=0; i<fNMod[0]; i++) {
489 cal = (AliITSCalibration*) calNoisySPD->At(i);
490 SetCalibrationModel(i, cal);
491 cal = (AliITSCalibration*) calDeadSPD->At(i);
492 SetSPDDeadModel(i, cal);
494 fTriggerConditions = calPitCond;
499 //______________________________________________________________________
500 Bool_t AliITSDetTypeRec::GetCalibrationSDD(Bool_t cacheStatus) {
501 // Get SDD calibration objects from OCDB
503 AliCDBEntry *entrySDD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSDD");
504 AliCDBEntry *entry2SDD = AliCDBManager::Instance()->Get("ITS/Calib/RespSDD");
505 AliCDBEntry *drSpSDD = AliCDBManager::Instance()->Get("ITS/Calib/DriftSpeedSDD");
506 AliCDBEntry *ddlMapSDD = AliCDBManager::Instance()->Get("ITS/Calib/DDLMapSDD");
507 // AliCDBEntry *mapASDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsAnodeSDD");
508 AliCDBEntry *mapTSDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsTimeSDD");
510 if(!entrySDD || !entry2SDD || !drSpSDD || !ddlMapSDD || !mapTSDD ){
511 AliFatal("SDD Calibration object retrieval failed! ");
517 TObjArray *calSDD = (TObjArray *)entrySDD->GetObject();
518 if(!cacheStatus)entrySDD->SetObject(NULL);
519 entrySDD->SetOwner(kTRUE);
521 AliITSresponseSDD *pSDD = (AliITSresponseSDD*)entry2SDD->GetObject();
522 if(!cacheStatus)entry2SDD->SetObject(NULL);
523 entry2SDD->SetOwner(kTRUE);
525 TObjArray *drSp = (TObjArray *)drSpSDD->GetObject();
526 if(!cacheStatus)drSpSDD->SetObject(NULL);
527 drSpSDD->SetOwner(kTRUE);
529 AliITSDDLModuleMapSDD *ddlsdd=(AliITSDDLModuleMapSDD*)ddlMapSDD->GetObject();
530 if(!cacheStatus)ddlMapSDD->SetObject(NULL);
531 ddlMapSDD->SetOwner(kTRUE);
533 // TObjArray *mapAn = (TObjArray *)mapASDD->GetObject();
534 // if(!cacheStatus)mapASDD->SetObject(NULL);
535 // mapASDD->SetOwner(kTRUE);
537 TObjArray *mapT = (TObjArray *)mapTSDD->GetObject();
538 if(!cacheStatus)mapTSDD->SetObject(NULL);
539 mapTSDD->SetOwner(kTRUE);
542 // DB entries are deleted. In this way metadeta objects are deleted as well
552 if ((!pSDD)||(!calSDD) || (!drSp) || (!ddlsdd) || (!mapT) ){
553 AliWarning("Can not get SDD calibration from calibration database !");
557 fNMod[1] = calSDD->GetEntries();
561 AliITSCalibration* cal;
564 Bool_t oldMapFormat=kFALSE;
565 TObject* objmap=(TObject*)mapT->At(0);
566 TString cname(objmap->ClassName());
567 if(cname.CompareTo("AliITSMapSDD")==0){
569 AliInfo("SDD Maps converted to new format");
571 for(Int_t iddl=0; iddl<AliITSDDLModuleMapSDD::GetNDDLs(); iddl++){
572 for(Int_t icar=0; icar<AliITSDDLModuleMapSDD::GetNModPerDDL();icar++){
573 Int_t iMod=fDDLMapSDD->GetModuleNumber(iddl,icar);
574 if(iMod==-1) continue;
575 Int_t i=iMod - fgkDefaultNModulesSPD;
576 cal = (AliITSCalibration*) calSDD->At(i);
579 for(Int_t iAnode=0;iAnode< ((AliITSCalibrationSDD*)cal)->NOfAnodes(); iAnode++){
580 if(((AliITSCalibrationSDD*)cal)->IsBadChannel(iAnode)) continue;
581 avegain+= ((AliITSCalibrationSDD*)cal)->GetChannelGain(iAnode);
584 AliITSDriftSpeedArraySDD* arr0 = (AliITSDriftSpeedArraySDD*) drSp->At(i0);
585 AliITSDriftSpeedArraySDD* arr1 = (AliITSDriftSpeedArraySDD*) drSp->At(i1);
587 AliITSCorrMapSDD* mt0 = 0;
588 AliITSCorrMapSDD* mt1 = 0;
590 AliITSMapSDD* oldmap0=(AliITSMapSDD*)mapT->At(i0);
591 AliITSMapSDD* oldmap1=(AliITSMapSDD*)mapT->At(i1);
592 mt0=oldmap0->ConvertToNewFormat();
593 mt1=oldmap1->ConvertToNewFormat();
595 mt0=(AliITSCorrMapSDD*)mapT->At(i0);
596 mt1=(AliITSCorrMapSDD*)mapT->At(i1);
598 cal->SetDriftSpeed(0,arr0);
599 cal->SetDriftSpeed(1,arr1);
602 SetCalibrationModel(iMod, cal);
605 if(nGdAnodes) fAveGainSDD=avegain/nGdAnodes;
610 //______________________________________________________________________
611 Bool_t AliITSDetTypeRec::GetCalibrationSSD(Bool_t cacheStatus) {
612 // Get SSD calibration objects from OCDB
613 // AliCDBEntry *entrySSD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSSD");
615 AliCDBEntry *entryNoiseSSD = AliCDBManager::Instance()->Get("ITS/Calib/NoiseSSD");
616 AliCDBEntry *entryGainSSD = AliCDBManager::Instance()->Get("ITS/Calib/GainSSD");
617 AliCDBEntry *entryBadChannelsSSD = AliCDBManager::Instance()->Get("ITS/Calib/BadChannelsSSD");
619 if(!entryNoiseSSD || !entryGainSSD || !entryBadChannelsSSD){
620 AliFatal("SSD Calibration object retrieval failed! ");
624 TObject *emptyssd = 0; TString ssdobjectname;
625 AliITSNoiseSSDv2 *noiseSSD = NULL;
626 emptyssd = (TObject *)entryNoiseSSD->GetObject();
627 ssdobjectname = emptyssd->GetName();
628 if(ssdobjectname=="TObjArray") {
629 TObjArray *noiseSSDOld = (TObjArray *)entryNoiseSSD->GetObject();
630 noiseSSD = new AliITSNoiseSSDv2();
631 ReadOldSSDNoise(noiseSSDOld, noiseSSD);
633 else if(ssdobjectname=="AliITSNoiseSSDv2")
634 noiseSSD = (AliITSNoiseSSDv2 *)entryNoiseSSD->GetObject();
635 if(!cacheStatus)entryNoiseSSD->SetObject(NULL);
636 entryNoiseSSD->SetOwner(kTRUE);
638 AliITSGainSSDv2 *gainSSD = NULL;;
639 emptyssd = (TObject *)entryGainSSD->GetObject();
640 ssdobjectname = emptyssd->GetName();
641 if(ssdobjectname=="Gain") {
642 TObjArray *gainSSDOld = (TObjArray *)entryGainSSD->GetObject();
643 gainSSD = new AliITSGainSSDv2();
644 ReadOldSSDGain(gainSSDOld, gainSSD);
646 else if(ssdobjectname=="AliITSGainSSDv2")
647 gainSSD = (AliITSGainSSDv2 *)entryGainSSD->GetObject();
648 if(!cacheStatus)entryGainSSD->SetObject(NULL);
649 entryGainSSD->SetOwner(kTRUE);
651 AliITSBadChannelsSSDv2 *badChannelsSSD = NULL;
652 emptyssd = (TObject *)entryBadChannelsSSD->GetObject();
653 ssdobjectname = emptyssd->GetName();
654 if(ssdobjectname=="TObjArray") {
655 TObjArray *badChannelsSSDOld = (TObjArray *)entryBadChannelsSSD->GetObject();
656 badChannelsSSD = new AliITSBadChannelsSSDv2();
657 ReadOldSSDBadChannels(badChannelsSSDOld, badChannelsSSD);
659 else if(ssdobjectname=="AliITSBadChannelsSSDv2")
660 badChannelsSSD = (AliITSBadChannelsSSDv2*)entryBadChannelsSSD->GetObject();
661 if(!cacheStatus)entryBadChannelsSSD->SetObject(NULL);
662 entryBadChannelsSSD->SetOwner(kTRUE);
664 // DB entries are deleted. In this way metadeta objects are deleted as well
666 delete entryNoiseSSD;
668 delete entryBadChannelsSSD;
671 if ((!noiseSSD)|| (!gainSSD)|| (!badChannelsSSD)) {
672 AliWarning("Can not get SSD calibration from calibration database !");
676 fSSDCalibration->SetNoise(noiseSSD);
677 fSSDCalibration->SetGain(gainSSD);
678 fSSDCalibration->SetBadChannels(badChannelsSSD);
679 //fSSDCalibration->FillBadChipMap();
684 //________________________________________________________________
685 void AliITSDetTypeRec::SetDefaultClusterFindersV2(Bool_t rawdata){
687 //Set defaults for cluster finder V2
690 Warning("SetDefaults","Null pointer to AliITSgeom !");
694 AliITSClusterFinder *clf;
696 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
699 if(!GetReconstructionModel(dettype)){
700 clf = new AliITSClusterFinderV2SPD(this);
702 if(!rawdata) clf->SetDigits(DigitsAddress(0));
703 SetReconstructionModel(dettype,clf);
709 if(!GetReconstructionModel(dettype)){
710 clf = new AliITSClusterFinderV2SDD(this);
712 if(!rawdata) clf->SetDigits(DigitsAddress(1));
713 SetReconstructionModel(dettype,clf);
720 if(!GetReconstructionModel(dettype)){
721 clf = new AliITSClusterFinderV2SSD(this);
723 if(!rawdata) clf->SetDigits(DigitsAddress(2));
724 SetReconstructionModel(dettype,clf);
731 //______________________________________________________________________
732 void AliITSDetTypeRec::MakeBranch(TTree* tree, Option_t* option){
734 //Creates branches for clusters and recpoints
735 Bool_t cR = (strstr(option,"R")!=0);
736 Bool_t cRF = (strstr(option,"RF")!=0);
740 if(cR) MakeBranchR(tree);
741 if(cRF) MakeBranchRF(tree);
745 //___________________________________________________________________
746 void AliITSDetTypeRec::ResetDigits(){
747 // Reset number of digits and the digits array for the ITS detector.
750 for(Int_t i=0;i<fgkNdettypes;i++){
754 //___________________________________________________________________
755 void AliITSDetTypeRec::ResetDigits(Int_t branch){
756 // Reset number of digits and the digits array for this branch.
758 if(fDigits->At(branch)) ((TClonesArray*)fDigits->At(branch))->Clear();
761 //__________________________________________________________________
762 void AliITSDetTypeRec::MakeBranchR(TTree *treeR, Option_t *opt){
764 //Creates tree branches for recpoints
766 // cont char *file File name where RecPoints branch is to be written
767 // to. If blank it write the SDigits to the same
768 // file in which the Hits were found.
773 // only one branch for rec points for all detector types
774 Bool_t oFast= (strstr(opt,"Fast")!=0);
776 Char_t detname[10] = "ITS";
780 sprintf(branchname,"%sRecPointsF",detname);
782 sprintf(branchname,"%sRecPoints",detname);
785 if(!fRecPoints)fRecPoints = new TClonesArray("AliITSRecPoint",1000);
787 MakeBranchInTree(treeR,branchname,0,&fRecPoints,buffsz,99);
789 //______________________________________________________________________
790 void AliITSDetTypeRec::SetTreeAddressR(TTree* const treeR){
791 // Set branch address for the Reconstructed points Trees.
793 // TTree *treeR Tree containing the RecPoints.
799 Char_t namedet[10]="ITS";
802 if(fRecPoints==0x0) fRecPoints = new TClonesArray("AliITSRecPoint",1000);
804 sprintf(branchname,"%sRecPoints",namedet);
805 branch = treeR->GetBranch(branchname);
807 branch->SetAddress(&fRecPoints);
810 sprintf(branchname,"%sRecPointsF",namedet);
811 branch = treeR->GetBranch(branchname);
813 branch->SetAddress(&fRecPoints);
817 //____________________________________________________________________
818 void AliITSDetTypeRec::AddRecPoint(const AliITSRecPoint &r){
819 // Add a reconstructed space point to the list
821 // const AliITSRecPoint &r RecPoint class to be added to the tree
822 // of reconstructed points TreeR.
828 TClonesArray &lrecp = *fRecPoints;
829 new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
832 //______________________________________________________________________
833 void AliITSDetTypeRec::DigitsToRecPoints(TTree *treeD,TTree *treeR,Int_t lastentry,Option_t *opt, Int_t optCluFind){
834 // cluster finding and reconstruction of space points
835 // the condition below will disappear when the geom class will be
836 // initialized for all versions - for the moment it is only for v5 !
837 // 7 is the SDD beam test version
839 // TTree *treeD Digits tree
840 // TTree *treeR Clusters tree
841 // Int_t lastentry Offset for module when not all of the modules
843 // Option_t *opt String indicating which ITS sub-detectors should
844 // be processed. If ="All" then all of the ITS
845 // sub detectors are processed.
847 const char *all = strstr(opt,"All");
848 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
851 SetDefaultClusterFindersV2();
852 AliDebug(1,"V2 cluster finder has been selected \n");
854 SetDefaultClusterFindersV2();
855 AliInfo("Cluster Finder Option not implemented, V2 cluster finder will be used \n");
859 // Reset Fast-OR fired map
860 ResetFastOrFiredMap();
862 if (all || det[0]) { // SPD present
863 // Get the FO signals for this event
864 AliRunLoader* runLoader = AliRunLoader::Instance();
865 AliITSLoader* itsLoader = (AliITSLoader*) runLoader->GetLoader("ITSLoader");
867 AliError("ITS loader is NULL.");
870 fFOSignals = (AliITSFOSignalsSPD*)itsLoader->TreeD()->GetUserInfo()->FindObject("AliITSFOSignalsSPD");
871 if(!fFOSignals) AliError("FO signals not retrieved");
877 AliITSClusterFinder *rec = 0;
878 Int_t id,module,first=0;
879 for(module=0;module<GetITSgeom()->GetIndexMax();module++){
880 id = GetITSgeom()->GetModuleType(module);
881 if (!all && !det[id]) continue;
882 if(det[id]) first = GetITSgeom()->GetStartDet(id);
883 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
884 TClonesArray *itsDigits = DigitsAddress(id);
886 AliFatal("The reconstruction class was not instanciated!");
887 ResetDigits(); // MvL: Not sure we neeed this when rereading anyways
889 treeD->GetEvent(lastentry+module);
892 treeD->GetEvent(lastentry+(module-first));
894 Int_t ndigits = itsDigits->GetEntriesFast();
895 if (ndigits>0 || id==0) { // for SPD we always want to call FindRawClusters (to process FO signals)
896 rec->SetDetTypeRec(this);
897 rec->SetDigits(DigitsAddress(id));
898 // rec->SetClusters(ClustersAddress(id));
899 rec->FindRawClusters(module);
905 // Remove PIT in-active chips from Fast-OR fired map
906 if (all || det[0]) { // SPD present
907 RemoveFastOrFiredInActive();
908 // here removing bits which have no associated clusters
909 RemoveFastOrFiredFromDead(GetFiredChipMap(treeR));
912 AliITSRecPointContainer* rpcont = AliITSRecPointContainer::Instance();
914 nClu[0]=rpcont->GetNClustersInLayer(1,treeR);
915 for(Int_t iLay=2; iLay<=6; iLay++) nClu[iLay-1]=rpcont->GetNClustersInLayerFast(iLay);
916 AliInfo(Form("Number of RecPoints in ITS Layers = %d %d %d %d %d %d",
917 nClu[0],nClu[1],nClu[2],nClu[3],nClu[4],nClu[5]));
919 //______________________________________________________________________
920 void AliITSDetTypeRec::DigitsToRecPoints(AliRawReader* rawReader,TTree *treeR,Option_t *opt){
921 // cluster finding and reconstruction of space points
922 // the condition below will disappear when the geom class will be
923 // initialized for all versions - for the moment it is only for v5 !
924 // 7 is the SDD beam test version
926 // AliRawReader *rawReader Pointer to the raw-data reader
927 // TTree *treeR Clusters tree
932 const char *all = strstr(opt,"All");
933 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
937 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
939 TClonesArray* array = rpc->UncheckedGetClusters(0);
940 TBranch *branch = treeR->Branch("ITSRecPoints",&array);
941 DigitsToRecPoints(rawReader,opt);
944 for(Int_t iModule=0;iModule<GetITSgeom()->GetIndexMax();iModule++){
945 id = GetITSgeom()->GetModuleType(iModule);
946 if (!all && !det[id]) continue;
947 array = rpc->UncheckedGetClusters(iModule);
949 AliDebug(1,Form("data for module %d missing!",iModule));
951 branch->SetAddress(&array);
953 nClusters+=array->GetEntriesFast();
958 AliITSRecPointContainer* rpcont = AliITSRecPointContainer::Instance();
960 nClu[0]=rpcont->GetNClustersInLayer(1,treeR);
961 for(Int_t iLay=2; iLay<=6; iLay++) nClu[iLay-1]=rpcont->GetNClustersInLayerFast(iLay);
962 AliInfo(Form("Number of RecPoints in ITS Layers = %d %d %d %d %d %d, Total = %d",
963 nClu[0],nClu[1],nClu[2],nClu[3],nClu[4],nClu[5],nClusters));
965 //______________________________________________________________________
966 void AliITSDetTypeRec::DigitsToRecPoints(AliRawReader* rawReader,Option_t *opt){
967 // cluster finding and reconstruction of space points
968 // the condition below will disappear when the geom class will be
969 // initialized for all versions - for the moment it is only for v5 !
970 // 7 is the SDD beam test version
972 // AliRawReader *rawReader Pointer to the raw-data reader
977 const char *all = strstr(opt,"All");
978 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
981 // Reset Fast-OR fired map
982 ResetFastOrFiredMap();
984 AliITSClusterFinder *rec = 0;
988 if (!all && !det[id]) continue;
989 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
991 AliFatal("The reconstruction class was not instantiated");
992 rec->SetDetTypeRec(this);
993 rec->RawdataToClusters(rawReader);
996 // Remove PIT in-active chips from Fast-OR fired map
997 if (all || det[0]) { // SPD present
998 RemoveFastOrFiredInActive();
999 // here removing bits which have no associated clusters
1000 RemoveFastOrFiredFromDead(GetFiredChipMap());
1004 //______________________________________________________________________
1005 void AliITSDetTypeRec::ReadOldSSDNoise(const TObjArray *array,
1006 AliITSNoiseSSDv2 *noiseSSD) {
1007 //Reads the old SSD calibration object and converts it to the new format
1008 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
1009 const Int_t fgkSSDPSIDESTRIPSPERMODULE = 768;
1011 Int_t gNMod = array->GetEntries();
1012 cout<<"Converting old calibration object for noise..."<<endl;
1015 Double_t noise = 0.0;
1016 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1017 AliITSNoiseSSD *noiseModule = (AliITSNoiseSSD*) (array->At(iModule));
1018 for(Int_t iStrip = 0; iStrip < fgkSSDSTRIPSPERMODULE; iStrip++) {
1019 noise = (iStrip < fgkSSDPSIDESTRIPSPERMODULE) ? noiseModule->GetNoiseP(iStrip) : noiseModule->GetNoiseN(1535 - iStrip);
1020 if(iStrip < fgkSSDPSIDESTRIPSPERMODULE)
1021 noiseSSD->AddNoiseP(iModule,iStrip,noise);
1022 if(iStrip >= fgkSSDPSIDESTRIPSPERMODULE)
1023 noiseSSD->AddNoiseN(iModule,1535 - iStrip,noise);
1025 }//loop over modules
1028 //______________________________________________________________________
1029 void AliITSDetTypeRec::ReadOldSSDBadChannels(const TObjArray *array,
1030 AliITSBadChannelsSSDv2 *badChannelsSSD) {
1031 //Reads the old SSD calibration object and converts it to the new format
1032 Int_t gNMod = array->GetEntries();
1033 cout<<"Converting old calibration object for bad channels..."<<endl;
1034 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1035 //for (Int_t iModule = 0; iModule < 1; iModule++) {
1036 AliITSBadChannelsSSD *bad = (AliITSBadChannelsSSD*) (array->At(iModule));
1037 TArrayI arrayPSide = bad->GetBadPChannelsList();
1038 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1039 badChannelsSSD->AddBadChannelP(iModule,
1041 (Char_t)arrayPSide.At(iPCounter));
1043 TArrayI arrayNSide = bad->GetBadNChannelsList();
1044 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1045 badChannelsSSD->AddBadChannelN(iModule,
1047 (Char_t)arrayNSide.At(iNCounter));
1049 }//loop over modules
1052 //______________________________________________________________________
1053 void AliITSDetTypeRec::ReadOldSSDGain(const TObjArray *array,
1054 AliITSGainSSDv2 *gainSSD) {
1055 //Reads the old SSD calibration object and converts it to the new format
1057 Int_t gNMod = array->GetEntries();
1058 cout<<"Converting old calibration object for gain..."<<endl;
1061 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1062 AliITSGainSSD *gainModule = (AliITSGainSSD*) (array->At(iModule));
1063 TArrayF arrayPSide = gainModule->GetGainP();
1064 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1065 gainSSD->AddGainP(iModule,
1067 arrayPSide.At(iPCounter));
1068 TArrayF arrayNSide = gainModule->GetGainN();
1069 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1070 gainSSD->AddGainN(iModule,
1072 arrayNSide.At(iNCounter));
1073 }//loop over modules
1075 //______________________________________________________________________
1076 void AliITSDetTypeRec::RemoveFastOrFiredInActive() {
1077 // Removes the chips that were in-active in the pixel trigger (from fast-or fired map)
1079 if (fTriggerConditions==NULL) {
1080 AliError("Pixel trigger conditions are missing.");
1086 while (fTriggerConditions->GetNextInActiveChip(eq,hs,chip)) {
1087 UInt_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1088 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1091 //______________________________________________________________________
1092 TBits AliITSDetTypeRec::GetFiredChipMap() const {
1095 // TBits of the fired chips
1098 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
1100 TBits isfiredchip(1200);
1102 AliITSsegmentationSPD *segSPD = (AliITSsegmentationSPD*)GetSegmentationModel(0);
1104 AliError("no segmentation model for SPD available, the fired chip map is empty. Exiting");
1109 for(Int_t imod =0; imod < fgkDefaultNModulesSPD; imod++){
1110 TClonesArray *array = rpc->UncheckedGetClusters(imod);
1111 if(!array) continue;
1112 Int_t nCluster = array->GetEntriesFast();
1115 AliITSRecPoint* cluster = (AliITSRecPoint*)array->UncheckedAt(nCluster);
1116 if (cluster->GetLayer()>1)continue;
1117 Float_t local[3]={-1,-1};
1118 local[1]=cluster->GetDetLocalX();
1119 local[0]=cluster->GetDetLocalZ();
1121 Int_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(imod);
1122 Int_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(imod);
1124 segSPD->LocalToDet(0.5,local[0],row,col);
1125 Int_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(imod,col);
1126 Int_t chipkey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1127 isfiredchip.SetBitNumber(chipkey,kTRUE);
1135 //______________________________________________________________________
1136 TBits AliITSDetTypeRec::GetFiredChipMap(TTree *treeR) const{
1138 // TBits of the fired chips
1140 TBits isfiredchip(1200);
1143 AliError("no treeR. fired chip map stays empty. Exiting.");
1147 AliITSRecPointContainer* rpcont=AliITSRecPointContainer::Instance();
1148 TClonesArray *recpoints = rpcont->FetchClusters(0,treeR);
1149 if(!rpcont->GetStatusOK() || !rpcont->IsSPDActive()){
1150 AliError("no clusters. fired chip map stays empty. Exiting.");
1154 AliITSsegmentationSPD *segSPD = (AliITSsegmentationSPD*)GetSegmentationModel(0);
1156 for(Int_t imod =0; imod < fgkDefaultNModulesSPD; imod++){
1157 recpoints = rpcont->UncheckedGetClusters(imod);
1158 Int_t nCluster = recpoints->GetEntriesFast();
1160 // loop over clusters
1162 AliITSRecPoint* cluster = (AliITSRecPoint*)recpoints->UncheckedAt(nCluster);
1163 if (cluster->GetLayer()>1)continue;
1164 Float_t local[3]={-1,-1};
1165 local[1]=cluster->GetDetLocalX();
1166 local[0]=cluster->GetDetLocalZ();
1168 Int_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(imod);
1169 Int_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(imod);
1171 segSPD->LocalToDet(0.5,local[0],row,col);
1172 Int_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(imod,col);
1173 Int_t chipkey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1174 isfiredchip.SetBitNumber(chipkey,kTRUE);
1180 //______________________________________________________________________
1181 void AliITSDetTypeRec::RemoveFastOrFiredFromDead(TBits firedchipmap){
1183 // resetting of the fast-or bit on cluster basis.
1184 // fast-or bits can be remnant from SPD ideal simulation (no dead channels)
1187 for(Int_t chipKey=0; chipKey<1200; chipKey++){
1188 // FO masked chips have been previously removed
1189 if(!fFastOrFiredMap.TestBitNumber(chipKey)) continue;
1190 if(!firedchipmap.TestBitNumber(chipKey)) {
1191 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1192 AliDebug(2,Form("removing bit in key %i \n ",chipKey));
1197 //______________________________________________________________________
1198 void AliITSDetTypeRec::SetFastOrFiredMapOnline(UInt_t eq, UInt_t hs, UInt_t chip) {
1199 // Set fast-or fired map for this chip
1200 Int_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1201 return SetFastOrFiredMap(chipKey);