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) {
329 fDigits = new TObjArray(fgkNdettypes);
334 for (i=0; i<fgkNdettypes; i++) {
335 digclass = GetDigitClassName(i);
336 fDigits->AddAt(new TClonesArray(digclass,1000),i);
337 if (fgkNdettypes==3) snprintf(branchname,29,"%sDigits%s",det[3],det[i]);
338 else snprintf(branchname,29,"%sDigits%d",det[3],i+1);
339 branch = treeD->GetBranch(branchname);
340 if (branch) branch->SetAddress(&((*fDigits)[i]));
345 //_______________________________________________________________________
346 TBranch* AliITSDetTypeRec::MakeBranchInTree(TTree* const tree,
347 const char* name, const char *classname,
348 void* address,Int_t size,Int_t splitlevel)
351 // Makes branch in given tree and diverts them to a separate file
357 Error("MakeBranchInTree","Making Branch %s Tree is NULL",name);
360 TBranch *branch = tree->GetBranch(name);
365 branch = tree->Branch(name,classname,address,size,splitlevel);
368 branch = tree->Bronch(name, "TClonesArray", address, size, splitlevel);
374 //____________________________________________________________________
375 void AliITSDetTypeRec::SetDefaults(){
377 //Set defaults for segmentation and response
380 Warning("SetDefaults","null pointer to AliITSgeomGeom !");
384 AliITSsegmentation* seg;
385 if(!GetCalibration()) {AliFatal("Exit");exit(0);}
387 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
389 seg = new AliITSsegmentationSPD();
390 SetSegmentationModel(dettype,seg);
391 SetDigitClassName(dettype,"AliITSdigitSPD");
394 seg = new AliITSsegmentationSDD();
395 if(fLoadOnlySPDCalib==kFALSE){
396 AliITSCalibrationSDD* cal=(AliITSCalibrationSDD*)GetCalibrationModel(fgkDefaultNModulesSPD+1);
397 if(cal->IsAMAt20MHz()){
398 seg->SetPadSize(seg->Dpz(0),20.);
399 seg->SetNPads(seg->Npz()/2,128);
402 SetSegmentationModel(dettype,seg);
403 SetDigitClassName(dettype,"AliITSdigitSDD");
406 AliITSsegmentationSSD* seg2 = new AliITSsegmentationSSD();
407 SetSegmentationModel(dettype,seg2);
408 SetDigitClassName(dettype,"AliITSdigitSSD");
412 //______________________________________________________________________
413 Bool_t AliITSDetTypeRec::GetCalibration() {
414 // Get Default calibration if a storage is not defined.
417 AliITSCalibration* cal = GetCalibrationModel(0);
423 // SetRunNumber((Int_t)AliCDBManager::Instance()->GetRun());
424 // Int_t run=GetRunNumber();
426 Bool_t cacheStatus = AliCDBManager::Instance()->GetCacheFlag();
427 if (fCalibration==0) {
428 fCalibration = new TObjArray(GetITSgeom()->GetIndexMax());
429 fCalibration->SetOwner(!cacheStatus);
430 fCalibration->Clear();
433 Bool_t retCode=GetCalibrationSPD(cacheStatus);
434 if(retCode==kFALSE) return kFALSE;
436 if(fLoadOnlySPDCalib==kFALSE){
437 retCode=GetCalibrationSDD(cacheStatus);
438 if(retCode==kFALSE) return kFALSE;
439 retCode=GetCalibrationSSD(cacheStatus);
440 if(retCode==kFALSE) return kFALSE;
443 AliInfo(Form("%i SPD, %i SDD and %i SSD in calibration database",
444 fNMod[0], fNMod[1], fNMod[2]));
447 //______________________________________________________________________
448 Bool_t AliITSDetTypeRec::GetCalibrationSPD(Bool_t cacheStatus) {
449 // Get SPD calibration objects from OCDB
450 // dead pixel are not used for local reconstruction
453 AliCDBEntry *noisySPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDNoisy");
454 AliCDBEntry *deadSPD = AliCDBManager::Instance()->Get("ITS/Calib/SPDDead");
455 AliCDBEntry *pitCond = AliCDBManager::Instance()->Get("TRIGGER/SPD/PITConditions");
456 if(!noisySPD || !deadSPD || !pitCond ){
457 AliFatal("SPD Calibration object retrieval failed! ");
461 TObjArray *calNoisySPD = (TObjArray*) noisySPD->GetObject();
462 if (!cacheStatus) noisySPD->SetObject(NULL);
463 noisySPD->SetOwner(kTRUE);
465 TObjArray *calDeadSPD = (TObjArray*) deadSPD->GetObject();
466 if (!cacheStatus) deadSPD->SetObject(NULL);
467 deadSPD->SetOwner(kTRUE);
469 AliITSTriggerConditions *calPitCond = (AliITSTriggerConditions*) pitCond->GetObject();
470 if (!cacheStatus) pitCond->SetObject(NULL);
471 pitCond->SetOwner(kTRUE);
478 if ((!calNoisySPD) || (!calDeadSPD) || (!calPitCond)){
479 AliWarning("Can not get SPD calibration from calibration database !");
483 fNMod[0] = calNoisySPD->GetEntries();
485 AliITSCalibration* cal;
486 for (Int_t i=0; i<fNMod[0]; i++) {
487 cal = (AliITSCalibration*) calNoisySPD->At(i);
488 SetCalibrationModel(i, cal);
489 cal = (AliITSCalibration*) calDeadSPD->At(i);
490 SetSPDDeadModel(i, cal);
492 fTriggerConditions = calPitCond;
497 //______________________________________________________________________
498 Bool_t AliITSDetTypeRec::GetCalibrationSDD(Bool_t cacheStatus) {
499 // Get SDD calibration objects from OCDB
501 AliCDBEntry *entrySDD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSDD");
502 AliCDBEntry *entry2SDD = AliCDBManager::Instance()->Get("ITS/Calib/RespSDD");
503 AliCDBEntry *drSpSDD = AliCDBManager::Instance()->Get("ITS/Calib/DriftSpeedSDD");
504 AliCDBEntry *ddlMapSDD = AliCDBManager::Instance()->Get("ITS/Calib/DDLMapSDD");
505 // AliCDBEntry *mapASDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsAnodeSDD");
506 AliCDBEntry *mapTSDD = AliCDBManager::Instance()->Get("ITS/Calib/MapsTimeSDD");
508 if(!entrySDD || !entry2SDD || !drSpSDD || !ddlMapSDD || !mapTSDD ){
509 AliFatal("SDD Calibration object retrieval failed! ");
515 TObjArray *calSDD = (TObjArray *)entrySDD->GetObject();
516 if(!cacheStatus)entrySDD->SetObject(NULL);
517 entrySDD->SetOwner(kTRUE);
519 AliITSresponseSDD *pSDD = (AliITSresponseSDD*)entry2SDD->GetObject();
520 if(!cacheStatus)entry2SDD->SetObject(NULL);
521 entry2SDD->SetOwner(kTRUE);
523 TObjArray *drSp = (TObjArray *)drSpSDD->GetObject();
524 if(!cacheStatus)drSpSDD->SetObject(NULL);
525 drSpSDD->SetOwner(kTRUE);
527 AliITSDDLModuleMapSDD *ddlsdd=(AliITSDDLModuleMapSDD*)ddlMapSDD->GetObject();
528 if(!cacheStatus)ddlMapSDD->SetObject(NULL);
529 ddlMapSDD->SetOwner(kTRUE);
531 // TObjArray *mapAn = (TObjArray *)mapASDD->GetObject();
532 // if(!cacheStatus)mapASDD->SetObject(NULL);
533 // mapASDD->SetOwner(kTRUE);
535 TObjArray *mapT = (TObjArray *)mapTSDD->GetObject();
536 if(!cacheStatus)mapTSDD->SetObject(NULL);
537 mapTSDD->SetOwner(kTRUE);
540 // DB entries are deleted. In this way metadeta objects are deleted as well
550 if ((!pSDD)||(!calSDD) || (!drSp) || (!ddlsdd) || (!mapT) ){
551 AliWarning("Can not get SDD calibration from calibration database !");
555 fNMod[1] = calSDD->GetEntries();
559 AliITSCalibration* cal;
562 Bool_t oldMapFormat=kFALSE;
563 TObject* objmap=(TObject*)mapT->At(0);
564 TString cname(objmap->ClassName());
565 if(cname.CompareTo("AliITSMapSDD")==0){
567 AliInfo("SDD Maps converted to new format");
569 for(Int_t iddl=0; iddl<AliITSDDLModuleMapSDD::GetNDDLs(); iddl++){
570 for(Int_t icar=0; icar<AliITSDDLModuleMapSDD::GetNModPerDDL();icar++){
571 Int_t iMod=fDDLMapSDD->GetModuleNumber(iddl,icar);
572 if(iMod==-1) continue;
573 Int_t i=iMod - fgkDefaultNModulesSPD;
574 cal = (AliITSCalibration*) calSDD->At(i);
577 for(Int_t iAnode=0;iAnode< ((AliITSCalibrationSDD*)cal)->NOfAnodes(); iAnode++){
578 if(((AliITSCalibrationSDD*)cal)->IsBadChannel(iAnode)) continue;
579 avegain+= ((AliITSCalibrationSDD*)cal)->GetChannelGain(iAnode);
582 AliITSDriftSpeedArraySDD* arr0 = (AliITSDriftSpeedArraySDD*) drSp->At(i0);
583 AliITSDriftSpeedArraySDD* arr1 = (AliITSDriftSpeedArraySDD*) drSp->At(i1);
585 AliITSCorrMapSDD* mt0 = 0;
586 AliITSCorrMapSDD* mt1 = 0;
588 AliITSMapSDD* oldmap0=(AliITSMapSDD*)mapT->At(i0);
589 AliITSMapSDD* oldmap1=(AliITSMapSDD*)mapT->At(i1);
590 mt0=oldmap0->ConvertToNewFormat();
591 mt1=oldmap1->ConvertToNewFormat();
593 mt0=(AliITSCorrMapSDD*)mapT->At(i0);
594 mt1=(AliITSCorrMapSDD*)mapT->At(i1);
596 cal->SetDriftSpeed(0,arr0);
597 cal->SetDriftSpeed(1,arr1);
600 SetCalibrationModel(iMod, cal);
603 if(nGdAnodes) fAveGainSDD=avegain/nGdAnodes;
608 //______________________________________________________________________
609 Bool_t AliITSDetTypeRec::GetCalibrationSSD(Bool_t cacheStatus) {
610 // Get SSD calibration objects from OCDB
611 // AliCDBEntry *entrySSD = AliCDBManager::Instance()->Get("ITS/Calib/CalibSSD");
613 AliCDBEntry *entryNoiseSSD = AliCDBManager::Instance()->Get("ITS/Calib/NoiseSSD");
614 AliCDBEntry *entryGainSSD = AliCDBManager::Instance()->Get("ITS/Calib/GainSSD");
615 AliCDBEntry *entryBadChannelsSSD = AliCDBManager::Instance()->Get("ITS/Calib/BadChannelsSSD");
617 if(!entryNoiseSSD || !entryGainSSD || !entryBadChannelsSSD){
618 AliFatal("SSD Calibration object retrieval failed! ");
622 TObject *emptyssd = 0; TString ssdobjectname;
623 AliITSNoiseSSDv2 *noiseSSD = NULL;
624 emptyssd = (TObject *)entryNoiseSSD->GetObject();
625 ssdobjectname = emptyssd->GetName();
626 if(ssdobjectname=="TObjArray") {
627 TObjArray *noiseSSDOld = (TObjArray *)entryNoiseSSD->GetObject();
628 noiseSSD = new AliITSNoiseSSDv2();
629 ReadOldSSDNoise(noiseSSDOld, noiseSSD);
631 else if(ssdobjectname=="AliITSNoiseSSDv2")
632 noiseSSD = (AliITSNoiseSSDv2 *)entryNoiseSSD->GetObject();
633 if(!cacheStatus)entryNoiseSSD->SetObject(NULL);
634 entryNoiseSSD->SetOwner(kTRUE);
636 AliITSGainSSDv2 *gainSSD = NULL;;
637 emptyssd = (TObject *)entryGainSSD->GetObject();
638 ssdobjectname = emptyssd->GetName();
639 if(ssdobjectname=="Gain") {
640 TObjArray *gainSSDOld = (TObjArray *)entryGainSSD->GetObject();
641 gainSSD = new AliITSGainSSDv2();
642 ReadOldSSDGain(gainSSDOld, gainSSD);
644 else if(ssdobjectname=="AliITSGainSSDv2")
645 gainSSD = (AliITSGainSSDv2 *)entryGainSSD->GetObject();
646 if(!cacheStatus)entryGainSSD->SetObject(NULL);
647 entryGainSSD->SetOwner(kTRUE);
649 AliITSBadChannelsSSDv2 *badChannelsSSD = NULL;
650 emptyssd = (TObject *)entryBadChannelsSSD->GetObject();
651 ssdobjectname = emptyssd->GetName();
652 if(ssdobjectname=="TObjArray") {
653 TObjArray *badChannelsSSDOld = (TObjArray *)entryBadChannelsSSD->GetObject();
654 badChannelsSSD = new AliITSBadChannelsSSDv2();
655 ReadOldSSDBadChannels(badChannelsSSDOld, badChannelsSSD);
657 else if(ssdobjectname=="AliITSBadChannelsSSDv2")
658 badChannelsSSD = (AliITSBadChannelsSSDv2*)entryBadChannelsSSD->GetObject();
659 if(!cacheStatus)entryBadChannelsSSD->SetObject(NULL);
660 entryBadChannelsSSD->SetOwner(kTRUE);
662 // DB entries are deleted. In this way metadeta objects are deleted as well
664 delete entryNoiseSSD;
666 delete entryBadChannelsSSD;
669 if ((!noiseSSD)|| (!gainSSD)|| (!badChannelsSSD)) {
670 AliWarning("Can not get SSD calibration from calibration database !");
674 fSSDCalibration->SetNoise(noiseSSD);
675 fSSDCalibration->SetGain(gainSSD);
676 fSSDCalibration->SetBadChannels(badChannelsSSD);
677 //fSSDCalibration->FillBadChipMap();
682 //________________________________________________________________
683 void AliITSDetTypeRec::SetDefaultClusterFindersV2(Bool_t rawdata){
685 //Set defaults for cluster finder V2
688 Warning("SetDefaults","Null pointer to AliITSgeom !");
692 AliITSClusterFinder *clf;
694 for(Int_t dettype=0;dettype<fgkNdettypes;dettype++){
697 if(!GetReconstructionModel(dettype)){
698 clf = new AliITSClusterFinderV2SPD(this);
700 if(!rawdata) clf->SetDigits(DigitsAddress(0));
701 SetReconstructionModel(dettype,clf);
707 if(!GetReconstructionModel(dettype)){
708 clf = new AliITSClusterFinderV2SDD(this);
710 if(!rawdata) clf->SetDigits(DigitsAddress(1));
711 SetReconstructionModel(dettype,clf);
718 if(!GetReconstructionModel(dettype)){
719 clf = new AliITSClusterFinderV2SSD(this);
721 if(!rawdata) clf->SetDigits(DigitsAddress(2));
722 SetReconstructionModel(dettype,clf);
729 //______________________________________________________________________
730 void AliITSDetTypeRec::MakeBranch(TTree* tree, Option_t* option){
732 //Creates branches for clusters and recpoints
733 Bool_t cR = (strstr(option,"R")!=0);
734 Bool_t cRF = (strstr(option,"RF")!=0);
738 if(cR) MakeBranchR(tree);
739 if(cRF) MakeBranchRF(tree);
743 //___________________________________________________________________
744 void AliITSDetTypeRec::ResetDigits(){
745 // Reset number of digits and the digits array for the ITS detector.
748 for(Int_t i=0;i<fgkNdettypes;i++){
752 //___________________________________________________________________
753 void AliITSDetTypeRec::ResetDigits(Int_t branch){
754 // Reset number of digits and the digits array for this branch.
756 if(fDigits->At(branch)) ((TClonesArray*)fDigits->At(branch))->Clear();
759 //__________________________________________________________________
760 void AliITSDetTypeRec::MakeBranchR(TTree *treeR, Option_t *opt){
762 //Creates tree branches for recpoints
764 // cont char *file File name where RecPoints branch is to be written
765 // to. If blank it write the SDigits to the same
766 // file in which the Hits were found.
771 // only one branch for rec points for all detector types
772 Bool_t oFast= (strstr(opt,"Fast")!=0);
774 Char_t detname[10] = "ITS";
778 snprintf(branchname,29,"%sRecPointsF",detname);
780 snprintf(branchname,29,"%sRecPoints",detname);
783 if(!fRecPoints)fRecPoints = new TClonesArray("AliITSRecPoint",1000);
785 MakeBranchInTree(treeR,branchname,0,&fRecPoints,buffsz,99);
787 //______________________________________________________________________
788 void AliITSDetTypeRec::SetTreeAddressR(TTree* const treeR){
789 // Set branch address for the Reconstructed points Trees.
791 // TTree *treeR Tree containing the RecPoints.
797 Char_t namedet[10]="ITS";
800 if(fRecPoints==0x0) fRecPoints = new TClonesArray("AliITSRecPoint",1000);
802 snprintf(branchname,29,"%sRecPoints",namedet);
803 branch = treeR->GetBranch(branchname);
805 branch->SetAddress(&fRecPoints);
808 snprintf(branchname,29,"%sRecPointsF",namedet);
809 branch = treeR->GetBranch(branchname);
811 branch->SetAddress(&fRecPoints);
815 //____________________________________________________________________
816 void AliITSDetTypeRec::AddRecPoint(const AliITSRecPoint &r){
817 // Add a reconstructed space point to the list
819 // const AliITSRecPoint &r RecPoint class to be added to the tree
820 // of reconstructed points TreeR.
826 TClonesArray &lrecp = *fRecPoints;
827 new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
830 //______________________________________________________________________
831 void AliITSDetTypeRec::DigitsToRecPoints(TTree *treeD,TTree *treeR,Int_t lastentry,Option_t *opt, Int_t optCluFind){
832 // cluster finding and reconstruction of space points
833 // the condition below will disappear when the geom class will be
834 // initialized for all versions - for the moment it is only for v5 !
835 // 7 is the SDD beam test version
837 // TTree *treeD Digits tree
838 // TTree *treeR Clusters tree
839 // Int_t lastentry Offset for module when not all of the modules
841 // Option_t *opt String indicating which ITS sub-detectors should
842 // be processed. If ="All" then all of the ITS
843 // sub detectors are processed.
845 const char *all = strstr(opt,"All");
846 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
849 SetDefaultClusterFindersV2();
850 AliDebug(1,"V2 cluster finder has been selected \n");
852 SetDefaultClusterFindersV2();
853 AliInfo("Cluster Finder Option not implemented, V2 cluster finder will be used \n");
857 // Reset Fast-OR fired map
858 ResetFastOrFiredMap();
860 if (all || det[0]) { // SPD present
861 // Get the FO signals for this event
862 AliRunLoader* runLoader = AliRunLoader::Instance();
863 AliITSLoader* itsLoader = (AliITSLoader*) runLoader->GetLoader("ITSLoader");
865 AliError("ITS loader is NULL.");
868 fFOSignals = (AliITSFOSignalsSPD*)itsLoader->TreeD()->GetUserInfo()->FindObject("AliITSFOSignalsSPD");
869 if(!fFOSignals) AliError("FO signals not retrieved");
875 AliITSClusterFinder *rec = 0;
876 Int_t id,module,first=0;
877 for(module=0;module<GetITSgeom()->GetIndexMax();module++){
878 id = GetITSgeom()->GetModuleType(module);
879 if (!all && !det[id]) continue;
880 if(det[id]) first = GetITSgeom()->GetStartDet(id);
881 rec = (AliITSClusterFinder*)GetReconstructionModel(id);
882 TClonesArray *itsDigits = DigitsAddress(id);
884 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");
994 rec->SetDetTypeRec(this);
995 rec->RawdataToClusters(rawReader);
998 // Remove PIT in-active chips from Fast-OR fired map
999 if (all || det[0]) { // SPD present
1000 RemoveFastOrFiredInActive();
1001 // here removing bits which have no associated clusters
1002 RemoveFastOrFiredFromDead(GetFiredChipMap());
1006 //______________________________________________________________________
1007 void AliITSDetTypeRec::ReadOldSSDNoise(const TObjArray *array,
1008 AliITSNoiseSSDv2 *noiseSSD) {
1009 //Reads the old SSD calibration object and converts it to the new format
1010 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
1011 const Int_t fgkSSDPSIDESTRIPSPERMODULE = 768;
1013 Int_t gNMod = array->GetEntries();
1014 cout<<"Converting old calibration object for noise..."<<endl;
1017 Double_t noise = 0.0;
1018 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1019 AliITSNoiseSSD *noiseModule = (AliITSNoiseSSD*) (array->At(iModule));
1020 for(Int_t iStrip = 0; iStrip < fgkSSDSTRIPSPERMODULE; iStrip++) {
1021 noise = (iStrip < fgkSSDPSIDESTRIPSPERMODULE) ? noiseModule->GetNoiseP(iStrip) : noiseModule->GetNoiseN(1535 - iStrip);
1022 if(iStrip < fgkSSDPSIDESTRIPSPERMODULE)
1023 noiseSSD->AddNoiseP(iModule,iStrip,noise);
1024 if(iStrip >= fgkSSDPSIDESTRIPSPERMODULE)
1025 noiseSSD->AddNoiseN(iModule,1535 - iStrip,noise);
1027 }//loop over modules
1030 //______________________________________________________________________
1031 void AliITSDetTypeRec::ReadOldSSDBadChannels(const TObjArray *array,
1032 AliITSBadChannelsSSDv2 *badChannelsSSD) {
1033 //Reads the old SSD calibration object and converts it to the new format
1034 Int_t gNMod = array->GetEntries();
1035 cout<<"Converting old calibration object for bad channels..."<<endl;
1036 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1037 //for (Int_t iModule = 0; iModule < 1; iModule++) {
1038 AliITSBadChannelsSSD *bad = (AliITSBadChannelsSSD*) (array->At(iModule));
1039 TArrayI arrayPSide = bad->GetBadPChannelsList();
1040 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1041 badChannelsSSD->AddBadChannelP(iModule,
1043 (Char_t)arrayPSide.At(iPCounter));
1045 TArrayI arrayNSide = bad->GetBadNChannelsList();
1046 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1047 badChannelsSSD->AddBadChannelN(iModule,
1049 (Char_t)arrayNSide.At(iNCounter));
1051 }//loop over modules
1054 //______________________________________________________________________
1055 void AliITSDetTypeRec::ReadOldSSDGain(const TObjArray *array,
1056 AliITSGainSSDv2 *gainSSD) {
1057 //Reads the old SSD calibration object and converts it to the new format
1059 Int_t gNMod = array->GetEntries();
1060 cout<<"Converting old calibration object for gain..."<<endl;
1063 for (Int_t iModule = 0; iModule < gNMod; iModule++) {
1064 AliITSGainSSD *gainModule = (AliITSGainSSD*) (array->At(iModule));
1065 TArrayF arrayPSide = gainModule->GetGainP();
1066 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
1067 gainSSD->AddGainP(iModule,
1069 arrayPSide.At(iPCounter));
1070 TArrayF arrayNSide = gainModule->GetGainN();
1071 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
1072 gainSSD->AddGainN(iModule,
1074 arrayNSide.At(iNCounter));
1075 }//loop over modules
1077 //______________________________________________________________________
1078 void AliITSDetTypeRec::RemoveFastOrFiredInActive() {
1079 // Removes the chips that were in-active in the pixel trigger (from fast-or fired map)
1081 if (fTriggerConditions==NULL) {
1082 AliError("Pixel trigger conditions are missing.");
1088 while (fTriggerConditions->GetNextInActiveChip(eq,hs,chip)) {
1089 UInt_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1090 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1093 //______________________________________________________________________
1094 TBits AliITSDetTypeRec::GetFiredChipMap() const {
1097 // TBits of the fired chips
1100 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
1102 TBits isfiredchip(1200);
1104 AliITSsegmentationSPD *segSPD = (AliITSsegmentationSPD*)GetSegmentationModel(0);
1106 AliError("no segmentation model for SPD available, the fired chip map is empty. Exiting");
1111 for(Int_t imod =0; imod < fgkDefaultNModulesSPD; imod++){
1112 TClonesArray *array = rpc->UncheckedGetClusters(imod);
1113 if(!array) continue;
1114 Int_t nCluster = array->GetEntriesFast();
1117 AliITSRecPoint* cluster = (AliITSRecPoint*)array->UncheckedAt(nCluster);
1118 if (cluster->GetLayer()>1)continue;
1119 Float_t local[3]={-1,-1};
1120 local[1]=cluster->GetDetLocalX();
1121 local[0]=cluster->GetDetLocalZ();
1123 Int_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(imod);
1124 Int_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(imod);
1126 segSPD->LocalToDet(0.5,local[0],row,col);
1127 Int_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(imod,col);
1128 Int_t chipkey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1129 isfiredchip.SetBitNumber(chipkey,kTRUE);
1137 //______________________________________________________________________
1138 TBits AliITSDetTypeRec::GetFiredChipMap(TTree *treeR) const{
1140 // TBits of the fired chips
1142 TBits isfiredchip(1200);
1145 AliError("no treeR. fired chip map stays empty. Exiting.");
1149 AliITSRecPointContainer* rpcont=AliITSRecPointContainer::Instance();
1150 TClonesArray *recpoints = NULL;
1151 rpcont->FetchClusters(0,treeR);
1152 if(!rpcont->GetStatusOK() || !rpcont->IsSPDActive()){
1153 AliError("no clusters. fired chip map stays empty. Exiting.");
1157 AliITSsegmentationSPD *segSPD = (AliITSsegmentationSPD*)GetSegmentationModel(0);
1159 for(Int_t imod =0; imod < fgkDefaultNModulesSPD; imod++){
1160 recpoints = rpcont->UncheckedGetClusters(imod);
1161 Int_t nCluster = recpoints->GetEntriesFast();
1163 // loop over clusters
1165 AliITSRecPoint* cluster = (AliITSRecPoint*)recpoints->UncheckedAt(nCluster);
1166 if (cluster->GetLayer()>1)continue;
1167 Float_t local[3]={-1,-1};
1168 local[1]=cluster->GetDetLocalX();
1169 local[0]=cluster->GetDetLocalZ();
1171 Int_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(imod);
1172 Int_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(imod);
1174 segSPD->LocalToDet(0.5,local[0],row,col);
1175 Int_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(imod,col);
1176 Int_t chipkey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1177 isfiredchip.SetBitNumber(chipkey,kTRUE);
1183 //______________________________________________________________________
1184 void AliITSDetTypeRec::RemoveFastOrFiredFromDead(TBits firedchipmap){
1186 // resetting of the fast-or bit on cluster basis.
1187 // fast-or bits can be remnant from SPD ideal simulation (no dead channels)
1190 for(Int_t chipKey=0; chipKey<1200; chipKey++){
1191 // FO masked chips have been previously removed
1192 if(!fFastOrFiredMap.TestBitNumber(chipKey)) continue;
1193 if(!firedchipmap.TestBitNumber(chipKey)) {
1194 fFastOrFiredMap.SetBitNumber(chipKey,kFALSE);
1195 AliDebug(2,Form("removing bit in key %i \n ",chipKey));
1200 //______________________________________________________________________
1201 void AliITSDetTypeRec::SetFastOrFiredMapOnline(UInt_t eq, UInt_t hs, UInt_t chip) {
1202 // Set fast-or fired map for this chip
1203 Int_t chipKey = AliITSRawStreamSPD::GetOfflineChipKeyFromOnline(eq,hs,chip);
1204 return SetFastOrFiredMap(chipKey);