In NtupleTrigger(): added loading mappig from CDB (now required)
[u/mrichter/AliRoot.git] / HLT / TPCLib / AliHLTTPCClusterFinder.cxx
CommitLineData
a38a7850 1// @(#) $Id$
4aa41877 2// Original: AliHLTClustFinderNew.cxx,v 1.29 2005/06/14 10:55:21 cvetan Exp
a38a7850 3
600e6a1b 4/**************************************************************************
9be2600f 5 * This file is property of and copyright by the ALICE HLT Project *
6 * ALICE Experiment at CERN, All rights reserved. *
600e6a1b 7 * *
9be2600f 8 * Primary Authors: Anders Vestbo <mailto:vestbo@fi.uib.no>, *
600e6a1b 9 * Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de> *
10 * Jochen Thaeder <mailto:thaeder@kip.uni-heidelberg.de> *
11 * for The ALICE Off-line Project. *
12 * *
13 * Permission to use, copy, modify and distribute this software and its *
14 * documentation strictly for non-commercial purposes is hereby granted *
15 * without fee, provided that the above copyright notice appears in all *
16 * copies and that both the copyright notice and this permission notice *
17 * appear in the supporting documentation. The authors make no claims *
18 * about the suitability of this software for any purpose. It is *
19 * provided "as is" without express or implied warranty. *
20 **************************************************************************/
a38a7850 21
600e6a1b 22/** @file AliHLTTPCClusterFinder.cxx
23 @author Anders Vestbo <mailto:vestbo@fi.uib.no>,
24 Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de>
25 Jochen Thaeder <mailto:thaeder@kip.uni-heidelberg.de>
26 @date
27 @brief Cluster Finder for the TPC
28*/
a38a7850 29
30#include "AliHLTTPCDigitReader.h"
a38a7850 31#include "AliHLTTPCRootTypes.h"
32#include "AliHLTTPCLogging.h"
33#include "AliHLTTPCClusterFinder.h"
34#include "AliHLTTPCDigitData.h"
35#include "AliHLTTPCTransform.h"
36#include "AliHLTTPCSpacePointData.h"
37#include "AliHLTTPCMemHandler.h"
46b33a24 38#include "AliHLTTPCPad.h"
a38a7850 39
40#if __GNUC__ >= 3
41using namespace std;
42#endif
43
44/** \class AliHLTTPCClusterFinder
a38a7850 45//
46// The current cluster finder for HLT
47// (Based on STAR L3)
48//
49// The cluster finder is initialized with the Init function,
50// providing the slice and patch information to work on.
51//
52// The input is a provided by the AliHLTTPCDigitReader class,
53// using the init() funktion, and the next() funktion in order
54// to get the next bin. Either packed or unpacked data can be
55// processed, dependent if one uses AliHLTTPCDigitReaderPacked
56// class or AliHLTTPCDigitReaderUnpacked class in the
57// Clusterfinder Component.
58// The resulting space points will be in the
59// array given by the SetOutputArray function.
60//
61// There are several setters which control the behaviour:
62//
63// - SetXYError(Float_t): set fixed error in XY direction
64// - SetZError(Float_t): set fixed error in Z direction
65// (used if errors are not calculated)
66// - SetDeconv(Bool_t): switch on/off deconvolution
67// - SetThreshold(UInt_t): set charge threshold for cluster
68// - SetMatchWidth(UInt_t): set the match distance in
69// time for sequences to be merged
70// - SetSTDOutput(Bool_t): switch on/off output about found clusters
71// - SetCalcErr(Bool_t): switch on/off calculation of
72// space point errors (or widths in raw system)
73// - SetRawSP(Bool_t): switch on/off convertion to raw system
74//
75//
76// Example Usage:
77//
78// AliHLTTPCFileHandler *file = new AliHLTTPCFileHandler();
79// file->SetAliInput(digitfile); //give some input file
80// for(int slice=0; slice<=35; slice++){
81// for(int patch=0; pat<6; pat++){
82// file->Init(slice,patch);
83// UInt_t ndigits=0;
84// UInt_t maxclusters=100000;
85// UInt_t pointsize = maxclusters*sizeof(AliHLTTPCSpacePointData);
86// AliHLTTPCSpacePointData *points = (AliHLTTPCSpacePointData*)memory->Allocate(pointsize);
87// AliHLTTPCDigitRowData *digits = (AliHLTTPCDigitRowData*)file->AliAltroDigits2Memory(ndigits,event);
88// AliHLTTPCClusterFinder *cf = new AliHLTTPCClusterFinder();
89// cf->SetMatchWidth(2);
90// cf->InitSlice( slice, patch, row[0], row[1], maxPoints );
91// cf->SetSTDOutput(kTRUE); //Some output to standard IO
92// cf->SetRawSP(kFALSE); //Convert space points to local system
93// cf->SetThreshold(5); //Threshold of cluster charge
94// cf->SetDeconv(kTRUE); //Deconv in pad and time direction
95// cf->SetCalcErr(kTRUE); //Calculate the errors of the spacepoints
96// cf->SetOutputArray(points); //Move the spacepoints to the array
97// cf->Read(iter->fPtr, iter->fSize ); //give the data to the cf
98// cf->ProcessDigits(); //process the rows given by init
99// Int_t npoints = cf->GetNumberOfClusters();
100// AliHLTTPCMemHandler *out= new AliHLTTPCMemHandler();
101// out->SetBinaryOutput(fname);
102// out->Memory2Binary(npoints,points); //store the spacepoints
103// out->CloseBinaryOutput();
104// delete out;
105// file->free();
106// delete cf;
107// }
108// }
a38a7850 109*/
110
111ClassImp(AliHLTTPCClusterFinder)
112
113AliHLTTPCClusterFinder::AliHLTTPCClusterFinder()
46b33a24 114 :
2a083ac4 115 fSpacePointData(NULL),
116 fDigitReader(NULL),
117 fPtr(NULL),
118 fSize(0),
46b33a24 119 fDeconvTime(kTRUE),
2a083ac4 120 fDeconvPad(kTRUE),
46b33a24 121 fStdout(kFALSE),
122 fCalcerr(kTRUE),
123 fRawSP(kFALSE),
124 fFirstRow(0),
125 fLastRow(0),
2a083ac4 126 fCurrentRow(0),
127 fCurrentSlice(0),
128 fCurrentPatch(0),
129 fMatch(1),
130 fThreshold(10),
131 fSignalThreshold(-1),
132 fNClusters(0),
133 fMaxNClusters(0),
134 fXYErr(0.2),
135 fZErr(0.3),
136 fOccupancyLimit(1.0)
a38a7850 137{
138 //constructor
46b33a24 139}
a38a7850 140
46b33a24 141AliHLTTPCClusterFinder::AliHLTTPCClusterFinder(const AliHLTTPCClusterFinder& src)
142 :
2a083ac4 143 fSpacePointData(NULL),
144 fDigitReader(src.fDigitReader),
145 fPtr(NULL),
146 fSize(src.fSize),
46b33a24 147 fDeconvTime(src.fDeconvTime),
2a083ac4 148 fDeconvPad(src.fDeconvPad),
46b33a24 149 fStdout(src.fStdout),
150 fCalcerr(src.fCalcerr),
151 fRawSP(src.fRawSP),
152 fFirstRow(src.fFirstRow),
153 fLastRow(src.fLastRow),
2a083ac4 154 fCurrentRow(src.fCurrentRow),
155 fCurrentSlice(src.fCurrentSlice),
156 fCurrentPatch(src.fCurrentPatch),
157 fMatch(src.fMatch),
158 fThreshold(src.fThreshold),
159 fSignalThreshold(src.fSignalThreshold),
160 fNClusters(src.fNClusters),
161 fMaxNClusters(src.fMaxNClusters),
162 fXYErr(src.fXYErr),
163 fZErr(src.fZErr),
164 fOccupancyLimit(src.fOccupancyLimit)
46b33a24 165{
166}
167
168AliHLTTPCClusterFinder& AliHLTTPCClusterFinder::operator=(const AliHLTTPCClusterFinder& src)
169{
170 fMatch=src.fMatch;
171 fThreshold=src.fThreshold;
84645eb0 172 fSignalThreshold=src.fSignalThreshold;
5235c3e9 173 fOccupancyLimit=src.fOccupancyLimit;
46b33a24 174 fXYErr=src.fXYErr;
175 fZErr=src.fZErr;
176 fDeconvPad=src.fDeconvPad;
177 fDeconvTime=src.fDeconvTime;
178 fStdout=src.fStdout;
179 fCalcerr=src.fCalcerr;
180 fRawSP=src.fRawSP;
181 fFirstRow=src.fFirstRow;
182 fLastRow=src.fLastRow;
183 fDigitReader=src.fDigitReader;
184 return (*this);
a38a7850 185}
186
187AliHLTTPCClusterFinder::~AliHLTTPCClusterFinder()
188{
189 //destructor
190}
191
192void AliHLTTPCClusterFinder::InitSlice(Int_t slice,Int_t patch,Int_t firstrow, Int_t lastrow,Int_t nmaxpoints)
193{
194 //init slice
195 fNClusters = 0;
196 fMaxNClusters = nmaxpoints;
197 fCurrentSlice = slice;
198 fCurrentPatch = patch;
199 fFirstRow = firstrow;
200 fLastRow = lastrow;
201}
202
203void AliHLTTPCClusterFinder::InitSlice(Int_t slice,Int_t patch,Int_t nmaxpoints)
204{
205 //init slice
206 fNClusters = 0;
207 fMaxNClusters = nmaxpoints;
208 fCurrentSlice = slice;
209 fCurrentPatch = patch;
210 fFirstRow=AliHLTTPCTransform::GetFirstRow(patch);
211 fLastRow=AliHLTTPCTransform::GetLastRow(patch);
212}
213
214void AliHLTTPCClusterFinder::SetOutputArray(AliHLTTPCSpacePointData *pt)
215{
216 //set pointer to output
217 fSpacePointData = pt;
218}
219
220void AliHLTTPCClusterFinder::Read(void* ptr,unsigned long size){
221 //set input pointer
222 fPtr = (UChar_t*)ptr;
223 fSize = size;
224}
225
226void AliHLTTPCClusterFinder::ProcessDigits()
227{
228 bool readValue = true;
229 Int_t newRow = 0;
230 Int_t rowOffset = 0;
46b33a24 231 UShort_t time=0,newTime=0;
232 UInt_t pad=0,newPad=0;
233 AliHLTTPCSignal_t charge=0;
a38a7850 234
235 fNClusters = 0;
236
237 // initialize block for reading packed data
db16520a 238 fDigitReader->InitBlock(fPtr,fSize,fFirstRow,fLastRow,fCurrentPatch,fCurrentSlice);
a38a7850 239 readValue = fDigitReader->Next();
240
5235c3e9 241 // Matthias 08.11.2006 the following return would cause termination without writing the
6c1a9d9e 242 // ClusterData and thus would block the component. I just want to have the commented line
5235c3e9 243 // here for information
244 //if (!readValue)return;
db16520a 245
a38a7850 246 pad = fDigitReader->GetPad();
247 time = fDigitReader->GetTime();
248 fCurrentRow = fDigitReader->GetRow();
249
250 if ( fCurrentPatch >= 2 ) // Outer sector, patches 2, 3, 4, 5
251 rowOffset = AliHLTTPCTransform::GetFirstRow( 2 );
252
253 fCurrentRow += rowOffset;
254
255 UInt_t lastpad = 123456789;
2a083ac4 256 const UInt_t kPadArraySize=5000;
257 const UInt_t kClusterListSize=10000;
84645eb0 258 AliClusterData *pad1[kPadArraySize]; //2 lists for internal memory=2pads
259 AliClusterData *pad2[kPadArraySize]; //2 lists for internal memory=2pads
260 AliClusterData clusterlist[kClusterListSize]; //Clusterlist
a38a7850 261
262 AliClusterData **currentPt; //List of pointers to the current pad
263 AliClusterData **previousPt; //List of pointers to the previous pad
264 currentPt = pad2;
265 previousPt = pad1;
266 UInt_t nprevious=0,ncurrent=0,ntotal=0;
267
46b33a24 268 /* quick implementation of baseline calculation and zero suppression
269 open a pad object for each pad and delete it after processing.
270 later a list of pad objects with base line history can be used
271 The whole thing only works if we really get unprocessed raw data, if
272 the data is already zero suppressed, there might be gaps in the time
273 bins.
274 */
275 Int_t gatingGridOffset=50;
276 AliHLTTPCPad baseline(gatingGridOffset, AliHLTTPCTransform::GetNTimeBins());
277 // just to make later conversion to a list of objects easier
84645eb0 278 AliHLTTPCPad* pCurrentPad=NULL;
279 if (fSignalThreshold>=0) {
280 pCurrentPad=&baseline;
281 baseline.SetThreshold(fSignalThreshold);
282 }
46b33a24 283
a38a7850 284 while ( readValue ){ // Reads through all digits in block
285
286 if(pad != lastpad){
287 //This is a new pad
288
289 //Switch the lists:
290 if(currentPt == pad2){
291 currentPt = pad1;
292 previousPt = pad2;
293 }
294 else {
295 currentPt = pad2;
296 previousPt = pad1;
297 }
298 nprevious = ncurrent;
299 ncurrent = 0;
300 if(pad != lastpad+1){
301 //this happens if there is a pad with no signal.
302 nprevious = ncurrent = 0;
303 }
304 lastpad = pad;
305 }
306
307 Bool_t newcluster = kTRUE;
308 UInt_t seqcharge=0,seqaverage=0,seqerror=0;
2a083ac4 309 AliHLTTPCSignal_t lastcharge=0;
310 UInt_t bLastWasFalling=0;
a38a7850 311 Int_t newbin=-1;
312
313
314 if(fDeconvTime){
315 redo: //This is a goto.
316
317 if(newbin > -1){
318 //bin = newbin;
319 newbin = -1;
320 }
321
322 lastcharge=0;
2a083ac4 323 bLastWasFalling = 0;
a38a7850 324 }
325
46b33a24 326 while(1){ //Loop over time bins of current pad
327 // read all the values for one pad at once to calculate the base line
328 if (pCurrentPad) {
329 if (!pCurrentPad->IsStarted()) {
330 //HLTDebug("reading data for pad %d, padrow %d", fDigitReader->GetPad(), fDigitReader->GetRow()+rowOffset);
331 pCurrentPad->SetID(fDigitReader->GetRow()+rowOffset,fDigitReader->GetPad());
332 if ((pCurrentPad->StartEvent())>=0) {
333 do {
334 if ((fDigitReader->GetRow()+rowOffset)!=pCurrentPad->GetRowNumber()) break;
335 if (fDigitReader->GetPad()!=pCurrentPad->GetPadNumber()) break;
336 pCurrentPad->SetRawData(fDigitReader->GetTime(), fDigitReader->GetSignal());
337 //HLTDebug("set raw data to pad: bin %d charge %d", fDigitReader->GetTime(), fDigitReader->GetSignal());
338 } while ((readValue = fDigitReader->Next())!=0);
339 }
340 pCurrentPad->CalculateBaseLine(AliHLTTPCTransform::GetNTimeBins()/2);
341 if (pCurrentPad->Next(kTRUE/*do zero suppression*/)==0) {
84645eb0 342 //HLTDebug("no data available after zero suppression");
46b33a24 343 pCurrentPad->StopEvent();
344 pCurrentPad->ResetHistory();
345 break;
346 }
347 time=pCurrentPad->GetCurrentPosition();
348 if (time>pCurrentPad->GetSize()) {
349 HLTError("invalid time bin for pad");
350 break;
351 }
352 }
353 }
a38a7850 354
46b33a24 355 if (pCurrentPad) {
5235c3e9 356 Float_t occupancy=pCurrentPad->GetOccupancy();
357 //HLTDebug("pad %d occupancy level: %f", pCurrentPad->GetPadNumber(), occupancy);
358 if ( occupancy < fOccupancyLimit ) {
359 charge = pCurrentPad->GetCorrectedData();
360 } else {
361 charge = 0;
362 //HLTDebug("ignoring pad %d with occupancy level %f", pCurrentPad->GetPadNumber(), occupancy);
363 }
46b33a24 364 } else {
5235c3e9 365 charge = fDigitReader->GetSignal();
46b33a24 366 }
367 //HLTDebug("get next charge value: position %d charge %d", time, charge);
a38a7850 368
738c049f 369
370 // CHARGE DEBUG
371 if (fDigitReader->GetRow() == 90){
372///// LOG(AliHLTTPCLog::kFatal,"AliHLTTPCClusterFinder::Row","row90") << "PAD=" << fDigitReader->GetPad() << " TIME=" << fDigitReader->GetTime()
373 // << " SIGNAL=" << fDigitReader->GetSignal() << ENDLOG;
374
375 }
376
a38a7850 377 if(time >= AliHLTTPCTransform::GetNTimeBins()){
46b33a24 378 HLTWarning("Timebin (%d) out of range (%d)", time, AliHLTTPCTransform::GetNTimeBins());
a38a7850 379 break;
380 }
46b33a24 381
382
a38a7850 383 //Get the current ADC-value
384 if(fDeconvTime){
385
386 //Check if the last pixel in the sequence is smaller than this
387 if(charge > lastcharge){
2a083ac4 388 if(bLastWasFalling){
a38a7850 389 newbin = 1;
390 break;
391 }
392 }
2a083ac4 393 else bLastWasFalling = 1; //last pixel was larger than this
a38a7850 394 lastcharge = charge;
395 }
396
397 //Sum the total charge of this sequence
398 seqcharge += charge;
399 seqaverage += time*charge;
400 seqerror += time*time*charge;
401
46b33a24 402 if (pCurrentPad) {
403
404 if((pCurrentPad->Next(kTRUE/*do zero suppression*/))==0) {
405 pCurrentPad->StopEvent();
406 pCurrentPad->ResetHistory();
407 if(readValue) {
408 newPad = fDigitReader->GetPad();
409 newTime = fDigitReader->GetTime();
410 newRow = fDigitReader->GetRow() + rowOffset;
411 }
412 break;
413 }
414
415 newPad=pCurrentPad->GetPadNumber();
416 newTime=pCurrentPad->GetCurrentPosition();
417 newRow=pCurrentPad->GetRowNumber();
418 } else {
a38a7850 419 readValue = fDigitReader->Next();
a38a7850 420 //Check where to stop:
421 if(!readValue) break; //No more value
422
423 newPad = fDigitReader->GetPad();
424 newTime = fDigitReader->GetTime();
425 newRow = fDigitReader->GetRow() + rowOffset;
46b33a24 426 }
a38a7850 427
428 if(newPad != pad)break; //new pad
429 if(newTime != time+1) break; //end of sequence
430
431 // pad = newpad; is equal
432 time = newTime;
433
434 }//end loop over sequence
435
46b33a24 436 //HLTDebug("ended time bin sequence loop: seqcharge=%d readValue=%d", seqcharge, readValue);
437 //HLTDebug("pad=%d newpad=%d current row=%d newrow=%d", pad, newPad, fCurrentRow, newRow);
438 if (seqcharge<=0) {
439 // with active zero suppression zero values are possible
440 continue;
441 }
a38a7850 442
443 //Calculate mean of sequence:
444 Int_t seqmean=0;
445 if(seqcharge)
446 seqmean = seqaverage/seqcharge;
447 else{
448 LOG(AliHLTTPCLog::kFatal,"AliHLTTPCClusterFinder::ProcessRow","Data")
449 <<"Error in data given to the cluster finder"<<ENDLOG;
450 seqmean = 1;
451 seqcharge = 1;
452 }
453
454 //Calculate mean in pad direction:
455 Int_t padmean = seqcharge*pad;
456 Int_t paderror = pad*padmean;
457
a38a7850 458 //Compare with results on previous pad:
84645eb0 459 for(UInt_t p=0; p<nprevious && p<kPadArraySize && ncurrent<kPadArraySize; p++){
a38a7850 460
461 //dont merge sequences on the same pad twice
462 if(previousPt[p]->fLastMergedPad==pad) continue;
463
464 Int_t difference = seqmean - previousPt[p]->fMean;
465 if(difference < -fMatch) break;
466
467 if(difference <= fMatch){ //There is a match here!!
468 AliClusterData *local = previousPt[p];
469
470 if(fDeconvPad){
471 if(seqcharge > local->fLastCharge){
472 if(local->fChargeFalling){ //The previous pad was falling
473 break; //create a new cluster
474 }
475 }
476 else local->fChargeFalling = 1;
477 local->fLastCharge = seqcharge;
478 }
479
480 //Don't create a new cluster, because we found a match
481 newcluster = kFALSE;
482
483 //Update cluster on current pad with the matching one:
484 local->fTotalCharge += seqcharge;
485 local->fPad += padmean;
486 local->fPad2 += paderror;
487 local->fTime += seqaverage;
488 local->fTime2 += seqerror;
489 local->fMean = seqmean;
490 local->fFlags++; //means we have more than one pad
491 local->fLastMergedPad = pad;
492
493 currentPt[ncurrent] = local;
494 ncurrent++;
495
496 break;
497 } //Checking for match at previous pad
498 } //Loop over results on previous pad.
499
84645eb0 500 if(newcluster && ncurrent<kPadArraySize){
a38a7850 501 //Start a new cluster. Add it to the clusterlist, and update
502 //the list of pointers to clusters in current pad.
503 //current pad will be previous pad on next pad.
504
505 //Add to the clusterlist:
506 AliClusterData *tmp = &clusterlist[ntotal];
507 tmp->fTotalCharge = seqcharge;
508 tmp->fPad = padmean;
509 tmp->fPad2 = paderror;
510 tmp->fTime = seqaverage;
511 tmp->fTime2 = seqerror;
512 tmp->fMean = seqmean;
513 tmp->fFlags = 0; //flags for single pad clusters
514 tmp->fLastMergedPad = pad;
515
516 if(fDeconvPad){
517 tmp->fChargeFalling = 0;
518 tmp->fLastCharge = seqcharge;
519 }
520
521 //Update list of pointers to previous pad:
522 currentPt[ncurrent] = &clusterlist[ntotal];
523 ntotal++;
524 ncurrent++;
525 }
526
527 if(fDeconvTime)
528 if(newbin >= 0) goto redo;
529
530 // to prevent endless loop
531 if(time >= AliHLTTPCTransform::GetNTimeBins()){
46b33a24 532 HLTWarning("Timebin (%d) out of range (%d)", time, AliHLTTPCTransform::GetNTimeBins());
a38a7850 533 break;
534 }
535
536
537 if(!readValue) break; //No more value
84645eb0 538
539 if (ntotal>=kClusterListSize || ncurrent>=kPadArraySize) {
540 HLTWarning("pad array size exceeded ntotal=%d ncurrent=%d, skip rest of the data", ntotal, ncurrent);
541 break;
542 }
a38a7850 543
544 if(fCurrentRow != newRow){
545 WriteClusters(ntotal,clusterlist);
546
547 lastpad = 123456789;
548
549 currentPt = pad2;
550 previousPt = pad1;
551 nprevious=0;
552 ncurrent=0;
553 ntotal=0;
554
555 fCurrentRow = newRow;
556 }
557
558 pad = newPad;
559 time = newTime;
db16520a 560
a38a7850 561 } // END while(readValue)
562
563 WriteClusters(ntotal,clusterlist);
564
46b33a24 565 HLTInfo("ClusterFinder found %d clusters in slice %d patch %d", fNClusters, fCurrentSlice, fCurrentPatch);
a38a7850 566
567} // ENDEND
568
569void AliHLTTPCClusterFinder::WriteClusters(Int_t nclusters,AliClusterData *list)
570{
571 //write cluster to output pointer
572 Int_t thisrow,thissector;
573 UInt_t counter = fNClusters;
574
575 for(int j=0; j<nclusters; j++)
576 {
577 if(!list[j].fFlags) continue; //discard single pad clusters
578 if(list[j].fTotalCharge < fThreshold) continue; //noise cluster
579
580 Float_t xyz[3];
581 Float_t fpad =(Float_t)list[j].fPad / list[j].fTotalCharge;
582 Float_t fpad2=fXYErr*fXYErr; //fixed given error
583 Float_t ftime =(Float_t)list[j].fTime / list[j].fTotalCharge;
584 Float_t ftime2=fZErr*fZErr; //fixed given error
585
738c049f 586
587
588
589
590
a38a7850 591 if(fCalcerr) { //calc the errors, otherwice take the fixed error
592 Int_t patch = AliHLTTPCTransform::GetPatch(fCurrentRow);
593 UInt_t q2=list[j].fTotalCharge*list[j].fTotalCharge;
594 Float_t sy2=list[j].fPad2 * list[j].fTotalCharge - list[j].fPad * list[j].fPad;
595 sy2/=q2;
596 if(sy2 < 0) {
597 LOG(AliHLTTPCLog::kError,"AliHLTTPCClusterFinder::WriteClusters","Cluster width")
598 <<"SigmaY2 negative "<<sy2<<" on row "<<fCurrentRow<<" "<<fpad<<" "<<ftime<<ENDLOG;
599 continue;
600 } else {
601 if(!fRawSP){
602 fpad2 = (sy2 + 1./12)*AliHLTTPCTransform::GetPadPitchWidth(patch)*AliHLTTPCTransform::GetPadPitchWidth(patch);
603 if(sy2 != 0){
604 fpad2*=0.108; //constants are from offline studies
605 if(patch<2)
606 fpad2*=2.07;
607 }
608 } else fpad2=sy2; //take the width not the error
609 }
610 Float_t sz2=list[j].fTime2*list[j].fTotalCharge - list[j].fTime*list[j].fTime;
611 sz2/=q2;
612 if(sz2 < 0){
613 LOG(AliHLTTPCLog::kError,"AliHLTTPCClusterFinder::WriteClusters","Cluster width")
614 <<"SigmaZ2 negative "<<sz2<<" on row "<<fCurrentRow<<" "<<fpad<<" "<<ftime<<ENDLOG;
615 continue;
616 } else {
617 if(!fRawSP){
618 ftime2 = (sz2 + 1./12)*AliHLTTPCTransform::GetZWidth()*AliHLTTPCTransform::GetZWidth();
619 if(sz2 != 0) {
620 ftime2 *= 0.169; //constants are from offline studies
621 if(patch<2)
622 ftime2 *= 1.77;
623 }
624 } else ftime2=sz2; //take the width, not the error
625 }
626 }
627 if(fStdout==kTRUE)
628 cout<<"WriteCluster: padrow "<<fCurrentRow<<" pad "<<fpad << " +- "<<fpad2<<" time "<<ftime<<" +- "<<ftime2<<" charge "<<list[j].fTotalCharge<<endl;
629
630 if(!fRawSP){
631 AliHLTTPCTransform::Slice2Sector(fCurrentSlice,fCurrentRow,thissector,thisrow);
632 AliHLTTPCTransform::Raw2Local(xyz,thissector,thisrow,fpad,ftime);
633
634 if(xyz[0]==0) LOG(AliHLTTPCLog::kError,"AliHLTTPCClustFinder","Cluster Finder")
635 <<AliHLTTPCLog::kDec<<"Zero cluster"<<ENDLOG;
636 if(fNClusters >= fMaxNClusters)
637 {
638 LOG(AliHLTTPCLog::kError,"AliHLTTPCClustFinder::WriteClusters","Cluster Finder")
639 <<AliHLTTPCLog::kDec<<"Too many clusters "<<fNClusters<<ENDLOG;
640 return;
641 }
642
643 fSpacePointData[counter].fX = xyz[0];
644 fSpacePointData[counter].fY = xyz[1];
645 fSpacePointData[counter].fZ = xyz[2];
646
647 } else {
648 fSpacePointData[counter].fX = fCurrentRow;
649 fSpacePointData[counter].fY = fpad;
650 fSpacePointData[counter].fZ = ftime;
651 }
652
653 fSpacePointData[counter].fCharge = list[j].fTotalCharge;
654 fSpacePointData[counter].fPadRow = fCurrentRow;
655 fSpacePointData[counter].fSigmaY2 = fpad2;
656 fSpacePointData[counter].fSigmaZ2 = ftime2;
657
44be0fde 658 fSpacePointData[counter].fUsed = kFALSE; // only used / set in AliHLTTPCDisplay
738c049f 659 fSpacePointData[counter].fTrackN = -1; // only used / set in AliHLTTPCDisplay
44be0fde 660
a38a7850 661 Int_t patch=fCurrentPatch;
662 if(patch==-1) patch=0; //never store negative patch number
663 fSpacePointData[counter].fID = counter
664 +((fCurrentSlice&0x7f)<<25)+((patch&0x7)<<22);//Uli
665
666#ifdef do_mc
667 Int_t trackID[3];
668 GetTrackID((Int_t)rint(fpad),(Int_t)rint(ftime),trackID);
669
670 fSpacePointData[counter].fTrackID[0] = trackID[0];
671 fSpacePointData[counter].fTrackID[1] = trackID[1];
672 fSpacePointData[counter].fTrackID[2] = trackID[2];
673
674 //cout<<"padrow "<<fCurrentRow<<" pad "<<(Int_t)rint(fpad)<<" time "<<(Int_t)rint(ftime)<<" Trackid "<<trackID[0]<<endl;
675#endif
676
677 fNClusters++;
678 counter++;
679 }
680}
681
682// STILL TO FIX ----------------------------------------------------------------------------
683
684#ifdef do_mc
685void AliHLTTPCClusterFinder::GetTrackID(Int_t pad,Int_t time,Int_t *trackID)
686{
687 //get mc id
688 AliHLTTPCDigitRowData *rowPt = (AliHLTTPCDigitRowData*)fDigitRowData;
689
690 trackID[0]=trackID[1]=trackID[2]=-2;
691 //cout<<"Looking for pad "<<pad<<" time "<<time<<endl;
692 for(Int_t i=fFirstRow; i<=fLastRow; i++){
693 if(rowPt->fRow < (UInt_t)fCurrentRow){
694 AliHLTTPCMemHandler::UpdateRowPointer(rowPt);
695 continue;
696 }
697 AliHLTTPCDigitData *digPt = (AliHLTTPCDigitData*)rowPt->fDigitData;
698 for(UInt_t j=0; j<rowPt->fNDigit; j++){
699 Int_t cpad = digPt[j].fPad;
700 Int_t ctime = digPt[j].fTime;
701 if(cpad != pad) continue;
702 if(ctime != time) continue;
703
704 trackID[0] = digPt[j].fTrackID[0];
705 trackID[1] = digPt[j].fTrackID[1];
706 trackID[2] = digPt[j].fTrackID[2];
707
708 //cout<<"Reading row "<<fCurrentRow<<" pad "<<cpad<<" time "<<ctime<<" trackID "<<digPt[j].fTrackID[0]<<endl;
709 break;
710 }
711 break;
712 }
713}
714#endif