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(),
80 fTriggerConditions(0),
90 fFastOrFiredMap(1200){
91 // Standard Constructor
99 fReconstruction = new TObjArray(fgkNdettypes);
100 fDigits = new TObjArray(fgkNdettypes);
101 for(Int_t i=0; i<3; i++){
104 fSSDCalibration=new AliITSCalibrationSSD();
105 fNMod = new Int_t [fgkNdettypes];
106 fNMod[0] = fgkDefaultNModulesSPD;
107 fNMod[1] = fgkDefaultNModulesSDD;
108 fNMod[2] = fgkDefaultNModulesSSD;
109 fRecPoints = new TClonesArray("AliITSRecPoint",3000);
115 //______________________________________________________________________
116 AliITSDetTypeRec::AliITSDetTypeRec(const AliITSDetTypeRec & rec):TObject(rec),
118 fITSgeom(rec.fITSgeom),
119 fReconstruction(rec.fReconstruction),
120 fSegmentation(rec.fSegmentation),
121 fCalibration(rec.fCalibration),
122 fSSDCalibration(rec.fSSDCalibration),
123 fSPDDead(rec.fSPDDead),
124 fSPDSparseDead(rec.fSPDSparseDead),
125 fTriggerConditions(rec.fTriggerConditions),
126 fDigits(rec.fDigits),
127 fFOSignals(rec.fFOSignals),
128 fDDLMapSDD(rec.fDDLMapSDD),
129 fRespSDD(rec.fRespSDD),
130 fAveGainSDD(rec.fAveGainSDD),
131 fRecPoints(rec.fRecPoints),
132 fNRecPoints(rec.fNRecPoints),
133 fFirstcall(rec.fFirstcall),
134 fLoadOnlySPDCalib(rec.fLoadOnlySPDCalib),
135 fFastOrFiredMap(rec.fFastOrFiredMap){
138 for(Int_t i=0; i<3; i++){
139 fkDigClassName[i]=rec.fkDigClassName[i]; // NB only copies Char_t*, so not so safe, but this code should never be reached anyways
142 //______________________________________________________________________
143 AliITSDetTypeRec& AliITSDetTypeRec::operator=(const AliITSDetTypeRec& source){
144 // Assignment operator.
145 this->~AliITSDetTypeRec();
146 new(this) AliITSDetTypeRec(source);
151 //_____________________________________________________________________
152 AliITSDetTypeRec::~AliITSDetTypeRec(){
156 fReconstruction->Delete();
157 delete fReconstruction;
161 fSegmentation->Delete();
162 delete fSegmentation;
166 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
167 fCalibration->Delete();
170 if(fRespSDD) delete fRespSDD;
171 if(fDDLMapSDD) delete fDDLMapSDD;
175 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
176 delete fSSDCalibration;
177 fSSDCalibration = NULL;
181 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
188 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
189 fSPDSparseDead->Delete();
190 delete fSPDSparseDead;
194 if(fTriggerConditions){
195 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
196 fTriggerConditions->Delete();
197 delete fTriggerConditions;
198 fTriggerConditions = 0;
207 fRecPoints->Delete();
213 if (fITSgeom) delete fITSgeom;
217 //___________________________________________________________________
218 void AliITSDetTypeRec::SetReconstructionModel(Int_t dettype,AliITSClusterFinder *clf){
220 //Set reconstruction model for detector type
222 if(fReconstruction==0) fReconstruction = new TObjArray(fgkNdettypes);
223 if(fReconstruction->At(dettype)!=0) delete fReconstruction->At(dettype);
224 fReconstruction->AddAt(clf,dettype);
226 //______________________________________________________________________
227 AliITSClusterFinder* AliITSDetTypeRec::GetReconstructionModel(Int_t dettype) const{
229 //Get reconstruction model for detector type
230 if(fReconstruction==0) {
231 Warning("GetReconstructionModel","fReconstruction is 0!");
234 return (AliITSClusterFinder*)fReconstruction->At(dettype);
237 //______________________________________________________________________
238 void AliITSDetTypeRec::SetSegmentationModel(Int_t dettype,AliITSsegmentation *seg){
240 //Set segmentation model for detector type
242 if(fSegmentation==0) fSegmentation = new TObjArray(fgkNdettypes);
243 if(fSegmentation->At(dettype)!=0) delete fSegmentation->At(dettype);
244 fSegmentation->AddAt(seg,dettype);
247 //______________________________________________________________________
248 AliITSsegmentation* AliITSDetTypeRec::GetSegmentationModel(Int_t dettype) const {
250 //Get segmentation model for detector type
252 if(fSegmentation==0) {
253 Warning("GetSegmentationModel","fSegmentation is 0!");
256 return (AliITSsegmentation*)fSegmentation->At(dettype);
259 //_______________________________________________________________________
260 void AliITSDetTypeRec::SetCalibrationModel(Int_t iMod, AliITSCalibration *cal){
262 //Set calibration (response) for the module iMod of type dettype
263 if (fCalibration==0) {
264 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
265 fCalibration->SetOwner(kTRUE);
266 fCalibration->Clear();
269 if (fCalibration->At(iMod) != 0)
270 delete (AliITSCalibration*) fCalibration->At(iMod);
271 fCalibration->AddAt(cal,iMod);
274 //_______________________________________________________________________
275 void AliITSDetTypeRec::SetSPDDeadModel(Int_t iMod, AliITSCalibration *cal){
277 //Set dead pixel info for the SPD module iMod
279 fSPDDead = new TObjArray(fgkDefaultNModulesSPD);
280 fSPDDead->SetOwner(kTRUE);
284 if (fSPDDead->At(iMod) != 0)
285 delete (AliITSCalibration*) fSPDDead->At(iMod);
286 fSPDDead->AddAt(cal,iMod);
288 //_______________________________________________________________________
289 void AliITSDetTypeRec::SetSPDSparseDeadModel(Int_t iMod, AliITSCalibration *cal){
291 //Set dead pixel info for the SPD ACTIVE module iMod
292 if (fSPDSparseDead==0) {
293 fSPDSparseDead = new TObjArray(fgkDefaultNModulesSPD);
294 fSPDSparseDead->SetOwner(kTRUE);
295 fSPDSparseDead->Clear();
298 if (fSPDSparseDead->At(iMod) != 0)
299 delete (AliITSCalibration*) fSPDSparseDead->At(iMod);
300 fSPDSparseDead->AddAt(cal,iMod);
302 //_______________________________________________________________________
303 AliITSCalibration* AliITSDetTypeRec::GetCalibrationModel(Int_t iMod) const {
305 //Get calibration model for module type
307 if(fCalibration==0) {
308 Warning("GetalibrationModel","fCalibration is 0!");
312 if(iMod<fgkDefaultNModulesSPD+fgkDefaultNModulesSDD){
313 return (AliITSCalibration*)fCalibration->At(iMod);
315 Int_t i=iMod-(fgkDefaultNModulesSPD+fgkDefaultNModulesSDD);
316 fSSDCalibration->SetModule(i);
317 return (AliITSCalibration*)fSSDCalibration;
321 //_______________________________________________________________________
322 AliITSCalibration* AliITSDetTypeRec::GetSPDDeadModel(Int_t iMod) const {
324 //Get SPD dead for module iMod
327 AliWarning("fSPDDead is 0!");
330 return (AliITSCalibration*)fSPDDead->At(iMod);
332 //_______________________________________________________________________
333 AliITSCalibration* AliITSDetTypeRec::GetSPDSparseDeadModel(Int_t iMod) const {
335 //Get SPD dead for module iMod
337 if(fSPDSparseDead==0) {
338 AliWarning("fSPDSparseDead is 0!");
341 return (AliITSCalibration*)fSPDSparseDead->At(iMod);
343 //_______________________________________________________________________
344 AliITSTriggerConditions* AliITSDetTypeRec::GetTriggerConditions() const {
345 //Get Pixel Trigger Conditions
346 if (fTriggerConditions==0) {
347 AliWarning("fTriggerConditions is 0!");
349 return fTriggerConditions;
351 //______________________________________________________________________
352 void AliITSDetTypeRec::SetTreeAddressD(TTree* const treeD){
353 // Set branch address for the tree of digits.
355 const char *det[4] = {"SPD","SDD","SSD","ITS"};
357 const Char_t* digclass;
362 if (fDigits == 0x0) {
363 fDigits = new TObjArray(fgkNdettypes);
368 for (i=0; i<fgkNdettypes; i++) {
369 digclass = GetDigitClassName(i);
370 fDigits->AddAt(new TClonesArray(digclass,1000),i);
371 if (fgkNdettypes==3) snprintf(branchname,29,"%sDigits%s",det[3],det[i]);
372 else snprintf(branchname,29,"%sDigits%d",det[3],i+1);
373 branch = treeD->GetBranch(branchname);
374 if (branch) branch->SetAddress(&((*fDigits)[i]));
379 //_______________________________________________________________________
380 TBranch* AliITSDetTypeRec::MakeBranchInTree(TTree* const tree,
381 const char* name, const char *classname,
382 void* address,Int_t size,Int_t splitlevel)
385 // Makes branch in given tree and diverts them to a separate file
391 Error("MakeBranchInTree","Making Branch %s Tree is NULL",name);
394 TBranch *branch = tree->GetBranch(name);
399 branch = tree->Branch(name,classname,address,size,splitlevel);
402 branch = tree->Bronch(name, "TClonesArray", address, size, splitlevel);
408 //____________________________________________________________________
409 void AliITSDetTypeRec::SetDefaults(){
411 //Set defaults for segmentation and response
414 Warning("SetDefaults","null pointer to AliITSgeomGeom !");
418 AliITSsegmentation* seg;
419 if(!GetCalibration()) {AliFatal("Exit");exit(0);}
421 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
423 seg = new AliITSsegmentationSPD();
424 SetSegmentationModel(dettype,seg);
425 SetDigitClassName(dettype,"AliITSdigitSPD");
428 seg = new AliITSsegmentationSDD();
429 if(fLoadOnlySPDCalib==kFALSE){
430 AliITSCalibrationSDD* cal=(AliITSCalibrationSDD*)GetCalibrationModel(fgkDefaultNModulesSPD+1);
431 if(cal->IsAMAt20MHz()){
432 seg->SetPadSize(seg->Dpz(0),20.);
433 seg->SetNPads(seg->Npz()/2,128);
436 SetSegmentationModel(dettype,seg);
437 SetDigitClassName(dettype,"AliITSdigitSDD");
440 AliITSsegmentationSSD* seg2 = new AliITSsegmentationSSD();
441 SetSegmentationModel(dettype,seg2);
442 SetDigitClassName(dettype,"AliITSdigitSSD");
446 //______________________________________________________________________
447 Bool_t AliITSDetTypeRec::GetCalibration() {
448 // Get Default calibration if a storage is not defined.
451 AliITSCalibration* cal = GetCalibrationModel(0);
457 // SetRunNumber((Int_t)AliCDBManager::Instance()->GetRun());
458 // Int_t run=GetRunNumber();
460 Bool_t cacheStatus = AliCDBManager::Instance()->GetCacheFlag();
461 if (fCalibration==0) {
462 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
463 fCalibration->SetOwner(!cacheStatus);
464 fCalibration->Clear();
467 Bool_t retCode=GetCalibrationSPD(cacheStatus);
468 if(retCode==kFALSE) return kFALSE;
470 if(fLoadOnlySPDCalib==kFALSE){
471 retCode=GetCalibrationSDD(cacheStatus);
472 if(retCode==kFALSE) return kFALSE;
473 retCode=GetCalibrationSSD(cacheStatus);
474 if(retCode==kFALSE) return kFALSE;
477 AliInfo(Form("%i SPD, %i SDD and %i SSD in calibration database",
478 fNMod[0], fNMod[1], fNMod[2]));
481 //______________________________________________________________________
482 Bool_t AliITSDetTypeRec::GetCalibrationSPD(Bool_t cacheStatus) {
483 // Get SPD calibration objects from OCDB
484 // dead pixel are not used for local reconstruction
487 AliCDBEntry *noisySPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDNoisy");
488 AliCDBEntry *deadSPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDDead");
489 AliCDBEntry *deadSparseSPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDSparseDead");
490 AliCDBEntry *pitCond = AliCDBManager::Instance()->Get("TRIGGER/SPD/PITConditions");
491 if(!noisySPD || !deadSPD || !pitCond ){
492 AliFatal("SPD Calibration object retrieval failed! ");
496 TObjArray *calNoisySPD = (TObjArray*) noisySPD->GetObject();
497 if (!cacheStatus) noisySPD->SetObject(NULL);
498 noisySPD->SetOwner(kTRUE);
500 TObjArray *calDeadSPD = (TObjArray*) deadSPD->GetObject();
501 if (!cacheStatus) deadSPD->SetObject(NULL);
502 deadSPD->SetOwner(kTRUE);
504 TObjArray *calSparseDeadSPD = (TObjArray*) deadSparseSPD->GetObject();
505 if (!cacheStatus) deadSparseSPD->SetObject(NULL);
506 deadSparseSPD->SetOwner(kTRUE);
509 AliITSTriggerConditions *calPitCond = (AliITSTriggerConditions*) pitCond->GetObject();
510 if (!cacheStatus) pitCond->SetObject(NULL);
511 pitCond->SetOwner(kTRUE);
516 delete deadSparseSPD;
519 if ((!calNoisySPD) || (!calDeadSPD) || (!calSparseDeadSPD) || (!calPitCond)){
520 AliWarning("Can not get SPD calibration from calibration database !");
524 fNMod[0] = calNoisySPD->GetEntries();
526 AliITSCalibration* cal;
527 for (Int_t i=0; i<fNMod[0]; i++) {
528 cal = (AliITSCalibration*) calNoisySPD->At(i);
529 SetCalibrationModel(i, cal);
530 cal = (AliITSCalibration*) calDeadSPD->At(i);
531 SetSPDDeadModel(i, cal);
532 cal = (AliITSCalibration*) calSparseDeadSPD->At(i);
533 SetSPDSparseDeadModel(i, cal);
535 fTriggerConditions = calPitCond;
540 //______________________________________________________________________
541 Bool_t AliITSDetTypeRec::GetCalibrationSDD(Bool_t cacheStatus) {
542 // Get SDD calibration objects from OCDB
544 AliCDBEntry *entrySDD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSDD");
545 AliCDBEntry *entry2SDD = AliCDBManager::Instance()->Get("ITS/Calib/RespSDD");
546 AliCDBEntry *drSpSDD = AliCDBManager::Instance()->Get("ITS/Calib/DriftSpeedSDD");
547 AliCDBEntry *ddlMapSDD = AliCDBManager::Instance()->Get("ITS/Calib/DDLMapSDD");
548 // AliCDBEntry *mapASDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsAnodeSDD");
549 AliCDBEntry *mapTSDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsTimeSDD");
551 if(!entrySDD || !entry2SDD || !drSpSDD || !ddlMapSDD || !mapTSDD ){
552 AliFatal("SDD Calibration object retrieval failed! ");
558 TObjArray *calSDD = (TObjArray *)entrySDD->GetObject();
559 if(!cacheStatus)entrySDD->SetObject(NULL);
560 entrySDD->SetOwner(kTRUE);
562 AliITSresponseSDD *pSDD = (AliITSresponseSDD*)entry2SDD->GetObject();
563 if(!cacheStatus)entry2SDD->SetObject(NULL);
564 entry2SDD->SetOwner(kTRUE);
566 TObjArray *drSp = (TObjArray *)drSpSDD->GetObject();
567 if(!cacheStatus)drSpSDD->SetObject(NULL);
568 drSpSDD->SetOwner(kTRUE);
570 AliITSDDLModuleMapSDD *ddlsdd=(AliITSDDLModuleMapSDD*)ddlMapSDD->GetObject();
571 if(!cacheStatus)ddlMapSDD->SetObject(NULL);
572 ddlMapSDD->SetOwner(kTRUE);
574 // TObjArray *mapAn = (TObjArray *)mapASDD->GetObject();
575 // if(!cacheStatus)mapASDD->SetObject(NULL);
576 // mapASDD->SetOwner(kTRUE);
578 TObjArray *mapT = (TObjArray *)mapTSDD->GetObject();
579 if(!cacheStatus)mapTSDD->SetObject(NULL);
580 mapTSDD->SetOwner(kTRUE);
583 // DB entries are deleted. In this way metadeta objects are deleted as well
593 if ((!pSDD)||(!calSDD) || (!drSp) || (!ddlsdd) || (!mapT) ){
594 AliWarning("Can not get SDD calibration from calibration database !");
598 fNMod[1] = calSDD->GetEntries();
602 AliITSCalibration* cal;
605 Bool_t oldMapFormat=kFALSE;
606 TObject* objmap=(TObject*)mapT->At(0);
607 TString cname(objmap->ClassName());
608 if(cname.CompareTo("AliITSMapSDD")==0){
610 AliInfo("SDD Maps converted to new format");
612 for(Int_t iddl=0; iddl<AliITSDDLModuleMapSDD::GetNDDLs(); iddl++){
613 for(Int_t icar=0; icar<AliITSDDLModuleMapSDD::GetNModPerDDL();icar++){
614 Int_t iMod=fDDLMapSDD->GetModuleNumber(iddl,icar);
615 if(iMod==-1) continue;
616 Int_t i=iMod - fgkDefaultNModulesSPD;
617 cal = (AliITSCalibration*) calSDD->At(i);
620 for(Int_t iAnode=0;iAnode< ((AliITSCalibrationSDD*)cal)->NOfAnodes(); iAnode++){
621 if(((AliITSCalibrationSDD*)cal)->IsBadChannel(iAnode)) continue;
622 avegain+= ((AliITSCalibrationSDD*)cal)->GetChannelGain(iAnode);
625 AliITSDriftSpeedArraySDD* arr0 = (AliITSDriftSpeedArraySDD*) drSp->At(i0);
626 AliITSDriftSpeedArraySDD* arr1 = (AliITSDriftSpeedArraySDD*) drSp->At(i1);
628 AliITSCorrMapSDD* mt0 = 0;
629 AliITSCorrMapSDD* mt1 = 0;
631 AliITSMapSDD* oldmap0=(AliITSMapSDD*)mapT->At(i0);
632 AliITSMapSDD* oldmap1=(AliITSMapSDD*)mapT->At(i1);
633 mt0=oldmap0->ConvertToNewFormat();
634 mt1=oldmap1->ConvertToNewFormat();
636 mt0=(AliITSCorrMapSDD*)mapT->At(i0);
637 mt1=(AliITSCorrMapSDD*)mapT->At(i1);
639 cal->SetDriftSpeed(0,arr0);
640 cal->SetDriftSpeed(1,arr1);
643 SetCalibrationModel(iMod, cal);
646 if(nGdAnodes) fAveGainSDD=avegain/nGdAnodes;
651 //______________________________________________________________________
652 Bool_t AliITSDetTypeRec::GetCalibrationSSD(Bool_t cacheStatus) {
653 // Get SSD calibration objects from OCDB
654 // AliCDBEntry *entrySSD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSSD");
656 AliCDBEntry *entryNoiseSSD = AliCDBManager::Instance()->Get("ITS/Calib/NoiseSSD");
657 AliCDBEntry *entryGainSSD = AliCDBManager::Instance()->Get("ITS/Calib/GainSSD");
658 AliCDBEntry *entryBadChannelsSSD = AliCDBManager::Instance()->Get("ITS/Calib/BadChannelsSSD");
660 if(!entryNoiseSSD || !entryGainSSD || !entryBadChannelsSSD){
661 AliFatal("SSD Calibration object retrieval failed! ");
665 TObject *emptyssd = 0; TString ssdobjectname;
666 AliITSNoiseSSDv2 *noiseSSD = NULL;
667 emptyssd = (TObject *)entryNoiseSSD->GetObject();
668 ssdobjectname = emptyssd->GetName();
669 if(ssdobjectname=="TObjArray") {
670 TObjArray *noiseSSDOld = (TObjArray *)entryNoiseSSD->GetObject();
671 noiseSSD = new AliITSNoiseSSDv2();
672 ReadOldSSDNoise(noiseSSDOld, noiseSSD);
674 else if(ssdobjectname=="AliITSNoiseSSDv2")
675 noiseSSD = (AliITSNoiseSSDv2 *)entryNoiseSSD->GetObject();
676 if(!cacheStatus)entryNoiseSSD->SetObject(NULL);
677 entryNoiseSSD->SetOwner(kTRUE);
679 AliITSGainSSDv2 *gainSSD = NULL;;
680 emptyssd = (TObject *)entryGainSSD->GetObject();
681 ssdobjectname = emptyssd->GetName();
682 if(ssdobjectname=="Gain") {
683 TObjArray *gainSSDOld = (TObjArray *)entryGainSSD->GetObject();
684 gainSSD = new AliITSGainSSDv2();
685 ReadOldSSDGain(gainSSDOld, gainSSD);
687 else if(ssdobjectname=="AliITSGainSSDv2")
688 gainSSD = (AliITSGainSSDv2 *)entryGainSSD->GetObject();
689 if(!cacheStatus)entryGainSSD->SetObject(NULL);
690 entryGainSSD->SetOwner(kTRUE);
692 AliITSBadChannelsSSDv2 *badChannelsSSD = NULL;
693 emptyssd = (TObject *)entryBadChannelsSSD->GetObject();
694 ssdobjectname = emptyssd->GetName();
695 if(ssdobjectname=="TObjArray") {
696 TObjArray *badChannelsSSDOld = (TObjArray *)entryBadChannelsSSD->GetObject();
697 badChannelsSSD = new AliITSBadChannelsSSDv2();
698 ReadOldSSDBadChannels(badChannelsSSDOld, badChannelsSSD);
700 else if(ssdobjectname=="AliITSBadChannelsSSDv2")
701 badChannelsSSD = (AliITSBadChannelsSSDv2*)entryBadChannelsSSD->GetObject();
702 if(!cacheStatus)entryBadChannelsSSD->SetObject(NULL);
703 entryBadChannelsSSD->SetOwner(kTRUE);
705 // DB entries are deleted. In this way metadeta objects are deleted as well
707 delete entryNoiseSSD;
709 delete entryBadChannelsSSD;
712 if ((!noiseSSD)|| (!gainSSD)|| (!badChannelsSSD)) {
713 AliWarning("Can not get SSD calibration from calibration database !");
717 fSSDCalibration->SetNoise(noiseSSD);
718 fSSDCalibration->SetGain(gainSSD);
719 fSSDCalibration->SetBadChannels(badChannelsSSD);
720 //fSSDCalibration->FillBadChipMap();
725 //________________________________________________________________
726 void AliITSDetTypeRec::SetDefaultClusterFindersV2(Bool_t rawdata){
728 //Set defaults for cluster finder V2
731 Warning("SetDefaults","Null pointer to AliITSgeom !");
735 AliITSClusterFinder *clf;
737 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
740 if(!GetReconstructionModel(dettype)){
741 clf = new AliITSClusterFinderV2SPD(this);
743 if(!rawdata) clf->SetDigits(DigitsAddress(0));
744 SetReconstructionModel(dettype,clf);
750 if(!GetReconstructionModel(dettype)){
751 clf = new AliITSClusterFinderV2SDD(this);
753 if(!rawdata) clf->SetDigits(DigitsAddress(1));
754 SetReconstructionModel(dettype,clf);
761 if(!GetReconstructionModel(dettype)){
762 clf = new AliITSClusterFinderV2SSD(this);
764 if(!rawdata) clf->SetDigits(DigitsAddress(2));
765 SetReconstructionModel(dettype,clf);
772 //______________________________________________________________________
773 void AliITSDetTypeRec::MakeBranch(TTree* tree, Option_t* option){
775 //Creates branches for clusters and recpoints
776 Bool_t cR = (strstr(option,"R")!=0);
777 Bool_t cRF = (strstr(option,"RF")!=0);
781 if(cR) MakeBranchR(tree);
782 if(cRF) MakeBranchRF(tree);
786 //___________________________________________________________________
787 void AliITSDetTypeRec::ResetDigits(){
788 // Reset number of digits and the digits array for the ITS detector.
791 for(Int_t i=0;i<fgkNdettypes;i++){
795 //___________________________________________________________________
796 void AliITSDetTypeRec::ResetDigits(Int_t branch){
797 // Reset number of digits and the digits array for this branch.
799 if(fDigits->At(branch)) ((TClonesArray*)fDigits->At(branch))->Clear();
802 //__________________________________________________________________
803 void AliITSDetTypeRec::MakeBranchR(TTree *treeR, Option_t *opt){
805 //Creates tree branches for recpoints
807 // cont char *file File name where RecPoints branch is to be written
808 // to. If blank it write the SDigits to the same
809 // file in which the Hits were found.
814 // only one branch for rec points for all detector types
815 Bool_t oFast= (strstr(opt,"Fast")!=0);
817 Char_t detname[10] = "ITS";
821 snprintf(branchname,29,"%sRecPointsF",detname);
823 snprintf(branchname,29,"%sRecPoints",detname);
826 if(!fRecPoints)fRecPoints = new TClonesArray("AliITSRecPoint",1000);
828 MakeBranchInTree(treeR,branchname,0,&fRecPoints,buffsz,99);
830 //______________________________________________________________________
831 void AliITSDetTypeRec::SetTreeAddressR(TTree* const treeR){
832 // Set branch address for the Reconstructed points Trees.
834 // TTree *treeR Tree containing the RecPoints.
840 Char_t namedet[10]="ITS";
843 if(fRecPoints==0x0) fRecPoints = new TClonesArray("AliITSRecPoint",1000);
845 snprintf(branchname,29,"%sRecPoints",namedet);
846 branch = treeR->GetBranch(branchname);
848 branch->SetAddress(&fRecPoints);
851 snprintf(branchname,29,"%sRecPointsF",namedet);
852 branch = treeR->GetBranch(branchname);
854 branch->SetAddress(&fRecPoints);
858 //____________________________________________________________________
859 void AliITSDetTypeRec::AddRecPoint(const AliITSRecPoint &r){
860 // Add a reconstructed space point to the list
862 // const AliITSRecPoint &r RecPoint class to be added to the tree
863 // of reconstructed points TreeR.
869 TClonesArray &lrecp = *fRecPoints;
870 new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
873 //______________________________________________________________________
874 void AliITSDetTypeRec::DigitsToRecPoints(TTree *treeD,TTree *treeR,Int_t lastentry,Option_t *opt, Int_t optCluFind){
875 // cluster finding and reconstruction of space points
876 // the condition below will disappear when the geom class will be
877 // initialized for all versions - for the moment it is only for v5 !
878 // 7 is the SDD beam test version
880 // TTree *treeD Digits tree
881 // TTree *treeR Clusters tree
882 // Int_t lastentry Offset for module when not all of the modules
884 // Option_t *opt String indicating which ITS sub-detectors should
885 // be processed. If ="All" then all of the ITS
886 // sub detectors are processed.
888 const char *all = strstr(opt,"All");
889 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
892 SetDefaultClusterFindersV2();
893 AliDebug(1,"V2 cluster finder has been selected \n");
895 SetDefaultClusterFindersV2();
896 AliInfo("Cluster Finder Option not implemented, V2 cluster finder will be used \n");
900 // Reset Fast-OR fired map
901 ResetFastOrFiredMap();
903 if (all || det[0]) { // SPD present
904 // Get the FO signals for this event
905 AliRunLoader* runLoader = AliRunLoader::Instance();
906 AliITSLoader* itsLoader = (AliITSLoader*) runLoader->GetLoader("ITSLoader");
908 AliError("ITS loader is NULL.");
911 fFOSignals = (AliITSFOSignalsSPD*)itsLoader->TreeD()->GetUserInfo()->FindObject("AliITSFOSignalsSPD");
912 if(!fFOSignals) AliError("FO signals not retrieved");
918 AliITSClusterFinder *rec = 0;
919 Int_t id,module,first=0;
920 for(module=0;module<GetITSgeom()->GetIndexMax();module++){
921 id = GetITSgeom()->GetModuleType(module);
922 if (!all && !det[id]) continue;
923 if(det[id]) first = GetITSgeom()->GetStartDet(id);
924 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
925 TClonesArray *itsDigits = DigitsAddress(id);
927 AliFatal("The reconstruction class was not instanciated!");
930 ResetDigits(); // MvL: Not sure we neeed this when rereading anyways
932 treeD->GetEvent(lastentry+module);
935 treeD->GetEvent(lastentry+(module-first));
937 Int_t ndigits = itsDigits->GetEntriesFast();
938 if (ndigits>0 || id==0) { // for SPD we always want to call FindRawClusters (to process FO signals)
939 rec->SetDetTypeRec(this);
940 rec->SetDigits(DigitsAddress(id));
941 // rec->SetClusters(ClustersAddress(id));
942 rec->FindRawClusters(module);
948 // Remove PIT in-active chips from Fast-OR fired map
949 if (all || det[0]) { // SPD present
950 RemoveFastOrFiredInActive();
951 // here removing bits which have no associated clusters
952 RemoveFastOrFiredFromDead(GetFiredChipMap(treeR));
955 AliITSRecPointContainer* rpcont = AliITSRecPointContainer::Instance();
957 nClu[0]=rpcont->GetNClustersInLayer(1,treeR);
958 for(Int_t iLay=2; iLay<=6; iLay++) nClu[iLay-1]=rpcont->GetNClustersInLayerFast(iLay);
959 AliInfo(Form("Number of RecPoints in ITS Layers = %d %d %d %d %d %d",
960 nClu[0],nClu[1],nClu[2],nClu[3],nClu[4],nClu[5]));
962 //______________________________________________________________________
963 void AliITSDetTypeRec::DigitsToRecPoints(AliRawReader* rawReader,TTree *treeR,Option_t *opt){
964 // cluster finding and reconstruction of space points
965 // the condition below will disappear when the geom class will be
966 // initialized for all versions - for the moment it is only for v5 !
967 // 7 is the SDD beam test version
969 // AliRawReader *rawReader Pointer to the raw-data reader
970 // TTree *treeR Clusters tree
975 const char *all = strstr(opt,"All");
976 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
980 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
982 TClonesArray* array = rpc->UncheckedGetClusters(0);
983 TBranch *branch = treeR->Branch("ITSRecPoints",&array);
984 DigitsToRecPoints(rawReader,opt);
987 for(Int_t iModule=0;iModule<GetITSgeom()->GetIndexMax();iModule++){
988 id = GetITSgeom()->GetModuleType(iModule);
989 if (!all && !det[id]) continue;
990 array = rpc->UncheckedGetClusters(iModule);
992 AliDebug(1,Form("data for module %d missing!",iModule));
994 branch->SetAddress(&array);
996 nClusters+=array->GetEntriesFast();
1001 AliITSRecPointContainer* rpcont = AliITSRecPointContainer::Instance();
1003 nClu[0]=rpcont->GetNClustersInLayer(1,treeR);
1004 for(Int_t iLay=2; iLay<=6; iLay++) nClu[iLay-1]=rpcont->GetNClustersInLayerFast(iLay);
1005 AliInfo(Form("Number of RecPoints in ITS Layers = %d %d %d %d %d %d, Total = %d",
1006 nClu[0],nClu[1],nClu[2],nClu[3],nClu[4],nClu[5],nClusters));
1008 //______________________________________________________________________
1009 void AliITSDetTypeRec::DigitsToRecPoints(AliRawReader* rawReader,Option_t *opt){
1010 // cluster finding and reconstruction of space points
1011 // the condition below will disappear when the geom class will be
1012 // initialized for all versions - for the moment it is only for v5 !
1013 // 7 is the SDD beam test version
1015 // AliRawReader *rawReader Pointer to the raw-data reader
1020 const char *all = strstr(opt,"All");
1021 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
1024 // Reset Fast-OR fired map
1025 ResetFastOrFiredMap();
1027 AliITSClusterFinder *rec = 0;
1030 for(id=0;id<3;id++){
1031 if (!all && !det[id]) continue;
1032 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
1034 AliFatal("The reconstruction class was not instantiated");
1037 rec->SetDetTypeRec(this);
1038 rec->RawdataToClusters(rawReader);
1041 // Remove PIT in-active chips from Fast-OR fired map
1042 if (all || det[0]) { // SPD present
1043 RemoveFastOrFiredInActive();
1044 // here removing bits which have no associated clusters
1045 RemoveFastOrFiredFromDead(GetFiredChipMap());
1049 //______________________________________________________________________
1050 void AliITSDetTypeRec::ReadOldSSDNoise(const TObjArray *array,
1051 AliITSNoiseSSDv2 *noiseSSD) {
1052 //Reads the old SSD calibration object and converts it to the new format
1053 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
1054 const Int_t fgkSSDPSIDESTRIPSPERMODULE = 768;
1056 Int_t gNMod = array->GetEntries();
1057 cout<<"Converting old calibration object for noise..."<<endl;
1060 Double_t noise = 0.0;
1061 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1062 AliITSNoiseSSD *noiseModule = (AliITSNoiseSSD*) (array->At(iModule));
1063 for(Int_t iStrip = 0; iStrip < fgkSSDSTRIPSPERMODULE; iStrip++) {
1064 noise = (iStrip < fgkSSDPSIDESTRIPSPERMODULE) ? noiseModule->GetNoiseP(iStrip) : noiseModule->GetNoiseN(1535 - iStrip);
1065 if(iStrip < fgkSSDPSIDESTRIPSPERMODULE)
1066 noiseSSD->AddNoiseP(iModule,iStrip,noise);
1067 if(iStrip >= fgkSSDPSIDESTRIPSPERMODULE)
1068 noiseSSD->AddNoiseN(iModule,1535 - iStrip,noise);
1070 }//loop over modules
1073 //______________________________________________________________________
1074 void AliITSDetTypeRec::ReadOldSSDBadChannels(const TObjArray *array,
1075 AliITSBadChannelsSSDv2 *badChannelsSSD) {
1076 //Reads the old SSD calibration object and converts it to the new format
1077 Int_t gNMod = array->GetEntries();
1078 cout<<"Converting old calibration object for bad channels..."<<endl;
1079 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1080 //for (Int_t iModule = 0; iModule < 1; iModule++) {
1081 AliITSBadChannelsSSD *bad = (AliITSBadChannelsSSD*) (array->At(iModule));
1082 TArrayI arrayPSide = bad->GetBadPChannelsList();
1083 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1084 badChannelsSSD->AddBadChannelP(iModule,
1086 (Char_t)arrayPSide.At(iPCounter));
1088 TArrayI arrayNSide = bad->GetBadNChannelsList();
1089 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1090 badChannelsSSD->AddBadChannelN(iModule,
1092 (Char_t)arrayNSide.At(iNCounter));
1094 }//loop over modules
1097 //______________________________________________________________________
1098 void AliITSDetTypeRec::ReadOldSSDGain(const TObjArray *array,
1099 AliITSGainSSDv2 *gainSSD) {
1100 //Reads the old SSD calibration object and converts it to the new format
1102 Int_t gNMod = array->GetEntries();
1103 cout<<"Converting old calibration object for gain..."<<endl;
1106 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1107 AliITSGainSSD *gainModule = (AliITSGainSSD*) (array->At(iModule));
1108 TArrayF arrayPSide = gainModule->GetGainP();
1109 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1110 gainSSD->AddGainP(iModule,
1112 arrayPSide.At(iPCounter));
1113 TArrayF arrayNSide = gainModule->GetGainN();
1114 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1115 gainSSD->AddGainN(iModule,
1117 arrayNSide.At(iNCounter));
1118 }//loop over modules
1120 //______________________________________________________________________
1121 void AliITSDetTypeRec::RemoveFastOrFiredInActive() {
1122 // Removes the chips that were in-active in the pixel trigger (from fast-or fired map)
1124 if (fTriggerConditions==NULL) {
1125 AliError("Pixel trigger conditions are missing.");
1131 while (fTriggerConditions->GetNextInActiveChip(eq,hs,chip)) {
1132 UInt_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1133 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1136 //______________________________________________________________________
1137 TBits AliITSDetTypeRec::GetFiredChipMap() const {
1140 // TBits of the fired chips
1143 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
1145 TBits isfiredchip(1200);
1147 AliITSsegmentationSPD *segSPD = (AliITSsegmentationSPD*)GetSegmentationModel(0);
1149 AliError("no segmentation model for SPD available, the fired chip map is empty. Exiting");
1154 for(Int_t imod =0; imod < fgkDefaultNModulesSPD; imod++){
1155 TClonesArray *array = rpc->UncheckedGetClusters(imod);
1156 if(!array) continue;
1157 Int_t nCluster = array->GetEntriesFast();
1160 AliITSRecPoint* cluster = (AliITSRecPoint*)array->UncheckedAt(nCluster);
1161 if (cluster->GetLayer()>1)continue;
1162 Float_t local[3]={-1,-1};
1163 local[1]=cluster->GetDetLocalX();
1164 local[0]=cluster->GetDetLocalZ();
1166 Int_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(imod);
1167 Int_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(imod);
1169 segSPD->LocalToDet(0.5,local[0],row,col);
1170 Int_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(imod,col);
1171 Int_t chipkey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1172 isfiredchip.SetBitNumber(chipkey,kTRUE);
1180 //______________________________________________________________________
1181 TBits AliITSDetTypeRec::GetFiredChipMap(TTree *treeR) const{
1183 // TBits of the fired chips
1185 TBits isfiredchip(1200);
1188 AliError("no treeR. fired chip map stays empty. Exiting.");
1192 AliITSRecPointContainer* rpcont=AliITSRecPointContainer::Instance();
1193 TClonesArray *recpoints = NULL;
1194 rpcont->FetchClusters(0,treeR);
1195 if(!rpcont->GetStatusOK() || !rpcont->IsSPDActive()){
1196 AliError("no clusters. fired chip map stays empty. Exiting.");
1200 AliITSsegmentationSPD *segSPD = (AliITSsegmentationSPD*)GetSegmentationModel(0);
1202 for(Int_t imod =0; imod < fgkDefaultNModulesSPD; imod++){
1203 recpoints = rpcont->UncheckedGetClusters(imod);
1204 Int_t nCluster = recpoints->GetEntriesFast();
1206 // loop over clusters
1208 AliITSRecPoint* cluster = (AliITSRecPoint*)recpoints->UncheckedAt(nCluster);
1209 if (cluster->GetLayer()>1)continue;
1210 Float_t local[3]={-1,-1};
1211 local[1]=cluster->GetDetLocalX();
1212 local[0]=cluster->GetDetLocalZ();
1214 Int_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(imod);
1215 Int_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(imod);
1217 segSPD->LocalToDet(0.5,local[0],row,col);
1218 Int_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(imod,col);
1219 Int_t chipkey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1220 isfiredchip.SetBitNumber(chipkey,kTRUE);
1226 //______________________________________________________________________
1227 void AliITSDetTypeRec::RemoveFastOrFiredFromDead(TBits firedchipmap){
1229 // resetting of the fast-or bit on cluster basis.
1230 // fast-or bits can be remnant from SPD ideal simulation (no dead channels)
1233 for(Int_t chipKey=0; chipKey<1200; chipKey++){
1234 // FO masked chips have been previously removed
1235 if(!fFastOrFiredMap.TestBitNumber(chipKey)) continue;
1236 if(!firedchipmap.TestBitNumber(chipKey)) {
1237 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1238 AliDebug(2,Form("removing bit in key %i \n ",chipKey));
1243 //______________________________________________________________________
1244 void AliITSDetTypeRec::SetFastOrFiredMapOnline(UInt_t eq, UInt_t hs, UInt_t chip) {
1245 // Set fast-or fired map for this chip
1246 Int_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1247 return SetFastOrFiredMap(chipKey);