1 /**************************************************************************
2 * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors 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 // Implementation of the ITS SDD fast clusterer class //
21 // Origin: Simone Capodicasa, Universita e INFN, capodica@to.infn.it //
23 ////////////////////////////////////////////////////////////////////////////
26 #include <TClonesArray.h>
31 #include "AliITSClusterFinderSDDfast.h"
32 #include "AliITSRecPoint.h"
33 #include "AliITSRecPointContainer.h"
34 #include "AliITSDetTypeRec.h"
35 #include "AliRawReader.h"
36 #include "AliITSRawStreamSDD.h"
37 #include "AliITSRawStreamSDDCompressed.h"
38 #include "AliITSCalibrationSDD.h"
39 #include "AliITSresponseSDD.h"
40 #include "AliITSDetTypeRec.h"
41 #include "AliITSReconstructor.h"
42 #include "AliITSsegmentationSDD.h"
43 #include "AliITSdigitSDD.h"
44 #include "AliITSgeomTGeo.h"
46 ClassImp(AliITSClusterFinderSDDfast)
48 AliITSClusterFinderSDDfast::AliITSClusterFinderSDDfast(AliITSDetTypeRec* dettyp):AliITSClusterFinder(dettyp),
56 fMaxDrTimeForTightCut(0.)
59 fNAnodes = GetSeg()->NpzHalf();
61 fNTimeBins = GetSeg()->Npx();
62 fNXbins = fNTimeBins+2;
63 AliDebug(2,Form("Cells in SDD cluster finder: Andoes=%d TimeBins=%d",fNAnodes,fNTimeBins));
64 SetPeakSelection(15.,30.,2000.);
65 fDDLBins.resize(kHybridsPerDDL);
69 //______________________________________________________________________
70 AliITSClusterFinderSDDfast::~AliITSClusterFinderSDDfast()
75 //______________________________________________________________________
76 void AliITSClusterFinderSDDfast::FindRawClusters(Int_t mod){
80 FindClustersSDD(fDigits);
84 //______________________________________________________________________
85 void AliITSClusterFinderSDDfast::FindClustersSDD(TClonesArray *digits){
87 std::vector<int> bins0;
88 std::vector<int> bins1;
89 const Int_t kMapDim=fNZbins*fNXbins/32;
92 for(Int_t j=0;j<kMapDim;++j){
95 AliITSCalibrationSDD* cal = (AliITSCalibrationSDD*)GetResp(fModule);
97 AliError(Form("Calibration object not present for SDD module %d\n",fModule));
102 Int_t i, ndigits=digits->GetEntriesFast();
103 for (i=0; i<ndigits; i++){
104 d=(AliITSdigitSDD*)digits->UncheckedAt(i);
105 Int_t ian=d->GetCoord1();
106 Int_t itb=d->GetCoord2();
107 Float_t gain=cal->GetChannelGain(ian)/fDetTypeRec->GetAverageGainSDD();
108 Float_t charge=d->GetSignal(); // returns expanded signal
109 // (10 bit, low threshold already added)
110 Float_t baseline = cal->GetBaseline(ian);
111 if(charge>baseline) charge-=baseline;
114 if(gain>0.){ // Bad channels have gain=0.
116 if(charge<cal->GetThresholdAnode(ian)) continue;
117 Int_t q=(Int_t)(charge+0.5);
120 Int_t iindex=y*fNZbins+z;
121 Float_t noise=cal->GetNoiseAfterElectronics(ian)*2.2; // applies zero suppression using the measured noise of each anode. Threshold values from ALICE-INT-1999-28 V10
124 bins0.push_back(iindex);
128 map0[iindex/32]|=(1<<(iindex%32));
135 bins1.push_back(iindex);
139 map1[iindex/32]|=(1<<(iindex%32));
144 FindClustersSDD(bins0, bins1, map0, map1, digits);
147 //______________________________________________________________________
148 void AliITSClusterFinderSDDfast::FindClustersSDD(std::vector<int>& bins0, std::vector<int>& bins1, const Int_t map0[], const Int_t map1[], TClonesArray *digits, TClonesArray *clusters, Int_t jitter){
150 static AliITSRecoParam *repa = NULL;
152 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
154 repa = AliITSRecoParam::GetHighFluxParam();
155 AliWarning("Using default AliITSRecoParam class");
158 const TGeoHMatrix *mT2L=AliITSgeomTGeo::GetTracking2LocalMatrix(fModule);
159 AliITSCalibrationSDD* cal = (AliITSCalibrationSDD*)GetResp(fModule);
161 AliError(Form("Calibration object not present for SDD module %d\n",fModule));
165 TClonesArray &cl=*clusters;
167 for (Int_t s=0; s<2; s++){
169 unsigned int binssize;
172 binssize=bins0.size();
177 binssize=bins1.size();
182 const Int_t rresto=fNZbins-1;
184 for(std::vector<int>::size_type i=0;i<binssize;i+=4){
188 Int_t resto=me%fNZbins;
190 Int_t idxs[1]={me+fNZbins};
191 for(Int_t k=0;k<1;++k)
192 if(map[idxs[k]/32]&(1<<(idxs[k]%32)))
193 for(std::vector<int>::size_type j=i+4;j<binssize;j+=4)
194 if(bins[j]==idxs[k]){
200 Int_t idxs[2]={me+1,me+fNZbins};
202 if(map[idxs[k]/32]&(1<<(idxs[k]%32)))
203 for(std::vector<int>::size_type j=i+4;j<binssize;j+=4)
204 if(bins[j]==idxs[k]){
210 for(std::vector<int>::size_type i=0;i<binssize;i+=4){
212 Int_t resto=me%fNZbins;
213 if(resto==fNZbins-1){
214 Int_t idxs[1]={me+fNZbins};
215 Int_t myid=bins[i+2];
216 for(Int_t l=0;l<1;++l){
217 if(map[idxs[l]/32]&(1<<(idxs[l]%32)))
218 for(std::vector<int>::size_type j=i+4;j<binssize;j+=4){
219 if(bins[j]==idxs[l]){
220 Int_t hisid=bins[j+2];
222 for(std::vector<int>::size_type k=2;k<binssize;k+=4)
232 Int_t idxs[2]={me+1,me+fNZbins};
233 Int_t myid=bins[i+2];
234 for(Int_t l=0;l<2;++l){
235 if(map[idxs[l]/32]&(1<<(idxs[l]%32)))
236 for(std::vector<int>::size_type j=i+4;j<binssize;j+=4){
237 if(bins[j]==idxs[l]){
238 Int_t hisid=bins[j+2];
240 for(std::vector<int>::size_type k=2;k<binssize;k+=4)
252 for(Int_t i=0;i<cid;++i)
253 for(Int_t j=0;j<12;++j)
256 Int_t kplab[cid][10];
257 for(Int_t i=0;i<cid;++i)
258 for(Int_t j=0;j<10;++j)
261 for(std::vector<int>::size_type i=0;i<binssize;i+=4){
263 Int_t me=bins[i+2]-1;
264 Int_t z=bins[i]%fNZbins;
265 Int_t x=bins[i]/fNZbins;
266 recp[me][0]+=q; //sumq
267 recp[me][1]+=z*q; //sumz
268 recp[me][2]+=x*q; //sumx
271 recp[me][3]+=z*z*q; //sigmaZ2
272 recp[me][4]+=x*x*q; //sigmaX2
281 recp[me][11]=bins[i];
284 if(recp[me][6]<z) recp[me][6]=z;
285 if(recp[me][7]>z) recp[me][7]=z;
286 if(recp[me][8]<x) recp[me][8]=x;
287 if(recp[me][9]>x) recp[me][9]=x;
290 recp[me][11]=bins[i];
296 for(Int_t ilab=0;ilab<10;++ilab)
297 kplab2[ilab]=kplab[me][ilab];
298 AliITSdigitSDD* d=(AliITSdigitSDD*)digits->UncheckedAt(bins[i+3]);
299 for (Int_t itrack=0;itrack<10;itrack++){
300 Int_t track = (d->GetTracks())[itrack];
302 AddLabel(kplab2, track);
305 for(Int_t ilab=0;ilab<10;++ilab)
306 kplab[me][ilab]=kplab2[ilab];
308 ++recp[me][5]; //nPiInClu
311 for(Int_t i=0;i<cid;++i){
312 if(recp[i][5]==0) continue;
313 if(recp[i][5]==1) continue;
315 Float_t q=recp[i][0];
317 Int_t clSizAnode=recp[i][6]-recp[i][7]+1;
318 Int_t clSizTb=recp[i][8]-recp[i][9]+1;
319 if(repa->GetUseSDDClusterSizeSelection()){
320 if(clSizTb==1) continue; // cut common mode noise spikes
321 if(clSizAnode>5) continue; // cut common mode noise spikes
322 if(clSizTb>10) continue; // cut clusters on noisy anodes
323 if(cal->IsAMAt20MHz() && clSizTb>8) continue; // cut clusters on noisy anodes
326 Float_t zz=(Float_t)recp[i][1]/q;
327 Float_t xx=(Float_t)recp[i][2]/q;
329 AliITSresponseSDD* rsdd = fDetTypeRec->GetResponseSDD();
331 Float_t zAnode=zz-0.5; // to have anode in range 0.-255. and centered on the mid of the pitch
332 Float_t timebin=xx-0.5; // to have time bin in range 0.-255. amd centered on the mid of the bin
334 if(s==1) zAnode+=fNAnodes; // right side has anodes from 256. to 511.
336 Float_t zdet=GetSeg()->GetLocalZFromAnode(zAnode);
337 Float_t driftTimeUncorr=GetSeg()->GetDriftTimeFromTb(timebin)+jitter*rsdd->GetCarlosRXClockPeriod();
338 Float_t driftTime=driftTimeUncorr-rsdd->GetTimeZero(fModule);
340 if(driftTime<fMaxDrTimeForTightCut && recp[i][10]<fCutOnPeakTight) continue;
341 if(recp[i][10]<fCutOnPeakLoose) continue;
343 Float_t driftSpeed=cal->GetDriftSpeedAtAnode(zAnode) + rsdd->GetDeltaVDrift(fModule,zAnode>255);
344 Float_t driftPathMicron=driftTime*driftSpeed;
345 const Double_t kMicronTocm=1.0e-4;
346 Float_t xdet=(driftPathMicron-GetSeg()->Dx())*kMicronTocm; // xdet is negative
347 if(s==0) xdet=-xdet; // left side has positive local x
349 if(repa->GetUseSDDCorrectionMaps()){
350 Float_t corrx=0, corrz=0;
351 cal->GetCorrections(zdet,xdet,corrz,corrx,GetSeg());
356 Double_t loc[3]={xdet,0.,zdet},trk[3]={0.,0.,0.};
357 mT2L->MasterToLocal(loc,trk);
361 q+=(driftTime*rsdd->GetADCvsDriftTime(fModule)); // correction for zero supp.
362 q/=rsdd->GetADCtokeV(fModule);
363 if(cal-> IsAMAt20MHz()) q*=2.; // account for 1/2 sampling freq.
364 if(q<repa->GetMinClusterChargeSDD()) continue; // remove noise clusters
367 Float_t hit[6]={xx,zz,recp[i][3],recp[i][4],q,0.};
369 Float_t hit[6]={xx, zz, 0.0030*0.0030, 0.0020*0.0020, q, 0.};
372 Int_t info[3]={clSizTb, clSizAnode, fNlayer[fModule]};
376 for(Int_t ilab=0;ilab<10;++ilab)
377 if(kplab[i][ilab]!=-2)
378 kplab2[ilab]=kplab[i][ilab];
383 if(fRawID2ClusID) kplab2[0]=fNClusters+1; // RS: store clID+1 as a reference to the cluster
384 for(Int_t ilab=1;ilab<10;++ilab)
387 if(digits) CheckLabels2(kplab2);
388 kplab2[3]=fNdet[fModule];
389 AliITSRecPoint cc(kplab2,hit,info);
391 cc.SetDriftTime(driftTimeUncorr);
393 cc.SetChargeRatio(recp[i][10]);
394 if(clusters) new (cl[nrp]) AliITSRecPoint(cc);
396 fDetTypeRec->AddRecPoint(cc);
402 AliDebug(2,Form("Clusters found on SDD module %d (unfolding %d) = %d\n",fModule,repa->GetUseUnfoldingInClusterFinderSDD(),nrp));
405 //______________________________________________________________________
406 void AliITSClusterFinderSDDfast::RawdataToClusters(AliRawReader* rawReader){
407 //------------------------------------------------------------
408 // This function creates ITS clusters from raw data
409 //------------------------------------------------------------
411 AliITSRawStream* inputSDD=AliITSRawStreamSDD::CreateRawStreamSDD(rawReader);
412 AliDebug(1,Form("%s is used",inputSDD->ClassName()));
414 AliITSDDLModuleMapSDD *ddlmap=(AliITSDDLModuleMapSDD*)fDetTypeRec->GetDDLModuleMapSDD();
415 inputSDD->SetDDLModuleMap(ddlmap);
416 for(Int_t iddl=0; iddl<AliITSDDLModuleMapSDD::GetNDDLs(); iddl++){
417 for(Int_t icar=0; icar<AliITSDDLModuleMapSDD::GetNModPerDDL();icar++){
418 Int_t iMod=ddlmap->GetModuleNumber(iddl,icar);
419 if(iMod==-1) continue;
420 AliITSCalibrationSDD* cal = (AliITSCalibrationSDD*)GetResp(iMod);
422 AliError(Form("Calibration object not present for SDD module %d\n",iMod));
425 Bool_t isZeroSupp=cal->GetZeroSupp();
427 for(Int_t iSid=0; iSid<2; iSid++) inputSDD->SetZeroSuppLowThreshold(iMod-240,iSid,cal->GetZSLowThreshold(iSid));
429 for(Int_t iSid=0; iSid<2; iSid++) inputSDD->SetZeroSuppLowThreshold(iMod-240,iSid,0);
433 FindClustersSDD(inputSDD);
437 void AliITSClusterFinderSDDfast::FindClustersSDD(AliITSRawStream* input){
439 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
440 Int_t nClustersSDD=0;
442 const Int_t kMapDim=fNZbins*fNXbins/32;
443 Int_t mapsDDL[kHybridsPerDDL][kMapDim];
444 for(Int_t i=0;i<kHybridsPerDDL;++i)
445 for(Int_t j=0;j<kMapDim;++j)
448 Int_t vectModId[kModulesPerDDL];
449 for(Int_t iMod=0; iMod<kModulesPerDDL; iMod++) vectModId[iMod]=-1;
450 // read raw data input stream
451 int countRW = 0; //RS
452 if (fRawID2ClusID) fRawID2ClusID->Reset(); //RS if array was provided, we shall store the rawID -> ClusterID
454 while (input->Next()) {
455 Int_t iModule = input->GetModuleID();
457 AliWarning(Form("Invalid SDD module number %d\n", iModule));
460 Int_t iCarlos=input->GetCarlosId();
461 Int_t iSide=input->GetChannel();
462 Int_t iHybrid=iCarlos*2+iSide;
464 if(input->IsCompletedModule()){
465 vectModId[iCarlos]=iModule; // store the module number
467 else if(input->IsCompletedDDL()){
468 // when all data from a DDL was read, search for clusters
469 Int_t jitter=input->GetJitter();
470 for(Int_t iMod=0;iMod<kModulesPerDDL;iMod++){
471 if(vectModId[iMod]>=0){
472 fModule = vectModId[iMod];
473 TClonesArray* clusters=rpc->UncheckedGetClusters(fModule);
474 std::vector<int> bins0;
475 std::vector<int> bins1;
476 bins0=fDDLBins[iMod*2];
477 bins1=fDDLBins[iMod*2+1];
480 for(Int_t j=0;j<kMapDim;++j){
483 for(Int_t i=0;i<kMapDim;++i){
484 map0[i]=mapsDDL[iMod*2][i];
485 map1[i]=mapsDDL[iMod*2+1][i];
488 FindClustersSDD(bins0, bins1, map0, map1, NULL, clusters,jitter);
490 Int_t nClusters = clusters->GetEntriesFast();
491 nClustersSDD += nClusters;
494 for (Int_t s=0; s<2; s++){
495 Int_t indexHyb=iMod*2+s;
496 for(std::vector<int>::size_type i=0;i<fDDLBins[indexHyb].size();++i)
497 fDDLBins[indexHyb].erase(fDDLBins[indexHyb].begin(),fDDLBins[indexHyb].end());
498 for(Int_t j=0;j<kMapDim;++j)
499 mapsDDL[indexHyb][j]=0;
503 else{ // fill the current digit into the bins array
504 if(iHybrid<0 || iHybrid>=kHybridsPerDDL){
505 AliWarning(Form("Invalid SDD hybrid number %d on module %d\n", iHybrid,iModule));
508 AliITSCalibrationSDD* cal=(AliITSCalibrationSDD*)GetResp(iModule);
510 AliError(Form("Calibration object not present for SDD module %d\n",iModule));
513 Float_t charge=input->GetSignal();
514 Int_t chan=input->GetCoord1()+fNAnodes*iSide;
515 Float_t gain=cal->GetChannelGain(chan)/fDetTypeRec->GetAverageGainSDD();;
516 Float_t baseline=cal->GetBaseline(chan);
517 if(charge>baseline) charge-=baseline;
519 if(gain>0.){ // Bad channels have gain=0
521 if(charge>=cal->GetThresholdAnode(chan)){
522 Int_t q=(Int_t)(charge+0.5);
523 Int_t iz = input->GetCoord1();
524 Int_t itb = input->GetCoord2();
528 Float_t noise=cal->GetNoiseAfterElectronics(ian)*2.2; // applies zero suppression using the measured noise of each anode. Threshold values from ALICE-INT-1999-28 V10
529 Int_t index=(itb+1)*fNZbins+(iz+1);
530 if((itb<fNTimeBins) && (iz<fNAnodes)){
531 if(q<noise) continue;
532 fDDLBins[iHybrid].push_back(index);
533 fDDLBins[iHybrid].push_back(q);
534 fDDLBins[iHybrid].push_back(0);
535 fDDLBins[iHybrid].push_back(countRW);
536 mapsDDL[iHybrid][index/32]|=(1<<(index%32));
539 AliWarning(Form("Invalid SDD cell: Anode=%d TimeBin=%d",iz,itb));
546 AliDebug(1,Form("found clusters in ITS SDD: %d", nClustersSDD));