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 **************************************************************************/
18 ////////////////////////////////////////////////////////////////////////
19 // This class defines the "Standard" reconstruction for the ITS //
22 ////////////////////////////////////////////////////////////////////////
23 #include "TObjArray.h"
26 #include "AliCDBManager.h"
27 #include "AliCDBEntry.h"
28 #include "AliITSClusterFinder.h"
29 #include "AliITSClusterFinderV2SPD.h"
30 #include "AliITSClusterFinderV2SDD.h"
31 #include "AliITSClusterFinderSDDfast.h"
32 #include "AliITSClusterFinderV2SSD.h"
33 #include "AliITSDetTypeRec.h"
34 #include "AliITSDDLModuleMapSDD.h"
35 #include "AliITSRecPoint.h"
36 #include "AliITSRecPointContainer.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"
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;
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 fSPDSparseDead(rec.fSPDSparseDead),
124 fTriggerConditions(rec.fTriggerConditions),
125 fDigits(rec.fDigits),
126 fFOSignals(rec.fFOSignals),
127 fDDLMapSDD(rec.fDDLMapSDD),
128 fRespSDD(rec.fRespSDD),
129 fAveGainSDD(rec.fAveGainSDD),
130 fRecPoints(rec.fRecPoints),
131 fNRecPoints(rec.fNRecPoints),
132 fFirstcall(rec.fFirstcall),
133 fLoadOnlySPDCalib(rec.fLoadOnlySPDCalib),
134 fFastOrFiredMap(rec.fFastOrFiredMap){
137 for(Int_t i=0; i<3; i++){
138 fkDigClassName[i]=rec.fkDigClassName[i]; // NB only copies Char_t*, so not so safe, but this code should never be reached anyways
141 //______________________________________________________________________
142 AliITSDetTypeRec& AliITSDetTypeRec::operator=(const AliITSDetTypeRec& source){
143 // Assignment operator.
144 this->~AliITSDetTypeRec();
145 new(this) AliITSDetTypeRec(source);
150 //_____________________________________________________________________
151 AliITSDetTypeRec::~AliITSDetTypeRec(){
155 fReconstruction->Delete();
156 delete fReconstruction;
159 fSegmentation->Delete();
160 delete fSegmentation;
163 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
164 fCalibration->Delete();
166 if(fRespSDD) delete fRespSDD;
167 if(fDDLMapSDD) delete fDDLMapSDD;
171 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
172 delete fSSDCalibration;
176 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
182 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
183 fSPDSparseDead->Delete();
184 delete fSPDSparseDead;
187 if(fTriggerConditions){
188 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
189 fTriggerConditions->Delete();
190 delete fTriggerConditions;
198 fRecPoints->Delete();
203 if (fITSgeom) delete fITSgeom;
207 //___________________________________________________________________
208 void AliITSDetTypeRec::SetReconstructionModel(Int_t dettype,AliITSClusterFinder *clf){
210 //Set reconstruction model for detector type
212 if(fReconstruction==0) fReconstruction = new TObjArray(fgkNdettypes);
213 if(fReconstruction->At(dettype)!=0) delete fReconstruction->At(dettype);
214 fReconstruction->AddAt(clf,dettype);
216 //______________________________________________________________________
217 AliITSClusterFinder* AliITSDetTypeRec::GetReconstructionModel(Int_t dettype) const{
219 //Get reconstruction model for detector type
220 if(fReconstruction==0) {
221 Warning("GetReconstructionModel","fReconstruction is 0!");
224 return (AliITSClusterFinder*)fReconstruction->At(dettype);
227 //______________________________________________________________________
228 void AliITSDetTypeRec::SetSegmentationModel(Int_t dettype,AliITSsegmentation *seg){
230 //Set segmentation model for detector type
232 if(fSegmentation==0) fSegmentation = new TObjArray(fgkNdettypes);
233 if(fSegmentation->At(dettype)!=0) delete fSegmentation->At(dettype);
234 fSegmentation->AddAt(seg,dettype);
237 //______________________________________________________________________
238 AliITSsegmentation* AliITSDetTypeRec::GetSegmentationModel(Int_t dettype) const {
240 //Get segmentation model for detector type
242 if(fSegmentation==0) {
243 Warning("GetSegmentationModel","fSegmentation is 0!");
246 return (AliITSsegmentation*)fSegmentation->At(dettype);
249 //_______________________________________________________________________
250 void AliITSDetTypeRec::SetCalibrationModel(Int_t iMod, AliITSCalibration *cal){
252 //Set calibration (response) for the module iMod of type dettype
253 if (fCalibration==0) {
254 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
255 fCalibration->SetOwner(kTRUE);
256 fCalibration->Clear();
259 if (fCalibration->At(iMod) != 0)
260 delete (AliITSCalibration*) fCalibration->At(iMod);
261 fCalibration->AddAt(cal,iMod);
264 //_______________________________________________________________________
265 void AliITSDetTypeRec::SetSPDDeadModel(Int_t iMod, AliITSCalibration *cal){
267 //Set dead pixel info for the SPD module iMod
269 fSPDDead = new TObjArray(fgkDefaultNModulesSPD);
270 fSPDDead->SetOwner(kTRUE);
274 if (fSPDDead->At(iMod) != 0)
275 delete (AliITSCalibration*) fSPDDead->At(iMod);
276 fSPDDead->AddAt(cal,iMod);
278 //_______________________________________________________________________
279 void AliITSDetTypeRec::SetSPDSparseDeadModel(Int_t iMod, AliITSCalibration *cal){
281 //Set dead pixel info for the SPD ACTIVE module iMod
282 if (fSPDSparseDead==0) {
283 fSPDSparseDead = new TObjArray(fgkDefaultNModulesSPD);
284 fSPDSparseDead->SetOwner(kTRUE);
285 fSPDSparseDead->Clear();
288 if (fSPDSparseDead->At(iMod) != 0)
289 delete (AliITSCalibration*) fSPDSparseDead->At(iMod);
290 fSPDSparseDead->AddAt(cal,iMod);
292 //_______________________________________________________________________
293 AliITSCalibration* AliITSDetTypeRec::GetCalibrationModel(Int_t iMod) const {
295 //Get calibration model for module type
297 if(fCalibration==0) {
298 Warning("GetalibrationModel","fCalibration is 0!");
302 if(iMod<fgkDefaultNModulesSPD+fgkDefaultNModulesSDD){
303 return (AliITSCalibration*)fCalibration->At(iMod);
305 Int_t i=iMod-(fgkDefaultNModulesSPD+fgkDefaultNModulesSDD);
306 fSSDCalibration->SetModule(i);
307 return (AliITSCalibration*)fSSDCalibration;
311 //_______________________________________________________________________
312 AliITSCalibration* AliITSDetTypeRec::GetSPDDeadModel(Int_t iMod) const {
314 //Get SPD dead for module iMod
317 AliWarning("fSPDDead is 0!");
320 return (AliITSCalibration*)fSPDDead->At(iMod);
322 //_______________________________________________________________________
323 AliITSCalibration* AliITSDetTypeRec::GetSPDSparseDeadModel(Int_t iMod) const {
325 //Get SPD dead for module iMod
327 if(fSPDSparseDead==0) {
328 AliWarning("fSPDSparseDead is 0!");
331 return (AliITSCalibration*)fSPDSparseDead->At(iMod);
333 //_______________________________________________________________________
334 AliITSTriggerConditions* AliITSDetTypeRec::GetTriggerConditions() const {
335 //Get Pixel Trigger Conditions
336 if (fTriggerConditions==0) {
337 AliWarning("fTriggerConditions is 0!");
339 return fTriggerConditions;
341 //______________________________________________________________________
342 void AliITSDetTypeRec::SetTreeAddressD(TTree* const treeD){
343 // Set branch address for the tree of digits.
345 const char *det[4] = {"SPD","SDD","SSD","ITS"};
347 const Char_t* digclass;
352 if (fDigits == 0x0) {
353 fDigits = new TObjArray(fgkNdettypes);
358 for (i=0; i<fgkNdettypes; i++) {
359 digclass = GetDigitClassName(i);
360 fDigits->AddAt(new TClonesArray(digclass,1000),i);
361 if (fgkNdettypes==3) snprintf(branchname,29,"%sDigits%s",det[3],det[i]);
362 else snprintf(branchname,29,"%sDigits%d",det[3],i+1);
363 branch = treeD->GetBranch(branchname);
364 if (branch) branch->SetAddress(&((*fDigits)[i]));
369 //_______________________________________________________________________
370 TBranch* AliITSDetTypeRec::MakeBranchInTree(TTree* const tree,
371 const char* name, const char *classname,
372 void* address,Int_t size,Int_t splitlevel)
375 // Makes branch in given tree and diverts them to a separate file
381 Error("MakeBranchInTree","Making Branch %s Tree is NULL",name);
384 TBranch *branch = tree->GetBranch(name);
389 branch = tree->Branch(name,classname,address,size,splitlevel);
392 branch = tree->Bronch(name, "TClonesArray", address, size, splitlevel);
398 //____________________________________________________________________
399 void AliITSDetTypeRec::SetDefaults(){
401 //Set defaults for segmentation and response
404 Warning("SetDefaults","null pointer to AliITSgeomGeom !");
408 AliITSsegmentation* seg;
409 if(!GetCalibration()) {AliFatal("Exit");exit(0);}
411 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
413 seg = new AliITSsegmentationSPD();
414 SetSegmentationModel(dettype,seg);
415 SetDigitClassName(dettype,"AliITSdigitSPD");
418 seg = new AliITSsegmentationSDD();
419 if(fLoadOnlySPDCalib==kFALSE){
420 AliITSCalibrationSDD* cal=(AliITSCalibrationSDD*)GetCalibrationModel(fgkDefaultNModulesSPD+1);
421 if(cal->IsAMAt20MHz()){
422 seg->SetPadSize(seg->Dpz(0),20.);
423 seg->SetNPads(seg->Npz()/2,128);
426 SetSegmentationModel(dettype,seg);
427 SetDigitClassName(dettype,"AliITSdigitSDD");
430 AliITSsegmentationSSD* seg2 = new AliITSsegmentationSSD();
431 SetSegmentationModel(dettype,seg2);
432 SetDigitClassName(dettype,"AliITSdigitSSD");
436 //______________________________________________________________________
437 Bool_t AliITSDetTypeRec::GetCalibration() {
438 // Get Default calibration if a storage is not defined.
441 AliITSCalibration* cal = GetCalibrationModel(0);
447 // SetRunNumber((Int_t)AliCDBManager::Instance()->GetRun());
448 // Int_t run=GetRunNumber();
450 Bool_t cacheStatus = AliCDBManager::Instance()->GetCacheFlag();
451 if (fCalibration==0) {
452 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
453 fCalibration->SetOwner(!cacheStatus);
454 fCalibration->Clear();
457 Bool_t retCode=GetCalibrationSPD(cacheStatus);
458 if(retCode==kFALSE) return kFALSE;
460 if(fLoadOnlySPDCalib==kFALSE){
461 retCode=GetCalibrationSDD(cacheStatus);
462 if(retCode==kFALSE) return kFALSE;
463 retCode=GetCalibrationSSD(cacheStatus);
464 if(retCode==kFALSE) return kFALSE;
467 AliInfo(Form("%i SPD, %i SDD and %i SSD in calibration database",
468 fNMod[0], fNMod[1], fNMod[2]));
471 //______________________________________________________________________
472 Bool_t AliITSDetTypeRec::GetCalibrationSPD(Bool_t cacheStatus) {
473 // Get SPD calibration objects from OCDB
474 // dead pixel are not used for local reconstruction
477 AliCDBEntry *noisySPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDNoisy");
478 AliCDBEntry *deadSPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDDead");
479 AliCDBEntry *deadSparseSPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDSparseDead");
480 AliCDBEntry *pitCond = AliCDBManager::Instance()->Get("TRIGGER/SPD/PITConditions");
481 if(!noisySPD || !deadSPD || !pitCond ){
482 AliFatal("SPD Calibration object retrieval failed! ");
486 TObjArray *calNoisySPD = (TObjArray*) noisySPD->GetObject();
487 if (!cacheStatus) noisySPD->SetObject(NULL);
488 noisySPD->SetOwner(kTRUE);
490 TObjArray *calDeadSPD = (TObjArray*) deadSPD->GetObject();
491 if (!cacheStatus) deadSPD->SetObject(NULL);
492 deadSPD->SetOwner(kTRUE);
494 TObjArray *calSparseDeadSPD = (TObjArray*) deadSparseSPD->GetObject();
495 if (!cacheStatus) deadSparseSPD->SetObject(NULL);
496 deadSparseSPD->SetOwner(kTRUE);
499 AliITSTriggerConditions *calPitCond = (AliITSTriggerConditions*) pitCond->GetObject();
500 if (!cacheStatus) pitCond->SetObject(NULL);
501 pitCond->SetOwner(kTRUE);
506 delete deadSparseSPD;
509 if ((!calNoisySPD) || (!calDeadSPD) || (!calSparseDeadSPD) || (!calPitCond)){
510 AliWarning("Can not get SPD calibration from calibration database !");
513 fNMod[0] = calNoisySPD->GetEntries();
515 AliITSCalibration* cal;
516 for (Int_t i=0; i<fNMod[0]; i++) {
517 cal = (AliITSCalibration*) calNoisySPD->At(i);
518 SetCalibrationModel(i, cal);
519 cal = (AliITSCalibration*) calDeadSPD->At(i);
520 SetSPDDeadModel(i, cal);
521 cal = (AliITSCalibration*) calSparseDeadSPD->At(i);
522 SetSPDSparseDeadModel(i, cal);
524 fTriggerConditions = calPitCond;
529 //______________________________________________________________________
530 Bool_t AliITSDetTypeRec::GetCalibrationSDD(Bool_t cacheStatus) {
531 // Get SDD calibration objects from OCDB
533 AliCDBEntry *entrySDD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSDD");
534 AliCDBEntry *entry2SDD = AliCDBManager::Instance()->Get("ITS/Calib/RespSDD");
535 AliCDBEntry *drSpSDD = AliCDBManager::Instance()->Get("ITS/Calib/DriftSpeedSDD");
536 AliCDBEntry *ddlMapSDD = AliCDBManager::Instance()->Get("ITS/Calib/DDLMapSDD");
537 // AliCDBEntry *mapASDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsAnodeSDD");
538 AliCDBEntry *mapTSDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsTimeSDD");
540 if(!entrySDD || !entry2SDD || !drSpSDD || !ddlMapSDD || !mapTSDD ){
541 AliFatal("SDD Calibration object retrieval failed! ");
547 TObjArray *calSDD = (TObjArray *)entrySDD->GetObject();
548 if(!cacheStatus)entrySDD->SetObject(NULL);
549 entrySDD->SetOwner(kTRUE);
551 AliITSresponseSDD *pSDD = (AliITSresponseSDD*)entry2SDD->GetObject();
552 if(!cacheStatus)entry2SDD->SetObject(NULL);
553 entry2SDD->SetOwner(kTRUE);
555 TObjArray *drSp = (TObjArray *)drSpSDD->GetObject();
556 if(!cacheStatus)drSpSDD->SetObject(NULL);
557 drSpSDD->SetOwner(kTRUE);
559 AliITSDDLModuleMapSDD *ddlsdd=(AliITSDDLModuleMapSDD*)ddlMapSDD->GetObject();
560 if(!cacheStatus)ddlMapSDD->SetObject(NULL);
561 ddlMapSDD->SetOwner(kTRUE);
563 // TObjArray *mapAn = (TObjArray *)mapASDD->GetObject();
564 // if(!cacheStatus)mapASDD->SetObject(NULL);
565 // mapASDD->SetOwner(kTRUE);
567 TObjArray *mapT = (TObjArray *)mapTSDD->GetObject();
568 if(!cacheStatus)mapTSDD->SetObject(NULL);
569 mapTSDD->SetOwner(kTRUE);
572 // DB entries are deleted. In this way metadeta objects are deleted as well
582 if ((!pSDD)||(!calSDD) || (!drSp) || (!ddlsdd) || (!mapT) ){
583 AliWarning("Can not get SDD calibration from calibration database !");
587 fNMod[1] = calSDD->GetEntries();
591 AliITSCalibration* cal;
594 Bool_t oldMapFormat=kFALSE;
595 TObject* objmap=(TObject*)mapT->At(0);
596 TString cname(objmap->ClassName());
597 if(cname.CompareTo("AliITSMapSDD")==0){
599 AliInfo("SDD Maps converted to new format");
601 for(Int_t iddl=0; iddl<AliITSDDLModuleMapSDD::GetNDDLs(); iddl++){
602 for(Int_t icar=0; icar<AliITSDDLModuleMapSDD::GetNModPerDDL();icar++){
603 Int_t iMod=fDDLMapSDD->GetModuleNumber(iddl,icar);
604 if(iMod==-1) continue;
605 Int_t i=iMod - fgkDefaultNModulesSPD;
606 cal = (AliITSCalibration*) calSDD->At(i);
609 for(Int_t iAnode=0;iAnode< ((AliITSCalibrationSDD*)cal)->NOfAnodes(); iAnode++){
610 if(((AliITSCalibrationSDD*)cal)->IsBadChannel(iAnode)) continue;
611 avegain+= ((AliITSCalibrationSDD*)cal)->GetChannelGain(iAnode);
614 AliITSDriftSpeedArraySDD* arr0 = (AliITSDriftSpeedArraySDD*) drSp->At(i0);
615 AliITSDriftSpeedArraySDD* arr1 = (AliITSDriftSpeedArraySDD*) drSp->At(i1);
617 AliITSCorrMapSDD* mt0 = 0;
618 AliITSCorrMapSDD* mt1 = 0;
620 AliITSMapSDD* oldmap0=(AliITSMapSDD*)mapT->At(i0);
621 AliITSMapSDD* oldmap1=(AliITSMapSDD*)mapT->At(i1);
622 mt0=oldmap0->ConvertToNewFormat();
623 mt1=oldmap1->ConvertToNewFormat();
625 mt0=(AliITSCorrMapSDD*)mapT->At(i0);
626 mt1=(AliITSCorrMapSDD*)mapT->At(i1);
628 cal->SetDriftSpeed(0,arr0);
629 cal->SetDriftSpeed(1,arr1);
632 SetCalibrationModel(iMod, cal);
635 if(nGdAnodes) fAveGainSDD=avegain/nGdAnodes;
640 //______________________________________________________________________
641 Bool_t AliITSDetTypeRec::GetCalibrationSSD(Bool_t cacheStatus) {
642 // Get SSD calibration objects from OCDB
643 // AliCDBEntry *entrySSD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSSD");
645 AliCDBEntry *entryNoiseSSD = AliCDBManager::Instance()->Get("ITS/Calib/NoiseSSD");
646 AliCDBEntry *entryGainSSD = AliCDBManager::Instance()->Get("ITS/Calib/GainSSD");
647 AliCDBEntry *entryBadChannelsSSD = AliCDBManager::Instance()->Get("ITS/Calib/BadChannelsSSD");
649 if(!entryNoiseSSD || !entryGainSSD || !entryBadChannelsSSD){
650 AliFatal("SSD Calibration object retrieval failed! ");
654 TObject *emptyssd = 0; TString ssdobjectname;
655 AliITSNoiseSSDv2 *noiseSSD = NULL;
656 emptyssd = (TObject *)entryNoiseSSD->GetObject();
657 ssdobjectname = emptyssd->GetName();
658 if(ssdobjectname=="TObjArray") {
659 TObjArray *noiseSSDOld = (TObjArray *)entryNoiseSSD->GetObject();
660 noiseSSD = new AliITSNoiseSSDv2();
661 ReadOldSSDNoise(noiseSSDOld, noiseSSD);
663 else if(ssdobjectname=="AliITSNoiseSSDv2")
664 noiseSSD = (AliITSNoiseSSDv2 *)entryNoiseSSD->GetObject();
665 if(!cacheStatus)entryNoiseSSD->SetObject(NULL);
666 entryNoiseSSD->SetOwner(kTRUE);
668 AliITSGainSSDv2 *gainSSD = NULL;;
669 emptyssd = (TObject *)entryGainSSD->GetObject();
670 ssdobjectname = emptyssd->GetName();
671 if(ssdobjectname=="Gain") {
672 TObjArray *gainSSDOld = (TObjArray *)entryGainSSD->GetObject();
673 gainSSD = new AliITSGainSSDv2();
674 ReadOldSSDGain(gainSSDOld, gainSSD);
676 else if(ssdobjectname=="AliITSGainSSDv2")
677 gainSSD = (AliITSGainSSDv2 *)entryGainSSD->GetObject();
678 if(!cacheStatus)entryGainSSD->SetObject(NULL);
679 entryGainSSD->SetOwner(kTRUE);
681 AliITSBadChannelsSSDv2 *badChannelsSSD = NULL;
682 emptyssd = (TObject *)entryBadChannelsSSD->GetObject();
683 ssdobjectname = emptyssd->GetName();
684 if(ssdobjectname=="TObjArray") {
685 TObjArray *badChannelsSSDOld = (TObjArray *)entryBadChannelsSSD->GetObject();
686 badChannelsSSD = new AliITSBadChannelsSSDv2();
687 ReadOldSSDBadChannels(badChannelsSSDOld, badChannelsSSD);
689 else if(ssdobjectname=="AliITSBadChannelsSSDv2")
690 badChannelsSSD = (AliITSBadChannelsSSDv2*)entryBadChannelsSSD->GetObject();
691 if(!cacheStatus)entryBadChannelsSSD->SetObject(NULL);
692 entryBadChannelsSSD->SetOwner(kTRUE);
694 // DB entries are deleted. In this way metadeta objects are deleted as well
696 delete entryNoiseSSD;
698 delete entryBadChannelsSSD;
701 if ((!noiseSSD)|| (!gainSSD)|| (!badChannelsSSD)) {
702 AliWarning("Can not get SSD calibration from calibration database !");
706 fSSDCalibration->SetNoise(noiseSSD);
707 fSSDCalibration->SetGain(gainSSD);
708 fSSDCalibration->SetBadChannels(badChannelsSSD);
709 //fSSDCalibration->FillBadChipMap();
714 //________________________________________________________________
715 void AliITSDetTypeRec::SetDefaultClusterFindersV2(Bool_t rawdata, Bool_t fastSDD){
717 //Set defaults for cluster finder V2
720 Warning("SetDefaults","Null pointer to AliITSgeom !");
724 AliITSClusterFinder *clf;
726 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
729 if(!GetReconstructionModel(dettype)){
730 clf = new AliITSClusterFinderV2SPD(this);
732 if(!rawdata) clf->SetDigits(DigitsAddress(0));
733 SetReconstructionModel(dettype,clf);
739 if(!GetReconstructionModel(dettype)){
741 clf = new AliITSClusterFinderSDDfast(this);
744 clf = new AliITSClusterFinderV2SDD(this);
747 if(!rawdata) clf->SetDigits(DigitsAddress(1));
748 SetReconstructionModel(dettype,clf);
755 if(!GetReconstructionModel(dettype)){
756 clf = new AliITSClusterFinderV2SSD(this);
758 if(!rawdata) clf->SetDigits(DigitsAddress(2));
759 SetReconstructionModel(dettype,clf);
766 //______________________________________________________________________
767 void AliITSDetTypeRec::MakeBranch(TTree* tree, Option_t* option){
769 //Creates branches for clusters and recpoints
770 Bool_t cR = (strstr(option,"R")!=0);
771 Bool_t cRF = (strstr(option,"RF")!=0);
775 if(cR) MakeBranchR(tree);
776 if(cRF) MakeBranchRF(tree);
780 //___________________________________________________________________
781 void AliITSDetTypeRec::ResetDigits(){
782 // Reset number of digits and the digits array for the ITS detector.
785 for(Int_t i=0;i<fgkNdettypes;i++){
789 //___________________________________________________________________
790 void AliITSDetTypeRec::ResetDigits(Int_t branch){
791 // Reset number of digits and the digits array for this branch.
793 if(fDigits->At(branch)) ((TClonesArray*)fDigits->At(branch))->Clear();
796 //__________________________________________________________________
797 void AliITSDetTypeRec::MakeBranchR(TTree *treeR, Option_t *opt){
799 //Creates tree branches for recpoints
801 // cont char *file File name where RecPoints branch is to be written
802 // to. If blank it write the SDigits to the same
803 // file in which the Hits were found.
808 // only one branch for rec points for all detector types
809 Bool_t oFast= (strstr(opt,"Fast")!=0);
811 Char_t detname[10] = "ITS";
815 snprintf(branchname,29,"%sRecPointsF",detname);
817 snprintf(branchname,29,"%sRecPoints",detname);
820 if(!fRecPoints)fRecPoints = new TClonesArray("AliITSRecPoint",1000);
822 MakeBranchInTree(treeR,branchname,0,&fRecPoints,buffsz,99);
824 //______________________________________________________________________
825 void AliITSDetTypeRec::SetTreeAddressR(TTree* const treeR){
826 // Set branch address for the Reconstructed points Trees.
828 // TTree *treeR Tree containing the RecPoints.
834 Char_t namedet[10]="ITS";
837 if(fRecPoints==0x0) fRecPoints = new TClonesArray("AliITSRecPoint",1000);
839 snprintf(branchname,29,"%sRecPoints",namedet);
840 branch = treeR->GetBranch(branchname);
842 branch->SetAddress(&fRecPoints);
845 snprintf(branchname,29,"%sRecPointsF",namedet);
846 branch = treeR->GetBranch(branchname);
848 branch->SetAddress(&fRecPoints);
852 //____________________________________________________________________
853 void AliITSDetTypeRec::AddRecPoint(const AliITSRecPoint &r){
854 // Add a reconstructed space point to the list
856 // const AliITSRecPoint &r RecPoint class to be added to the tree
857 // of reconstructed points TreeR.
862 TClonesArray &lrecp = *fRecPoints;
863 new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
866 //______________________________________________________________________
867 void AliITSDetTypeRec::DigitsToRecPoints(TTree *treeD,TTree *treeR,Int_t lastentry,Option_t *opt, Int_t optCluFind){
868 // cluster finding and reconstruction of space points
869 // the condition below will disappear when the geom class will be
870 // initialized for all versions - for the moment it is only for v5 !
871 // 7 is the SDD beam test version
873 // TTree *treeD Digits tree
874 // TTree *treeR Clusters tree
875 // Int_t lastentry Offset for module when not all of the modules
877 // Option_t *opt String indicating which ITS sub-detectors should
878 // be processed. If ="All" then all of the ITS
879 // sub detectors are processed.
881 const char *all = strstr(opt,"All");
882 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
885 SetDefaultClusterFindersV2();
886 AliDebug(1,"V2 cluster finder has been selected \n");
888 SetDefaultClusterFindersV2(kFALSE,kTRUE);
889 AliDebug(1,"SPD and SSD V2 Cluster Finder - SDD fast Cluster Finder \n");
893 // Reset Fast-OR fired map
894 ResetFastOrFiredMap();
896 if (all || det[0]) { // SPD present
897 // Get the FO signals for this event
898 AliRunLoader* runLoader = AliRunLoader::Instance();
899 AliITSLoader* itsLoader = (AliITSLoader*) runLoader->GetLoader("ITSLoader");
901 AliError("ITS loader is NULL.");
904 fFOSignals = (AliITSFOSignalsSPD*)itsLoader->TreeD()->GetUserInfo()->FindObject("AliITSFOSignalsSPD");
905 if(!fFOSignals) AliError("FO signals not retrieved");
911 AliITSClusterFinder *rec = 0;
912 Int_t id,module,first=0;
913 for(module=0;module<GetITSgeom()->GetIndexMax();module++){
914 id = GetITSgeom()->GetModuleType(module);
915 if (!all && !det[id]) continue;
916 if(det[id]) first = GetITSgeom()->GetStartDet(id);
917 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
918 TClonesArray *itsDigits = DigitsAddress(id);
920 AliFatal("The reconstruction class was not instanciated!");
923 ResetDigits(); // MvL: Not sure we neeed this when rereading anyways
925 treeD->GetEvent(lastentry+module);
928 treeD->GetEvent(lastentry+(module-first));
930 Int_t ndigits = itsDigits->GetEntriesFast();
931 if (ndigits>0 || id==0) { // for SPD we always want to call FindRawClusters (to process FO signals)
932 rec->SetDetTypeRec(this);
933 rec->SetDigits(DigitsAddress(id));
934 // rec->SetClusters(ClustersAddress(id));
935 rec->FindRawClusters(module);
941 // Remove PIT in-active chips from Fast-OR fired map
942 if (all || det[0]) { // SPD present
943 RemoveFastOrFiredInActive();
944 // here removing bits which have no associated clusters
945 RemoveFastOrFiredFromDead(GetFiredChipMap(treeR));
948 AliITSRecPointContainer* rpcont = AliITSRecPointContainer::Instance();
950 nClu[0]=rpcont->GetNClustersInLayer(1,treeR);
951 for(Int_t iLay=2; iLay<=6; iLay++) nClu[iLay-1]=rpcont->GetNClustersInLayerFast(iLay);
952 AliInfo(Form("Number of RecPoints in ITS Layers = %d %d %d %d %d %d",
953 nClu[0],nClu[1],nClu[2],nClu[3],nClu[4],nClu[5]));
955 //______________________________________________________________________
956 void AliITSDetTypeRec::DigitsToRecPoints(AliRawReader* rawReader,TTree *treeR,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 // TTree *treeR Clusters tree
968 const char *all = strstr(opt,"All");
969 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
973 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
975 TClonesArray* array = rpc->UncheckedGetClusters(0);
976 TBranch *branch = treeR->Branch("ITSRecPoints",&array);
977 DigitsToRecPoints(rawReader,opt);
980 for(Int_t iModule=0;iModule<GetITSgeom()->GetIndexMax();iModule++){
981 id = GetITSgeom()->GetModuleType(iModule);
982 if (!all && !det[id]) continue;
983 array = rpc->UncheckedGetClusters(iModule);
985 AliDebug(1,Form("data for module %d missing!",iModule));
987 branch->SetAddress(&array);
989 nClusters+=array->GetEntriesFast();
994 AliITSRecPointContainer* rpcont = AliITSRecPointContainer::Instance();
996 nClu[0]=rpcont->GetNClustersInLayer(1,treeR);
997 for(Int_t iLay=2; iLay<=6; iLay++) nClu[iLay-1]=rpcont->GetNClustersInLayerFast(iLay);
998 AliInfo(Form("Number of RecPoints in ITS Layers = %d %d %d %d %d %d, Total = %d",
999 nClu[0],nClu[1],nClu[2],nClu[3],nClu[4],nClu[5],nClusters));
1001 //______________________________________________________________________
1002 void AliITSDetTypeRec::DigitsToRecPoints(AliRawReader* rawReader,Option_t *opt){
1003 // cluster finding and reconstruction of space points
1004 // the condition below will disappear when the geom class will be
1005 // initialized for all versions - for the moment it is only for v5 !
1006 // 7 is the SDD beam test version
1008 // AliRawReader *rawReader Pointer to the raw-data reader
1013 const char *all = strstr(opt,"All");
1014 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
1017 // Reset Fast-OR fired map
1018 ResetFastOrFiredMap();
1020 AliITSClusterFinder *rec = 0;
1023 for(id=0;id<3;id++){
1024 if (!all && !det[id]) continue;
1025 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
1027 AliFatal("The reconstruction class was not instantiated");
1030 rec->SetDetTypeRec(this);
1031 rec->RawdataToClusters(rawReader);
1034 // Remove PIT in-active chips from Fast-OR fired map
1035 if (all || det[0]) { // SPD present
1036 RemoveFastOrFiredInActive();
1037 // here removing bits which have no associated clusters
1038 RemoveFastOrFiredFromDead(GetFiredChipMap());
1042 //______________________________________________________________________
1043 void AliITSDetTypeRec::ReadOldSSDNoise(const TObjArray *array,
1044 AliITSNoiseSSDv2 *noiseSSD) {
1045 //Reads the old SSD calibration object and converts it to the new format
1046 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
1047 const Int_t fgkSSDPSIDESTRIPSPERMODULE = 768;
1049 Int_t gNMod = array->GetEntries();
1050 AliInfo("Converting old calibration object for noise...\n");
1053 Double_t noise = 0.0;
1054 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1055 AliITSNoiseSSD *noiseModule = (AliITSNoiseSSD*) (array->At(iModule));
1056 for(Int_t iStrip = 0; iStrip < fgkSSDSTRIPSPERMODULE; iStrip++) {
1057 noise = (iStrip < fgkSSDPSIDESTRIPSPERMODULE) ? noiseModule->GetNoiseP(iStrip) : noiseModule->GetNoiseN(1535 - iStrip);
1058 if(iStrip < fgkSSDPSIDESTRIPSPERMODULE)
1059 noiseSSD->AddNoiseP(iModule,iStrip,noise);
1060 if(iStrip >= fgkSSDPSIDESTRIPSPERMODULE)
1061 noiseSSD->AddNoiseN(iModule,1535 - iStrip,noise);
1063 }//loop over modules
1066 //______________________________________________________________________
1067 void AliITSDetTypeRec::ReadOldSSDBadChannels(const TObjArray *array,
1068 AliITSBadChannelsSSDv2 *badChannelsSSD) {
1069 //Reads the old SSD calibration object and converts it to the new format
1070 Int_t gNMod = array->GetEntries();
1071 AliInfo("Converting old calibration object for bad channels...");
1072 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1073 //for (Int_t iModule = 0; iModule < 1; iModule++) {
1074 AliITSBadChannelsSSD *bad = (AliITSBadChannelsSSD*) (array->At(iModule));
1075 TArrayI arrayPSide = bad->GetBadPChannelsList();
1076 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1077 badChannelsSSD->AddBadChannelP(iModule,
1079 (Char_t)arrayPSide.At(iPCounter));
1081 TArrayI arrayNSide = bad->GetBadNChannelsList();
1082 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1083 badChannelsSSD->AddBadChannelN(iModule,
1085 (Char_t)arrayNSide.At(iNCounter));
1087 }//loop over modules
1090 //______________________________________________________________________
1091 void AliITSDetTypeRec::ReadOldSSDGain(const TObjArray *array,
1092 AliITSGainSSDv2 *gainSSD) {
1093 //Reads the old SSD calibration object and converts it to the new format
1095 Int_t gNMod = array->GetEntries();
1096 AliInfo("Converting old calibration object for gain...\n");
1099 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1100 AliITSGainSSD *gainModule = (AliITSGainSSD*) (array->At(iModule));
1101 TArrayF arrayPSide = gainModule->GetGainP();
1102 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1103 gainSSD->AddGainP(iModule,
1105 arrayPSide.At(iPCounter));
1106 TArrayF arrayNSide = gainModule->GetGainN();
1107 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1108 gainSSD->AddGainN(iModule,
1110 arrayNSide.At(iNCounter));
1111 }//loop over modules
1113 //______________________________________________________________________
1114 void AliITSDetTypeRec::RemoveFastOrFiredInActive() {
1115 // Removes the chips that were in-active in the pixel trigger (from fast-or fired map)
1117 if (fTriggerConditions==NULL) {
1118 AliError("Pixel trigger conditions are missing.");
1124 while (fTriggerConditions->GetNextInActiveChip(eq,hs,chip)) {
1125 UInt_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1126 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1129 //______________________________________________________________________
1130 TBits AliITSDetTypeRec::GetFiredChipMap() const {
1133 // TBits of the fired chips
1136 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
1138 TBits isfiredchip(1200);
1140 AliITSsegmentationSPD *segSPD = (AliITSsegmentationSPD*)GetSegmentationModel(0);
1142 AliError("no segmentation model for SPD available, the fired chip map is empty. Exiting");
1147 for(Int_t imod =0; imod < fgkDefaultNModulesSPD; imod++){
1148 TClonesArray *array = rpc->UncheckedGetClusters(imod);
1149 if(!array) continue;
1150 Int_t nCluster = array->GetEntriesFast();
1153 AliITSRecPoint* cluster = (AliITSRecPoint*)array->UncheckedAt(nCluster);
1154 if (cluster->GetLayer()>1)continue;
1155 Float_t local[3]={-1,-1};
1156 local[1]=cluster->GetDetLocalX();
1157 local[0]=cluster->GetDetLocalZ();
1159 Int_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(imod);
1160 Int_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(imod);
1162 segSPD->LocalToDet(0.5,local[0],row,col);
1163 Int_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(imod,col);
1164 Int_t chipkey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1165 isfiredchip.SetBitNumber(chipkey,kTRUE);
1173 //______________________________________________________________________
1174 TBits AliITSDetTypeRec::GetFiredChipMap(TTree *treeR) const{
1176 // TBits of the fired chips
1178 TBits isfiredchip(1200);
1181 AliError("no treeR. fired chip map stays empty. Exiting.");
1185 AliITSRecPointContainer* rpcont=AliITSRecPointContainer::Instance();
1186 TClonesArray *recpoints = NULL;
1187 rpcont->FetchClusters(0,treeR);
1188 if(!rpcont->GetStatusOK() || !rpcont->IsSPDActive()){
1189 AliError("no clusters. fired chip map stays empty. Exiting.");
1193 AliITSsegmentationSPD *segSPD = (AliITSsegmentationSPD*)GetSegmentationModel(0);
1195 for(Int_t imod =0; imod < fgkDefaultNModulesSPD; imod++){
1196 recpoints = rpcont->UncheckedGetClusters(imod);
1197 Int_t nCluster = recpoints->GetEntriesFast();
1199 // loop over clusters
1201 AliITSRecPoint* cluster = (AliITSRecPoint*)recpoints->UncheckedAt(nCluster);
1202 if (cluster->GetLayer()>1)continue;
1203 Float_t local[3]={-1,-1};
1204 local[1]=cluster->GetDetLocalX();
1205 local[0]=cluster->GetDetLocalZ();
1207 Int_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(imod);
1208 Int_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(imod);
1210 segSPD->LocalToDet(0.5,local[0],row,col);
1211 Int_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(imod,col);
1212 Int_t chipkey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1213 isfiredchip.SetBitNumber(chipkey,kTRUE);
1219 //______________________________________________________________________
1220 void AliITSDetTypeRec::RemoveFastOrFiredFromDead(TBits firedchipmap){
1222 // resetting of the fast-or bit on cluster basis.
1223 // fast-or bits can be remnant from SPD ideal simulation (no dead channels)
1226 for(Int_t chipKey=0; chipKey<1200; chipKey++){
1227 // FO masked chips have been previously removed
1228 if(!fFastOrFiredMap.TestBitNumber(chipKey)) continue;
1229 if(!firedchipmap.TestBitNumber(chipKey)) {
1230 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1231 AliDebug(2,Form("removing bit in key %i \n ",chipKey));
1236 //______________________________________________________________________
1237 void AliITSDetTypeRec::SetFastOrFiredMapOnline(UInt_t eq, UInt_t hs, UInt_t chip) {
1238 // Set fast-or fired map for this chip
1239 Int_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1240 return SetFastOrFiredMap(chipKey);