]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITSClusterFinderV2SSD.cxx
Changes to set the number of layers in trackerUpgrade and trackU from the reconstructor.
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderV2SSD.cxx
CommitLineData
04366a57 1/**************************************************************************
b42cfa25 2 * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
04366a57 3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
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 **************************************************************************/
308b5ea4 15
16/* $Id$ */
17
04366a57 18////////////////////////////////////////////////////////////////////////////
19// Implementation of the ITS clusterer V2 class //
20// //
21// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
d695268b 22// Last revision: 13-05-09 Enrico Fragiacomo //
23// enrico.fragiacomo@ts.infn.it //
04366a57 24// //
25///////////////////////////////////////////////////////////////////////////
26
66b89079 27#include "AliITSClusterFinderV2SSD.h"
28
308b5ea4 29#include <Riostream.h>
66b89079 30#include <TGeoGlobalMagField.h>
04366a57 31
66b89079 32#include "AliLog.h"
33#include "AliMagF.h"
00a7cc50 34#include "AliITSRecPoint.h"
01ef1bd4 35#include "AliITSRecPointContainer.h"
1f3e997f 36#include "AliITSgeomTGeo.h"
7d62fb64 37#include "AliITSDetTypeRec.h"
04366a57 38#include "AliRawReader.h"
39#include "AliITSRawStreamSSD.h"
04366a57 40#include <TClonesArray.h>
138df073 41#include <TCollection.h>
04366a57 42#include "AliITSdigitSSD.h"
a86176e3 43#include "AliITSReconstructor.h"
3a4139a2 44#include "AliITSCalibrationSSD.h"
8be4e1b1 45#include "AliITSsegmentationSSD.h"
04366a57 46
308b5ea4 47Short_t *AliITSClusterFinderV2SSD::fgPairs = 0x0;
48Int_t AliITSClusterFinderV2SSD::fgPairsSize = 0;
d695268b 49const Float_t AliITSClusterFinderV2SSD::fgkThreshold = 5.;
50
42ed6062 51const Float_t AliITSClusterFinderV2SSD::fgkCosmic2008StripShifts[16][9] =
52 {{-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35}, // DDL 512
53 {-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35}, // DDL 513
54 {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 514
55 {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 515
56 { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00}, // DDL 516
57 { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00}, // DDL 517
58 {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 518
59 {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 519
60 {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.25,-0.15}, // DDL 520
61 {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 521
62 {-0.10,-0.10,-0.10,-0.40,-0.40,-0.40,-0.10,-0.10,-0.45}, // DDL 522
63 {-0.10,-0.10,-0.10,-0.35,-0.35,-0.35,-0.10,-0.35,-0.50}, // DDL 523
64 { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00}, // DDL 524
65 { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00}, // DDL 525
66 { 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35}, // DDL 526
67 { 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45}}; // DDL 527
c157c94e 68
04366a57 69ClassImp(AliITSClusterFinderV2SSD)
70
71
66b89079 72 AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(AliITSDetTypeRec* dettyp):AliITSClusterFinder(dettyp),fLastSSD1(AliITSgeomTGeo::GetModuleIndex(6,1,1)-1), fLorentzShiftP(0), fLorentzShiftN(0)
d036ccd3 73{
8be4e1b1 74//Default constructor
66b89079 75 static AliITSRecoParam *repa = NULL;
76 if(!repa){
77 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
78 if(!repa){
79 repa = AliITSRecoParam::GetHighFluxParam();
80 AliWarning("Using default AliITSRecoParam class");
81 }
82 }
04366a57 83
66b89079 84 if (repa->GetCorrectLorentzAngleSSD()) {
85 AliMagF* field = dynamic_cast<AliMagF*>(TGeoGlobalMagField::Instance()->GetField());
045be90c 86 if (field == 0) {
66b89079 87 AliError("Cannot get magnetic field from TGeoGlobalMagField");
045be90c 88 }
89 else {
90 Float_t Bfield = field->SolenoidField();
91 // NB: spatial shift has opposite sign for lay 5 and 6, but strip numbering also changes direction, so no sign-change
92 // Shift due to ExB on drift N-side, units: strip width
93 fLorentzShiftP = -repa->GetTanLorentzAngleElectronsSSD() * 150.e-4/95.e-4 * Bfield / 5.0;
94 // Shift due to ExB on drift P-side, units: strip width
95 fLorentzShiftN = -repa->GetTanLorentzAngleHolesSSD() * 150.e-4/95.e-4 * Bfield / 5.0;
96 AliDebug(1,Form("Bfield %f Lorentz Shift P-side %f N-side %f",Bfield,fLorentzShiftN,fLorentzShiftP));
97 }
66b89079 98 }
04366a57 99}
100
308b5ea4 101//______________________________________________________________________
66b89079 102AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(const AliITSClusterFinderV2SSD &cf) : AliITSClusterFinder(cf), fLastSSD1(cf.fLastSSD1), fLorentzShiftP(cf.fLorentzShiftP), fLorentzShiftN(cf.fLorentzShiftN)
308b5ea4 103{
104 // Copy constructor
308b5ea4 105}
106
107//______________________________________________________________________
108AliITSClusterFinderV2SSD& AliITSClusterFinderV2SSD::operator=(const AliITSClusterFinderV2SSD& cf ){
109 // Assignment operator
110
111 this->~AliITSClusterFinderV2SSD();
112 new(this) AliITSClusterFinderV2SSD(cf);
113 return *this;
114}
115
04366a57 116
117void AliITSClusterFinderV2SSD::FindRawClusters(Int_t mod){
118
119 //Find clusters V2
120 SetModule(mod);
121 FindClustersSSD(fDigits);
122
123}
124
125void AliITSClusterFinderV2SSD::FindClustersSSD(TClonesArray *alldigits) {
126 //------------------------------------------------------------
127 // Actual SSD cluster finder
128 //------------------------------------------------------------
d695268b 129 Int_t smaxall=alldigits->GetEntriesFast();
130 if (smaxall==0) return;
a86176e3 131
d695268b 132
133 //---------------------------------------
134 // load recoparam and calibration
135 //
136 static AliITSRecoParam *repa = NULL;
a86176e3 137 if(!repa){
138 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
139 if(!repa){
ed446fa3 140 repa = AliITSRecoParam::GetHighFluxParam();
a86176e3 141 AliWarning("Using default AliITSRecoParam class");
142 }
143 }
144
3a4139a2 145 AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)GetResp(fModule);
146 Float_t gain=0;
d695268b 147 Float_t noise=0;
148 //---------------------------------------
3a4139a2 149
d695268b 150
151 //------------------------------------
152 // fill the digits array with zero-suppression condition
153 // Signal is converted in KeV
154 //
f7b30404 155 TObjArray digits;
04366a57 156 for (Int_t i=0;i<smaxall; i++){
157 AliITSdigitSSD *d=(AliITSdigitSSD*)alldigits->UncheckedAt(i);
3a4139a2 158
d695268b 159 if(d->IsSideP()) noise = cal->GetNoiseP(d->GetStripNumber());
160 else noise = cal->GetNoiseN(d->GetStripNumber());
161 if (d->GetSignal()<3.*noise) continue;
162
3a4139a2 163 if(d->IsSideP()) gain = cal->GetGainP(d->GetStripNumber());
164 else gain = cal->GetGainN(d->GetStripNumber());
165
d695268b 166 Float_t q=gain*d->GetSignal(); //
3a4139a2 167 q=cal->ADCToKeV(q); // converts the charge in KeV from ADC units
6490e1c5 168 d->SetSignal(Int_t(q));
3a4139a2 169
f7b30404 170 digits.AddLast(d);
04366a57 171 }
f7b30404 172 Int_t smax = digits.GetEntriesFast();
04366a57 173 if (smax==0) return;
d695268b 174 //------------------------------------
175
04366a57 176
177 const Int_t kMax=1000;
178 Int_t np=0, nn=0;
179 Ali1Dcluster pos[kMax], neg[kMax];
180 Float_t y=0., q=0., qmax=0.;
d695268b 181 Int_t lab[4]={-2,-2,-2,-2};
182 Bool_t flag5 = 0;
04366a57 183
d695268b 184 /*
185 cout<<"-----------------------------"<<endl;
186 cout<<"this is module "<<fModule;
187 cout<<endl;
188 cout<<endl;
189 */
66b89079 190 Int_t layer = 4;
191 if (fModule>fLastSSD1)
192 layer = 5;
d695268b 193
194 //--------------------------------------------------------
195 // start 1D-clustering from the first digit in the digits array
196 //
f7b30404 197 AliITSdigitSSD *d=(AliITSdigitSSD*)digits.UncheckedAt(0);
04366a57 198 q += d->GetSignal();
199 y += d->GetCoord2()*d->GetSignal();
200 qmax=d->GetSignal();
201 lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
d695268b 202
203 if(d->IsSideP()) {
204 noise = cal->GetNoiseP(d->GetStripNumber());
205 gain = cal->GetGainP(d->GetStripNumber());
206 }
207 else {
208 noise = cal->GetNoiseN(d->GetStripNumber());
209 gain = cal->GetGainN(d->GetStripNumber());
210 }
211 noise*=gain;
212 noise=cal->ADCToKeV(noise); // converts noise in KeV from ADC units
213
214 if(qmax>fgkThreshold*noise) flag5=1; // seed for the cluster
215
216 /*
217 cout<<d->GetSignal()<<" "<<noise<<" "<<flag5<<" "<<
218 d->GetCoord1()<<" "<<d->GetCoord2()<<endl;
219 */
220
04366a57 221 Int_t curr=d->GetCoord2();
222 Int_t flag=d->GetCoord1();
d695268b 223
224 // Note: the first side which will be processed is supposed to be the
225 // P-side which is neg
04366a57 226 Int_t *n=&nn;
227 Ali1Dcluster *c=neg;
d695268b 228 if(flag) {n=&np; c=pos;} // in case we have only Nstrips (P was bad!)
229
04366a57 230 Int_t nd=1;
231 Int_t milab[10];
232 for (Int_t ilab=0;ilab<10;ilab++){
233 milab[ilab]=-2;
234 }
235 milab[0]=d->GetTrack(0); milab[1]=d->GetTrack(1); milab[2]=d->GetTrack(2);
236
d695268b 237
238 //----------------------------------------------------------
239 // search for neighboring digits
240 //
04366a57 241 for (Int_t s=1; s<smax; s++) {
f7b30404 242 d=(AliITSdigitSSD*)digits.UncheckedAt(s);
04366a57 243 Int_t strip=d->GetCoord2();
d695268b 244
245 // if digits is not a neighbour or side did not change
246 // and at least one of the previous digits met the seed condition
247 // then creates a new 1D cluster
248 if ( ( ((strip-curr) > 1) || (flag!=d->GetCoord1()) ) ) {
249
250 if(flag5) {
251 //cout<<"here1"<<endl;
66b89079 252 Float_t dLorentz = 0;
253 if (!flag) { // P-side is neg clust
254 dLorentz = fLorentzShiftN;
255 }
256 else { // N-side is p clust
257 dLorentz = fLorentzShiftP;
258 }
259 c[*n].SetY(y/q+dLorentz);
04366a57 260 c[*n].SetQ(q);
261 c[*n].SetNd(nd);
262 CheckLabels2(milab);
263 c[*n].SetLabels(milab);
a86176e3 264
265 if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
d695268b 266 // Note: fUseUnfoldingInClusterFinderSSD=kFALSE by default in RecoParam
267
a86176e3 268 //Split suspiciously big cluster
269 if (nd>4&&nd<25) {
66b89079 270 c[*n].SetY(y/q-0.25*nd+dLorentz);
a86176e3 271 c[*n].SetQ(0.5*q);
272 (*n)++;
273 if (*n==kMax) {
274 Error("FindClustersSSD","Too many 1D clusters !");
275 return;
276 }
66b89079 277 c[*n].SetY(y/q+0.25*nd+dLorentz);
a86176e3 278 c[*n].SetQ(0.5*q);
279 c[*n].SetNd(nd);
280 c[*n].SetLabels(milab);
281 }
282
283 } // unfolding is on
284
04366a57 285 (*n)++;
286 if (*n==kMax) {
287 Error("FindClustersSSD","Too many 1D clusters !");
288 return;
289 }
d695268b 290
291 } // flag5 set
292
293 // reset everything
04366a57 294 y=q=qmax=0.;
295 nd=0;
d695268b 296 flag5=0;
04366a57 297 lab[0]=lab[1]=lab[2]=-2;
d695268b 298 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
299
300 // if side changed from P to N, switch to pos 1D clusters
301 // (if for some reason the side changed from N to P then do the opposite)
302 if (flag!=d->GetCoord1())
303 { if(!flag) {n=&np; c=pos;} else {n=&nn; c=neg;} }
304
305 } // end create new 1D cluster from previous neighboring digits
306
307 // continues adding digits to the previous cluster
308 // or start a new one
04366a57 309 flag=d->GetCoord1();
310 q += d->GetSignal();
311 y += d->GetCoord2()*d->GetSignal();
312 nd++;
d695268b 313
314 if(d->IsSideP()) {
315 noise = cal->GetNoiseP(d->GetStripNumber());
316 gain = cal->GetGainP(d->GetStripNumber());
317 }
318 else {
319 noise = cal->GetNoiseN(d->GetStripNumber());
320 gain = cal->GetGainN(d->GetStripNumber());
321 }
322 noise*=gain;
323 noise=cal->ADCToKeV(noise); // converts the charge in KeV from ADC units
324
325 if(d->GetSignal()>fgkThreshold*noise) flag5=1;
326
327 /*
328 cout<<d->GetSignal()<<" "<<noise<<" "<<flag5<<" "<<
329 d->GetCoord1()<<" "<<d->GetCoord2()<<endl;
330 */
331
04366a57 332 if (d->GetSignal()>qmax) {
333 qmax=d->GetSignal();
334 lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
335 }
336 for (Int_t ilab=0;ilab<10;ilab++) {
337 if (d->GetTrack(ilab)>=0) AddLabel(milab, (d->GetTrack(ilab)));
338 }
339 curr=strip;
a86176e3 340
d695268b 341
342 } // loop over digits, no more digits in the digits array
343
344
345 // add the last 1D cluster
346 if(flag5) {
347
348 // cout<<"here2"<<endl;
66b89079 349 Float_t dLorentz = 0;
350 if (!flag) { // P-side is neg clust
351 dLorentz = fLorentzShiftN;
352 }
353 else { // N-side is p clust
354 dLorentz = fLorentzShiftP;
355 }
356
357 c[*n].SetY(y/q + dLorentz);
d695268b 358 c[*n].SetQ(q);
359 c[*n].SetNd(nd);
360 c[*n].SetLabels(lab);
a86176e3 361
d695268b 362 if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
363
364 //Split suspiciously big cluster
365 if (nd>4 && nd<25) {
66b89079 366 c[*n].SetY(y/q-0.25*nd + dLorentz);
d695268b 367 c[*n].SetQ(0.5*q);
368 (*n)++;
369 if (*n==kMax) {
370 Error("FindClustersSSD","Too many 1D clusters !");
371 return;
372 }
66b89079 373 c[*n].SetY(y/q+0.25*nd + dLorentz);
d695268b 374 c[*n].SetQ(0.5*q);
375 c[*n].SetNd(nd);
376 c[*n].SetLabels(lab);
a86176e3 377 }
d695268b 378 } // unfolding is on
379
380 (*n)++;
381 if (*n==kMax) {
382 Error("FindClustersSSD","Too many 1D clusters !");
383 return;
a86176e3 384 }
04366a57 385
d695268b 386 } // if flag5 last 1D cluster added
387
388
389 //------------------------------------------------------
390 // call FindClustersSSD to pair neg and pos 1D clusters
391 // and create recpoints from the crosses
392 // Note1: neg are Pside and pos are Nside!!
393 // Note2: if there are no Pside digits nn=0 (bad strips!!) (same for Nside)
394 //
395 // cout<<nn<<" Pside and "<<np<<" Nside clusters"<<endl;
138df073 396
397 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
398 if (nn*np > 0) {
399 TClonesArray* clusters = rpc->UncheckedGetClusters(fModule);
400 clusters->Clear();
401 FindClustersSSD(neg, nn, pos, np, clusters);
402 TIter itr(clusters);
403 AliITSRecPoint *irp;
404 while ((irp = (AliITSRecPoint*)itr.Next())) fDetTypeRec->AddRecPoint(*irp);
405 }
d695268b 406 //-----------------------------------------------------
04366a57 407}
408
409
01ef1bd4 410void AliITSClusterFinderV2SSD::RawdataToClusters(AliRawReader* rawReader){
04366a57 411
66b89079 412 //------------------------------------------------------------
04366a57 413 // This function creates ITS clusters from raw data
414 //------------------------------------------------------------
415 rawReader->Reset();
416 AliITSRawStreamSSD inputSSD(rawReader);
01ef1bd4 417 FindClustersSSD(&inputSSD);
04366a57 418
419}
420
138df073 421
01ef1bd4 422void AliITSClusterFinderV2SSD::FindClustersSSD(AliITSRawStreamSSD* input)
04366a57 423{
424 //------------------------------------------------------------
425 // Actual SSD cluster finder for raw data
426 //------------------------------------------------------------
a86176e3 427
01ef1bd4 428 AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
a86176e3 429 static AliITSRecoParam *repa = NULL;
430 if(!repa){
431 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
432 if(!repa){
ed446fa3 433 repa = AliITSRecoParam::GetHighFluxParam();
a86176e3 434 AliWarning("Using default AliITSRecoParam class");
435 }
436 }
04366a57 437 Int_t nClustersSSD = 0;
138df073 438 const Int_t kNADC = 12;
439 const Int_t kMaxADCClusters = 1000;
440
441 Int_t strips[kNADC][2][kMaxADCClusters][2]; // [ADC],[side],[istrip], [0]=istrip [1]=signal
442 Int_t nStrips[kNADC][2];
443
444 for( int i=0; i<kNADC; i++ ){
445 nStrips[i][0] = 0;
446 nStrips[i][1] = 0;
447 }
448
449 Int_t ddl = -1;
450 Int_t ad = -1;
451
452 //*
453 //* Loop over modules DDL+AD
454 //*
455
04366a57 456 while (kTRUE) {
138df073 457
458 bool next = input->Next();
308b5ea4 459
138df073 460 //*
461 //* Continue if corrupted input
462 //*
463
464 if( (!next)&&(input->flag) ){
465 AliWarning(Form("ITSClustersFinderSSD: Corrupted data: warning from RawReader"));
466 continue;
308b5ea4 467 }
138df073 468
469 Int_t newDDL = input->GetDDL();
470 Int_t newAD = input->GetAD();
471
472 if( next ){
473 if( newDDL<0 || newDDL>15 ){
474 AliWarning(Form("ITSClustersFinderSSD: Corrupted data: wrong DDL number (%d)",newDDL));
475 continue;
476 }
308b5ea4 477
138df073 478 if( newAD<1 || newAD>9 ){
479 AliWarning(Form("ITSClustersFinderSSD: Corrupted data: wrong AD number (%d)",newAD));
480 continue;
481 }
308b5ea4 482 }
42ed6062 483
138df073 484 bool newModule = ( !next || ddl!= newDDL || ad!=newAD );
42ed6062 485
138df073 486 if( newModule && ddl>=0 && ad>=0 ){
bc4dd89a 487
138df073 488 //*
489 //* Reconstruct the previous block of 12 modules --- actual clusterfinder
490 //*
491 //cout<<endl;
492 for( int adc = 0; adc<kNADC; adc++ ){
b42cfa25 493
138df073 494 //* 1D clusterfinder
bc4dd89a 495
138df073 496 Ali1Dcluster clusters1D[2][kMaxADCClusters]; // per ADC, per side
497 Int_t nClusters1D[2] = {0,0};
498 //int nstat[2] = {0,0};
499 fModule = AliITSRawStreamSSD::GetModuleNumber(ddl, (ad - 1) * 12 + adc );
b42cfa25 500
138df073 501 if( fModule<0 ){
502// AliWarning(Form("ITSClustersFinderSSD: Corrupted data: module (ddl %d ad %d adc %d) not found in the map",ddl,ad,adc));
503//CM channels are always present even everything is suppressed
504 continue;
505 }
66b89079 506
507 Int_t layer = 4;
508 if (fModule>fLastSSD1)
509 layer = 5;
a64f9843 510
138df073 511 AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)fDetTypeRec->GetCalibrationModel(fModule);
512 if( !cal ){
513 AliWarning(Form("ITSClustersFinderSSD: No calibration found for module (ddl %d ad %d adc %d)",ddl,ad,adc));
514 continue;
515 }
308b5ea4 516
138df073 517 Float_t dStrip = 0;
518
519 if( repa->GetUseCosmicRunShiftsSSD()) { // Special condition for 2007/2008 cosmic data
520 dStrip = fgkCosmic2008StripShifts[ddl][ad-1];
521 if (TMath::Abs(dStrip) > 1.5){
522 AliWarning(Form("Indexing error in Cosmic calibration: ddl = %d, dStrip %f\n",ddl,dStrip));
523 dStrip = 0;
524 }
308b5ea4 525 }
b42cfa25 526
138df073 527 for( int side=0; side<=1; side++ ){
528
529 Int_t lab[3]={-2,-2,-2};
530 Float_t q = 0.;
531 Float_t y = 0.;
532 Int_t nDigits = 0;
533 Int_t ostrip = -2;
534 Bool_t snFlag = 0;
66b89079 535
536 Float_t dLorentz = 0;
537 if (side==0) { // P-side is neg clust
538 dLorentz = fLorentzShiftN;
539 }
540 else { // N-side is pos clust
541 dLorentz = fLorentzShiftP;
542 }
b42cfa25 543
138df073 544 Int_t n = nStrips[adc][side];
545 for( int istr = 0; istr<n+1; istr++ ){
b42cfa25 546
138df073 547 bool stripOK = 1;
548 Int_t strip=0;
549 Float_t signal=0.0, noise=0.0, gain=0.0;
308b5ea4 550
138df073 551 if( istr<n ){
552 strip = strips[adc][side][istr][0];
553 signal = strips[adc][side][istr][1];
a86176e3 554
138df073 555 //cout<<"strip "<<adc<<" / "<<side<<": "<<strip<<endl;
556
557 if( cal ){
558 noise = side ?cal->GetNoiseN(strip) :cal->GetNoiseP(strip);
559 gain = side ?cal->GetGainN(strip) :cal->GetGainP(strip);
560 stripOK = ( noise>=1. && signal>=3.0*noise
561 //&& !cal->IsPChannelBad(strip)
562 );
563 }
564 } else stripOK = 0; // end of data
565
566 bool newCluster = ( (abs(strip-ostrip)!=1) || !stripOK );
567
568 if( newCluster ){
569
570 //* Store the previous cluster
571
572 if( nDigits>0 && q>0 && snFlag ){
573
574 if (nClusters1D[side] >= kMaxADCClusters-1 ) {
575 AliWarning("HLT ClustersFinderSSD: Too many 1D clusters !");
576 }else {
577
578 Ali1Dcluster &cluster = clusters1D[side][nClusters1D[side]++];
66b89079 579 cluster.SetY( y / q + dStrip + dLorentz);
138df073 580 cluster.SetQ(q);
581 cluster.SetNd(nDigits);
582 cluster.SetLabels(lab);
583 //cout<<"cluster 1D side "<<side<<": y= "<<y<<" q= "<<q<<" d="<<dStrip<<" Y="<<cluster.GetY()<<endl;
584 //Split suspiciously big cluster
585
586 if( repa->GetUseUnfoldingInClusterFinderSSD()
587 && nDigits > 4 && nDigits < 25
588 ){
66b89079 589 cluster.SetY(y/q + dStrip - 0.25*nDigits + dLorentz);
138df073 590 cluster.SetQ(0.5*q);
591 Ali1Dcluster& cluster2 = clusters1D[side][nClusters1D[side]++];
66b89079 592 cluster2.SetY(y/q + dStrip + 0.25*nDigits + dLorentz);
138df073 593 cluster2.SetQ(0.5*q);
594 cluster2.SetNd(nDigits);
595 cluster2.SetLabels(lab);
596 } // unfolding is on
a86176e3 597 }
308b5ea4 598 }
138df073 599 y = q = 0.;
600 nDigits = 0;
601 snFlag = 0;
602
603 } //* End store the previous cluster
604
605 if( stripOK ){ // add new signal to the cluster
606// signal = (Int_t) (signal * gain); // signal is corrected for gain
607 if( signal>fgkThreshold*noise) snFlag = 1;
608 signal = signal * gain; // signal is corrected for gain
609// if( cal ) signal = (Int_t) cal->ADCToKeV( signal ); // signal is converted in KeV
610 if( cal ) signal = cal->ADCToKeV( signal ); // signal is converted in KeV
611 q += signal; // add digit to current cluster
612 y += strip * signal;
613 nDigits++;
614 //nstat[side]++;
615 ostrip = strip;
a64f9843 616
308b5ea4 617 }
138df073 618 } //* end loop over strips
b42cfa25 619
138df073 620 } //* end loop over ADC sides
308b5ea4 621
a64f9843 622
138df073 623 //* 2D clusterfinder
624 if( nClusters1D[0] && nClusters1D[1] && fModule>=0 ){
625 TClonesArray* clusters = rpc->UncheckedGetClusters(fModule);
626 FindClustersSSD( clusters1D[0], nClusters1D[0], clusters1D[1], nClusters1D[1], clusters);
627 Int_t nClustersn = clusters->GetEntriesFast();
628 nClustersSSD += nClustersn;
308b5ea4 629 }
630
138df073 631 //cout<<"SG: "<<ddl<<" "<<ad<<" "<<adc<<": strips "<<nstat[0]<<"+"<<nstat[1]<<", clusters 1D= "<<nClusters1D[0]<<" + "<<nClusters1D[1]<<", 2D= "<<clusters.size()<<endl;
a64f9843 632
138df073 633 }//* end loop over adc
634
635 }//* end of reconstruction of previous block of 12 modules
636
637 if( newModule ){
638
639 //*
640 //* Clean up arrays and set new module
641 //*
642
643 for( int i=0; i<kNADC; i++ ){
644 nStrips[i][0] = 0;
645 nStrips[i][1] = 0;
646 }
647 ddl = newDDL;
648 ad = newAD;
649 }
650
a64f9843 651
138df073 652 //*
653 //* Exit main loop when there is no more input
654 //*
04366a57 655
138df073 656 if( !next ) break;
657
658 //*
659 //* Fill the current strip information
660 //*
661
662 Int_t adc = input->GetADC();
663 if( adc<0 || adc>=kNADC+2 || (adc>5&&adc<8) ){
664 AliWarning(Form("HLT ClustersFinderSSD: Corrupted data: wrong adc number (%d)", adc));
665 continue;
666 }
a64f9843 667
138df073 668 if( adc>7 ) adc-= 2; // shift ADC numbers 8-13 to 6-11
669
670 Bool_t side = input->GetSideFlag();
671 Int_t strip = input->GetStrip();
672 Int_t signal = input->GetSignal();
673
a86176e3 674
138df073 675 //cout<<"SSD: "<<ddl<<" "<<ad<<" "<<adc<<" "<<side<<" "<<strip<<" : "<<signal<<endl;
04366a57 676
138df073 677 if( strip>767 ){
678 AliWarning(Form("HLT ClustersFinderSSD: Corrupted data: wrong strip number (ddl %d ad %d adc %d side %d, strip %d",
679 ddl, ad, adc, side,strip) );
680 continue;
681 }
682 if (strip < 0) continue;
683
684 int &n = nStrips[adc][side];
685 if( n >0 ){
686 Int_t oldStrip = strips[adc][side][n-1][0];
687
688 if( strip==oldStrip ){
689 AliWarning(Form("HLT ClustersFinderSSD: Corrupted data: duplicated signal: ddl %d ad %d adc %d, side %d, strip %d",
690 ddl, ad, adc, side, strip ));
691 continue;
692 }
693 }
694 strips[adc][side][n][0] = strip;
695 strips[adc][side][n][1] = signal;
696 n++;
3a4139a2 697
138df073 698 //cout<<"SSD: "<<input->GetDDL()<<" "<<input->GetAD()<<" "
699 //<<input->GetADC()<<" "<<input->GetSideFlag()<<" "<<((int)input->GetStrip())<<" "<<strip<<" : "<<input->GetSignal()<<endl;
04366a57 700
138df073 701 } //* End main loop over the input
308b5ea4 702
3b1d8321 703 AliDebug(1,Form("found clusters in ITS SSD: %d", nClustersSSD));
04366a57 704}
705
138df073 706
04366a57 707void AliITSClusterFinderV2SSD::
708FindClustersSSD(Ali1Dcluster* neg, Int_t nn,
709 Ali1Dcluster* pos, Int_t np,
710 TClonesArray *clusters) {
711 //------------------------------------------------------------
712 // Actual SSD cluster finder
713 //------------------------------------------------------------
b4704be3 714
715 const TGeoHMatrix *mT2L=AliITSgeomTGeo::GetTracking2LocalMatrix(fModule);
716
7101948c 717 //---------------------------------------
718 // load recoparam
719 //
720 static AliITSRecoParam *repa = NULL;
721 if(!repa){
722 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
723 if(!repa){
724 repa = AliITSRecoParam::GetHighFluxParam();
725 AliWarning("Using default AliITSRecoParam class");
726 }
727 }
728
138df073 729// TClonesArray &cl=*clusters;
d036ccd3 730
443e42aa 731 AliITSsegmentationSSD *seg = static_cast<AliITSsegmentationSSD*>(fDetTypeRec->GetSegmentationModel(2));
8be4e1b1 732 if (fModule>fLastSSD1)
733 seg->SetLayer(6);
734 else
735 seg->SetLayer(5);
736
737 Float_t hwSSD = seg->Dx()*1e-4/2;
738 Float_t hlSSD = seg->Dz()*1e-4/2;
739
04366a57 740 Int_t idet=fNdet[fModule];
741 Int_t ncl=0;
d036ccd3 742
04366a57 743 //
8be4e1b1 744 Int_t *cnegative = new Int_t[np];
745 Int_t *cused1 = new Int_t[np];
746 Int_t *negativepair = new Int_t[10*np];
747 Int_t *cpositive = new Int_t[nn];
748 Int_t *cused2 = new Int_t[nn];
749 Int_t *positivepair = new Int_t[10*nn];
750 for (Int_t i=0;i<np;i++) {cnegative[i]=0; cused1[i]=0;}
751 for (Int_t i=0;i<nn;i++) {cpositive[i]=0; cused2[i]=0;}
752 for (Int_t i=0;i<10*np;i++) {negativepair[i]=0;}
753 for (Int_t i=0;i<10*nn;i++) {positivepair[i]=0;}
c157c94e 754
308b5ea4 755 if ((np*nn) > fgPairsSize) {
d036ccd3 756
308b5ea4 757 if (fgPairs) delete [] fgPairs;
758 fgPairsSize = 4*np*nn;
759 fgPairs = new Short_t[fgPairsSize];
c157c94e 760 }
308b5ea4 761 memset(fgPairs,0,sizeof(Short_t)*np*nn);
762
04366a57 763 //
764 // find available pairs
765 //
138df073 766 Int_t ncross = 0;
04366a57 767 for (Int_t i=0; i<np; i++) {
d036ccd3 768 Float_t yp=pos[i].GetY();
a64f9843 769 if ( (pos[i].GetQ()>0) && (pos[i].GetQ()<3) ) continue;
04366a57 770 for (Int_t j=0; j<nn; j++) {
a64f9843 771 if ( (neg[j].GetQ()>0) && (neg[j].GetQ()<3) ) continue;
d036ccd3 772 Float_t yn=neg[j].GetY();
773
8be4e1b1 774 Float_t xt, zt;
775 seg->GetPadCxz(yn, yp, xt, zt);
d695268b 776 //cout<<yn<<" "<<yp<<" "<<xt<<" "<<zt<<endl;
d036ccd3 777
5846520e 778 if (TMath::Abs(xt)<hwSSD)
779 if (TMath::Abs(zt)<hlSSD) {
710f576f 780 Int_t in = i*10+cnegative[i];
781 Int_t ip = j*10+cpositive[j];
782 if ((in < 10*np) && (ip < 10*nn)) {
783 negativepair[in] =j; //index
784 positivepair[ip] =i;
785 cnegative[i]++; //counters
138df073 786 cpositive[j]++;
787 ncross++;
710f576f 788 fgPairs[i*nn+j]=100;
789 }
790 else
791 AliError(Form("Index out of range: ip=%d, in=%d",ip,in));
04366a57 792 }
793 }
794 }
308b5ea4 795
eb595a70 796 if (!ncross) {
797 delete [] cnegative;
798 delete [] cused1;
799 delete [] negativepair;
800 delete [] cpositive;
801 delete [] cused2;
802 delete [] positivepair;
803 return;
804 }
138df073 805//why not to allocate memorey here? if(!clusters) clusters = new TClonesArray("AliITSRecPoint", ncross);
806
5846520e 807 /* //
308b5ea4 808 // try to recover points out of but close to the module boundaries
04366a57 809 //
810 for (Int_t i=0; i<np; i++) {
d036ccd3 811 Float_t yp=pos[i].GetY();
a64f9843 812 if ( (pos[i].GetQ()>0) && (pos[i].GetQ()<3) ) continue;
04366a57 813 for (Int_t j=0; j<nn; j++) {
a64f9843 814 if ( (neg[j].GetQ()>0) && (neg[j].GetQ()<3) ) continue;
308b5ea4 815 // if both 1Dclusters have an other cross continue
04366a57 816 if (cpositive[j]&&cnegative[i]) continue;
8be4e1b1 817 Float_t yn=neg[j].GetY();
d036ccd3 818
8be4e1b1 819 Float_t xt, zt;
820 seg->GetPadCxz(yn, yp, xt, zt);
d036ccd3 821
8be4e1b1 822 if (TMath::Abs(xt)<hwSSD+0.1)
823 if (TMath::Abs(zt)<hlSSD+0.15) {
308b5ea4 824 // tag 1Dcluster (eventually will produce low quality recpoint)
04366a57 825 if (cnegative[i]==0) pos[i].SetNd(100); // not available pair
826 if (cpositive[j]==0) neg[j].SetNd(100); // not available pair
710f576f 827 Int_t in = i*10+cnegative[i];
828 Int_t ip = j*10+cpositive[j];
829 if ((in < 10*np) && (ip < 10*nn)) {
830 negativepair[in] =j; //index
831 positivepair[ip] =i;
832 cnegative[i]++; //counters
833 cpositive[j]++;
834 fgPairs[i*nn+j]=100;
835 }
836 else
837 AliError(Form("Index out of range: ip=%d, in=%d",ip,in));
04366a57 838 }
839 }
840 }
5846520e 841 */
308b5ea4 842
04366a57 843 //
d695268b 844 Float_t lp[6];
04366a57 845 Int_t milab[10];
846 Double_t ratio;
847
b42cfa25 848
a86176e3 849 if(repa->GetUseChargeMatchingInClusterFinderSSD()==kTRUE) {
308b5ea4 850
851
a86176e3 852 //
853 // sign gold tracks
854 //
855 for (Int_t ip=0;ip<np;ip++){
8be4e1b1 856 Float_t xbest=1000,zbest=1000,qbest=0;
a86176e3 857 //
858 // select gold clusters
859 if ( (cnegative[ip]==1) && cpositive[negativepair[10*ip]]==1){
d036ccd3 860 Float_t yp=pos[ip].GetY();
a86176e3 861 Int_t j = negativepair[10*ip];
a64f9843 862
863 if( (pos[ip].GetQ()==0) && (neg[j].GetQ() ==0) ) {
864 // both bad, hence continue;
865 // mark both as used (to avoid recover at the end)
866 cused1[ip]++;
867 cused2[j]++;
868 continue;
869 }
870
a86176e3 871 ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
d695268b 872 //cout<<"ratio="<<ratio<<endl;
a64f9843 873
874 // charge matching (note that if posQ or negQ is 0 -> ratio=1 and the following condition is met
2069484c 875 if (TMath::Abs(ratio)>0.2) continue; // note: 0.2=3xsigma_ratio calculated in cosmics tests
a64f9843 876
a86176e3 877 //
8be4e1b1 878 Float_t yn=neg[j].GetY();
879
880 Float_t xt, zt;
881 seg->GetPadCxz(yn, yp, xt, zt);
882
883 xbest=xt; zbest=zt;
2069484c 884
2069484c 885
a86176e3 886 qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
d695268b 887 if( (pos[ip].GetQ()==0)||(neg[j].GetQ()==0)) qbest*=2; // in case of bad strips on one side keep all charge from the other one
a86176e3 888
a86176e3 889 {
8be4e1b1 890 Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
a86176e3 891 mT2L->MasterToLocal(loc,trk);
892 lp[0]=trk[1];
893 lp[1]=trk[2];
00a7cc50 894 }
a86176e3 895 lp[4]=qbest; //Q
896 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
897 for (Int_t ilab=0;ilab<3;ilab++){
898 milab[ilab] = pos[ip].GetLabel(ilab);
899 milab[ilab+3] = neg[j].GetLabel(ilab);
00a7cc50 900 }
04366a57 901 //
a86176e3 902 CheckLabels2(milab);
903 milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
904 Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
d695268b 905
7eb157d7 906 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
907 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
d695268b 908 // out-of-diagonal element of covariance matrix
909 if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
910 else if ( (info[0]>1) && (info[1]>1) ) {
7eb157d7 911 lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2
912 lp[3]=0.0065; // 0.08*0.08; //SigmaZ2
913 lp[5]=-6.48e-05;
d695268b 914 }
915 else {
7eb157d7 916 lp[2]=4.80e-06; // 0.00219*0.00219
917 lp[3]=0.0093; // 0.0964*0.0964;
918 if (info[0]==1) {
919 lp[5]=-0.00014;
920 }
921 else {
922 lp[2]=2.79e-06; // 0.0017*0.0017;
923 lp[3]=0.00935; // 0.967*0.967;
924 lp[5]=-4.32e-05;
925 }
d695268b 926 }
927
a86176e3 928 AliITSRecPoint * cl2;
138df073 929 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
a86176e3 930
138df073 931 cl2->SetChargeRatio(ratio);
932 cl2->SetType(1);
933 fgPairs[ip*nn+j]=1;
a64f9843 934
138df073 935 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
936 cl2->SetType(2);
937 fgPairs[ip*nn+j]=2;
938 }
a64f9843 939
138df073 940 if(pos[ip].GetQ()==0) cl2->SetType(3);
941 if(neg[j].GetQ()==0) cl2->SetType(4);
a64f9843 942
138df073 943 cused1[ip]++;
944 cused2[j]++;
a64f9843 945
138df073 946 ncl++;
a86176e3 947 }
948 }
949
950 for (Int_t ip=0;ip<np;ip++){
8be4e1b1 951 Float_t xbest=1000,zbest=1000,qbest=0;
a86176e3 952 //
953 //
954 // select "silber" cluster
955 if ( cnegative[ip]==1 && cpositive[negativepair[10*ip]]==2){
956 Int_t in = negativepair[10*ip];
957 Int_t ip2 = positivepair[10*in];
958 if (ip2==ip) ip2 = positivepair[10*in+1];
959 Float_t pcharge = pos[ip].GetQ()+pos[ip2].GetQ();
a64f9843 960
d695268b 961
962
963 ratio = (pcharge-neg[in].GetQ())/(pcharge+neg[in].GetQ());
964 if ( (TMath::Abs(ratio)<0.2) && (pcharge!=0) ) {
965 //if ( (TMath::Abs(pcharge-neg[in].GetQ())<30) && (pcharge!=0) ) { //
a64f9843 966
a86176e3 967 //
968 // add first pair
a64f9843 969 if ( (fgPairs[ip*nn+in]==100)&&(pos[ip].GetQ() ) ) { //
970
d036ccd3 971 Float_t yp=pos[ip].GetY();
8be4e1b1 972 Float_t yn=neg[in].GetY();
973
974 Float_t xt, zt;
975 seg->GetPadCxz(yn, yp, xt, zt);
976
977 xbest=xt; zbest=zt;
2069484c 978
a86176e3 979 qbest =pos[ip].GetQ();
8be4e1b1 980 Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
a64f9843 981 mT2L->MasterToLocal(loc,trk);
982 lp[0]=trk[1];
983 lp[1]=trk[2];
984
a86176e3 985 lp[4]=qbest; //Q
986 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
987 for (Int_t ilab=0;ilab<3;ilab++){
988 milab[ilab] = pos[ip].GetLabel(ilab);
989 milab[ilab+3] = neg[in].GetLabel(ilab);
990 }
991 //
992 CheckLabels2(milab);
993 ratio = (pos[ip].GetQ()-neg[in].GetQ())/(pos[ip].GetQ()+neg[in].GetQ());
994 milab[3]=(((ip<<10) + in)<<10) + idet; // pos|neg|det
995 Int_t info[3] = {pos[ip].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
996
7eb157d7 997 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
998 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
d695268b 999 // out-of-diagonal element of covariance matrix
1000 if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1001 else if ( (info[0]>1) && (info[1]>1) ) {
7eb157d7 1002 lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2
1003 lp[3]=0.0065; // 0.08*0.08; //SigmaZ2
1004 lp[5]=-6.48e-05;
d695268b 1005 }
1006 else {
7eb157d7 1007 lp[2]=4.80e-06; // 0.00219*0.00219
1008 lp[3]=0.0093; // 0.0964*0.0964;
1009 if (info[0]==1) {
1010 lp[5]=-0.00014;
1011 }
1012 else {
1013 lp[2]=2.79e-06; // 0.0017*0.0017;
1014 lp[3]=0.00935; // 0.967*0.967;
1015 lp[5]=-4.32e-05;
1016 }
d695268b 1017 }
1018
7eb157d7 1019 AliITSRecPoint * cl2;
138df073 1020 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
a86176e3 1021 cl2->SetChargeRatio(ratio);
1022 cl2->SetType(5);
1023 fgPairs[ip*nn+in] = 5;
1024 if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster
1025 cl2->SetType(6);
1026 fgPairs[ip*nn+in] = 6;
1027 }
a86176e3 1028 ncl++;
04366a57 1029 }
04366a57 1030
a64f9843 1031
04366a57 1032 //
a86176e3 1033 // add second pair
00a7cc50 1034
a86176e3 1035 // if (!(cused1[ip2] || cused2[in])){ //
a64f9843 1036 if ( (fgPairs[ip2*nn+in]==100) && (pos[ip2].GetQ()) ) {
1037
d036ccd3 1038 Float_t yp=pos[ip2].GetY();
8be4e1b1 1039 Float_t yn=neg[in].GetY();
1040
1041 Float_t xt, zt;
1042 seg->GetPadCxz(yn, yp, xt, zt);
1043
1044 xbest=xt; zbest=zt;
2069484c 1045
a86176e3 1046 qbest =pos[ip2].GetQ();
a64f9843 1047
8be4e1b1 1048 Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
a64f9843 1049 mT2L->MasterToLocal(loc,trk);
1050 lp[0]=trk[1];
1051 lp[1]=trk[2];
1052
a86176e3 1053 lp[4]=qbest; //Q
1054 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1055 for (Int_t ilab=0;ilab<3;ilab++){
1056 milab[ilab] = pos[ip2].GetLabel(ilab);
1057 milab[ilab+3] = neg[in].GetLabel(ilab);
00a7cc50 1058 }
a86176e3 1059 //
1060 CheckLabels2(milab);
1061 ratio = (pos[ip2].GetQ()-neg[in].GetQ())/(pos[ip2].GetQ()+neg[in].GetQ());
1062 milab[3]=(((ip2<<10) + in)<<10) + idet; // pos|neg|det
1063 Int_t info[3] = {pos[ip2].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
308b5ea4 1064
7eb157d7 1065 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
1066 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
d695268b 1067 // out-of-diagonal element of covariance matrix
1068 if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1069 else if ( (info[0]>1) && (info[1]>1) ) {
7eb157d7 1070 lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2
1071 lp[3]=0.0065; // 0.08*0.08; //SigmaZ2
1072 lp[5]=-6.48e-05;
d695268b 1073 }
1074 else {
7eb157d7 1075 lp[2]=4.80e-06; // 0.00219*0.00219
1076 lp[3]=0.0093; // 0.0964*0.0964;
1077 if (info[0]==1) {
1078 lp[5]=-0.00014;
1079 }
1080 else {
1081 lp[2]=2.79e-06; // 0.0017*0.0017;
1082 lp[3]=0.00935; // 0.967*0.967;
1083 lp[5]=-4.32e-05;
1084 }
d695268b 1085 }
1086
a86176e3 1087 AliITSRecPoint * cl2;
138df073 1088 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
a86176e3 1089
a86176e3 1090 cl2->SetChargeRatio(ratio);
1091 cl2->SetType(5);
1092 fgPairs[ip2*nn+in] =5;
1093 if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster
1094 cl2->SetType(6);
1095 fgPairs[ip2*nn+in] =6;
1096 }
a86176e3 1097 ncl++;
a64f9843 1098 }
1099
a86176e3 1100 cused1[ip]++;
1101 cused1[ip2]++;
1102 cused2[in]++;
a64f9843 1103
1104 } // charge matching condition
1105
1106 } // 2 Pside cross 1 Nside
1107 } // loop over Pside clusters
a86176e3 1108
a64f9843 1109
1110
1111 //
1112 for (Int_t jn=0;jn<nn;jn++){
1113 if (cused2[jn]) continue;
8be4e1b1 1114 Float_t xbest=1000,zbest=1000,qbest=0;
a64f9843 1115 // select "silber" cluster
1116 if ( cpositive[jn]==1 && cnegative[positivepair[10*jn]]==2){
1117 Int_t ip = positivepair[10*jn];
1118 Int_t jn2 = negativepair[10*ip];
1119 if (jn2==jn) jn2 = negativepair[10*ip+1];
1120 Float_t pcharge = neg[jn].GetQ()+neg[jn2].GetQ();
04366a57 1121 //
a64f9843 1122
d695268b 1123
1124 ratio = (pcharge-pos[ip].GetQ())/(pcharge+pos[ip].GetQ());
1125 if ( (TMath::Abs(ratio)<0.2) && (pcharge!=0) ) {
1126
1127 /*
2069484c 1128 if ( (TMath::Abs(pcharge-pos[ip].GetQ())<30) && // charge matching
a64f9843 1129 (pcharge!=0) ) { // reject combinations of bad strips
d695268b 1130 */
1131
1132
a64f9843 1133 //
1134 // add first pair
1135 // if (!(cused1[ip]||cused2[jn])){
1136 if ( (fgPairs[ip*nn+jn]==100) && (neg[jn].GetQ()) ) { //
1137
d036ccd3 1138 Float_t yn=neg[jn].GetY();
8be4e1b1 1139 Float_t yp=pos[ip].GetY();
2069484c 1140
8be4e1b1 1141 Float_t xt, zt;
1142 seg->GetPadCxz(yn, yp, xt, zt);
1143
1144 xbest=xt; zbest=zt;
2069484c 1145
a64f9843 1146 qbest =neg[jn].GetQ();
d036ccd3 1147
a64f9843 1148 {
8be4e1b1 1149 Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
a64f9843 1150 mT2L->MasterToLocal(loc,trk);
1151 lp[0]=trk[1];
1152 lp[1]=trk[2];
b4704be3 1153 }
04366a57 1154
1155 lp[4]=qbest; //Q
1156 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1157 for (Int_t ilab=0;ilab<3;ilab++){
1158 milab[ilab] = pos[ip].GetLabel(ilab);
1159 milab[ilab+3] = neg[jn].GetLabel(ilab);
1160 }
1161 //
1162 CheckLabels2(milab);
1163 ratio = (pos[ip].GetQ()-neg[jn].GetQ())/(pos[ip].GetQ()+neg[jn].GetQ());
1164 milab[3]=(((ip<<10) + jn)<<10) + idet; // pos|neg|det
1165 Int_t info[3] = {pos[ip].GetNd(),neg[jn].GetNd(),fNlayer[fModule]};
1166
7eb157d7 1167 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
1168 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
d695268b 1169 // out-of-diagonal element of covariance matrix
1170 if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1171 else if ( (info[0]>1) && (info[1]>1) ) {
7eb157d7 1172 lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2
1173 lp[3]=0.0065; // 0.08*0.08; //SigmaZ2
1174 lp[5]=-6.48e-05;
d695268b 1175 }
1176 else {
7eb157d7 1177 lp[2]=4.80e-06; // 0.00219*0.00219
1178 lp[3]=0.0093; // 0.0964*0.0964;
1179 if (info[0]==1) {
1180 lp[5]=-0.00014;
1181 }
1182 else {
1183 lp[2]=2.79e-06; // 0.0017*0.0017;
1184 lp[3]=0.00935; // 0.967*0.967;
1185 lp[5]=-4.32e-05;
1186 }
d695268b 1187 }
1188
00a7cc50 1189 AliITSRecPoint * cl2;
138df073 1190 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
308b5ea4 1191
00a7cc50 1192 cl2->SetChargeRatio(ratio);
1193 cl2->SetType(7);
308b5ea4 1194 fgPairs[ip*nn+jn] =7;
00a7cc50 1195 if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster
1196 cl2->SetType(8);
308b5ea4 1197 fgPairs[ip*nn+jn]=8;
00a7cc50 1198 }
04366a57 1199 ncl++;
04366a57 1200 }
1201 //
1202 // add second pair
1203 // if (!(cused1[ip]||cused2[jn2])){
a64f9843 1204 if ( (fgPairs[ip*nn+jn2]==100)&&(neg[jn2].GetQ() ) ) { //
1205
d036ccd3 1206 Float_t yn=neg[jn2].GetY();
1207 Double_t yp=pos[ip].GetY();
2069484c 1208
8be4e1b1 1209 Float_t xt, zt;
1210 seg->GetPadCxz(yn, yp, xt, zt);
1211
1212 xbest=xt; zbest=zt;
2069484c 1213
04366a57 1214 qbest =neg[jn2].GetQ();
d036ccd3 1215
b4704be3 1216 {
8be4e1b1 1217 Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
b4704be3 1218 mT2L->MasterToLocal(loc,trk);
1219 lp[0]=trk[1];
1220 lp[1]=trk[2];
1221 }
d695268b 1222
04366a57 1223 lp[4]=qbest; //Q
1224 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1225 for (Int_t ilab=0;ilab<3;ilab++){
1226 milab[ilab] = pos[ip].GetLabel(ilab);
1227 milab[ilab+3] = neg[jn2].GetLabel(ilab);
1228 }
1229 //
1230 CheckLabels2(milab);
1231 ratio = (pos[ip].GetQ()-neg[jn2].GetQ())/(pos[ip].GetQ()+neg[jn2].GetQ());
1232 milab[3]=(((ip<<10) + jn2)<<10) + idet; // pos|neg|det
1233 Int_t info[3] = {pos[ip].GetNd(),neg[jn2].GetNd(),fNlayer[fModule]};
d695268b 1234
7eb157d7 1235 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
1236 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
d695268b 1237 // out-of-diagonal element of covariance matrix
1238 if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1239 else if ( (info[0]>1) && (info[1]>1) ) {
7eb157d7 1240 lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2
1241 lp[3]=0.0065; // 0.08*0.08; //SigmaZ2
1242 lp[5]=-6.48e-05;
d695268b 1243 }
1244 else {
7eb157d7 1245 lp[2]=4.80e-06; // 0.00219*0.00219
1246 lp[3]=0.0093; // 0.0964*0.0964;
1247 if (info[0]==1) {
1248 lp[5]=-0.00014;
1249 }
1250 else {
1251 lp[2]=2.79e-06; // 0.0017*0.0017;
1252 lp[3]=0.00935; // 0.967*0.967;
1253 lp[5]=-4.32e-05;
1254 }
d695268b 1255 }
1256
00a7cc50 1257 AliITSRecPoint * cl2;
138df073 1258 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
308b5ea4 1259
308b5ea4 1260
00a7cc50 1261 cl2->SetChargeRatio(ratio);
308b5ea4 1262 fgPairs[ip*nn+jn2]=7;
00a7cc50 1263 cl2->SetType(7);
1264 if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster
1265 cl2->SetType(8);
308b5ea4 1266 fgPairs[ip*nn+jn2]=8;
00a7cc50 1267 }
04366a57 1268 ncl++;
04366a57 1269 }
1270 cused1[ip]++;
1271 cused2[jn]++;
1272 cused2[jn2]++;
308b5ea4 1273
a64f9843 1274 } // charge matching condition
1275
1276 } // 2 Nside cross 1 Pside
1277 } // loop over Pside clusters
1278
1279
1280
1281 for (Int_t ip=0;ip<np;ip++){
d695268b 1282
1283 if(cused1[ip]) continue;
1284
1285
8be4e1b1 1286 Float_t xbest=1000,zbest=1000,qbest=0;
04366a57 1287 //
a64f9843 1288 // 2x2 clusters
04366a57 1289 //
d695268b 1290 if ( (cnegative[ip]==2) && cpositive[negativepair[10*ip]]==2){
1291 Float_t minchargediff =4.;
1292 Float_t minchargeratio =0.2;
1293
1294 Int_t j=-1;
1295 for (Int_t di=0;di<cnegative[ip];di++){
1296 Int_t jc = negativepair[ip*10+di];
1297 Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1298 ratio = (pos[ip].GetQ()-neg[jc].GetQ())/(pos[ip].GetQ()+neg[jc].GetQ());
1299 //if (TMath::Abs(chargedif)<minchargediff){
1300 if (TMath::Abs(ratio)<0.2){
1301 j =jc;
1302 minchargediff = TMath::Abs(chargedif);
1303 minchargeratio = TMath::Abs(ratio);
1304 }
1305 }
1306 if (j<0) continue; // not proper cluster
1307
1308
1309 Int_t count =0;
1310 for (Int_t di=0;di<cnegative[ip];di++){
1311 Int_t jc = negativepair[ip*10+di];
1312 Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1313 if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1314 }
1315 if (count>1) continue; // more than one "proper" cluster for positive
1316 //
1317
1318 count =0;
1319 for (Int_t dj=0;dj<cpositive[j];dj++){
1320 Int_t ic = positivepair[j*10+dj];
1321 Float_t chargedif = pos[ic].GetQ()-neg[j].GetQ();
1322 if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1323 }
1324 if (count>1) continue; // more than one "proper" cluster for negative
1325
1326 Int_t jp = 0;
1327
1328 count =0;
1329 for (Int_t dj=0;dj<cnegative[jp];dj++){
1330 Int_t ic = positivepair[jp*10+dj];
1331 Float_t chargedif = pos[ic].GetQ()-neg[jp].GetQ();
1332 if (TMath::Abs(chargedif)<minchargediff+4.) count++;
1333 }
1334 if (count>1) continue;
1335 if (fgPairs[ip*nn+j]<100) continue;
1336 //
1337
1338
1339
1340 //almost gold clusters
1341 Float_t yp=pos[ip].GetY();
1342 Float_t yn=neg[j].GetY();
1343 Float_t xt, zt;
1344 seg->GetPadCxz(yn, yp, xt, zt);
1345 xbest=xt; zbest=zt;
1346 qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
1347 {
1348 Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
1349 mT2L->MasterToLocal(loc,trk);
1350 lp[0]=trk[1];
1351 lp[1]=trk[2];
1352 }
1353 lp[4]=qbest; //Q
1354 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1355 for (Int_t ilab=0;ilab<3;ilab++){
1356 milab[ilab] = pos[ip].GetLabel(ilab);
1357 milab[ilab+3] = neg[j].GetLabel(ilab);
1358 }
1359 //
1360 CheckLabels2(milab);
1361 if ((neg[j].GetQ()==0)&&(pos[ip].GetQ()==0)) continue; // reject crosses of bad strips!!
1362 ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
1363 milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
1364 Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
1365
7eb157d7 1366 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
1367 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
d695268b 1368 // out-of-diagonal element of covariance matrix
1369 if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1370 else if ( (info[0]>1) && (info[1]>1) ) {
7eb157d7 1371 lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2
1372 lp[3]=0.0065; // 0.08*0.08; //SigmaZ2
1373 lp[5]=-6.48e-05;
d695268b 1374 }
1375 else {
7eb157d7 1376 lp[2]=4.80e-06; // 0.00219*0.00219
1377 lp[3]=0.0093; // 0.0964*0.0964;
1378 if (info[0]==1) {
1379 lp[5]=-0.00014;
1380 }
1381 else {
1382 lp[2]=2.79e-06; // 0.0017*0.0017;
1383 lp[3]=0.00935; // 0.967*0.967;
1384 lp[5]=-4.32e-05;
1385 }
d695268b 1386 }
1387
1388 AliITSRecPoint * cl2;
138df073 1389 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
d695268b 1390
1391 cl2->SetChargeRatio(ratio);
1392 cl2->SetType(10);
1393 fgPairs[ip*nn+j]=10;
1394 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
1395 cl2->SetType(11);
1396 fgPairs[ip*nn+j]=11;
1397 }
1398 cused1[ip]++;
1399 cused2[j]++;
d695268b 1400 ncl++;
1401
1402 } // 2X2
1403 } // loop over Pside 1Dclusters
1404
1405
1406
1407 for (Int_t ip=0;ip<np;ip++){
1408
1409 if(cused1[ip]) continue;
1410
1411
1412 Float_t xbest=1000,zbest=1000,qbest=0;
1413 //
1414 // manyxmany clusters
1415 //
a64f9843 1416 if ( (cnegative[ip]<5) && cpositive[negativepair[10*ip]]<5){
1417 Float_t minchargediff =4.;
1418 Int_t j=-1;
1419 for (Int_t di=0;di<cnegative[ip];di++){
1420 Int_t jc = negativepair[ip*10+di];
1421 Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1422 if (TMath::Abs(chargedif)<minchargediff){
1423 j =jc;
1424 minchargediff = TMath::Abs(chargedif);
1425 }
00a7cc50 1426 }
a64f9843 1427 if (j<0) continue; // not proper cluster
1428
1429 Int_t count =0;
1430 for (Int_t di=0;di<cnegative[ip];di++){
1431 Int_t jc = negativepair[ip*10+di];
1432 Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1433 if (TMath::Abs(chargedif)<minchargediff+3.) count++;
00a7cc50 1434 }
a64f9843 1435 if (count>1) continue; // more than one "proper" cluster for positive
1436 //
00a7cc50 1437
a64f9843 1438 count =0;
1439 for (Int_t dj=0;dj<cpositive[j];dj++){
1440 Int_t ic = positivepair[j*10+dj];
1441 Float_t chargedif = pos[ic].GetQ()-neg[j].GetQ();
1442 if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1443 }
1444 if (count>1) continue; // more than one "proper" cluster for negative
1445
1446 Int_t jp = 0;
1447
1448 count =0;
1449 for (Int_t dj=0;dj<cnegative[jp];dj++){
1450 Int_t ic = positivepair[jp*10+dj];
1451 Float_t chargedif = pos[ic].GetQ()-neg[jp].GetQ();
1452 if (TMath::Abs(chargedif)<minchargediff+4.) count++;
1453 }
1454 if (count>1) continue;
1455 if (fgPairs[ip*nn+j]<100) continue;
1456 //
1457
1458 //almost gold clusters
d036ccd3 1459 Float_t yp=pos[ip].GetY();
1460 Float_t yn=neg[j].GetY();
1461
8be4e1b1 1462
1463 Float_t xt, zt;
1464 seg->GetPadCxz(yn, yp, xt, zt);
1465
1466 xbest=xt; zbest=zt;
1467
a64f9843 1468 qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
d036ccd3 1469
a64f9843 1470 {
8be4e1b1 1471 Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
a64f9843 1472 mT2L->MasterToLocal(loc,trk);
1473 lp[0]=trk[1];
1474 lp[1]=trk[2];
1475 }
a64f9843 1476 lp[4]=qbest; //Q
1477 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1478 for (Int_t ilab=0;ilab<3;ilab++){
1479 milab[ilab] = pos[ip].GetLabel(ilab);
1480 milab[ilab+3] = neg[j].GetLabel(ilab);
1481 }
1482 //
1483 CheckLabels2(milab);
db6e54cd 1484 if ((neg[j].GetQ()==0)&&(pos[ip].GetQ()==0)) continue; // reject crosses of bad strips!!
a64f9843 1485 ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
1486 milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
1487 Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
d695268b 1488
7eb157d7 1489 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
1490 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
d695268b 1491 // out-of-diagonal element of covariance matrix
1492 if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1493 else if ( (info[0]>1) && (info[1]>1) ) {
7eb157d7 1494 lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2
1495 lp[3]=0.0065; // 0.08*0.08; //SigmaZ2
1496 lp[5]=-6.48e-05;
d695268b 1497 }
1498 else {
7eb157d7 1499 lp[2]=4.80e-06; // 0.00219*0.00219
1500 lp[3]=0.0093; // 0.0964*0.0964;
1501 if (info[0]==1) {
1502 lp[5]=-0.00014;
1503 }
1504 else {
1505 lp[2]=2.79e-06; // 0.0017*0.0017;
1506 lp[3]=0.00935; // 0.967*0.967;
1507 lp[5]=-4.32e-05;
1508 }
d695268b 1509 }
1510
a64f9843 1511 AliITSRecPoint * cl2;
138df073 1512 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
d036ccd3 1513
a64f9843 1514 cl2->SetChargeRatio(ratio);
d695268b 1515 cl2->SetType(12);
1516 fgPairs[ip*nn+j]=12;
a64f9843 1517 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
d695268b 1518 cl2->SetType(13);
1519 fgPairs[ip*nn+j]=13;
a64f9843 1520 }
1521 cused1[ip]++;
1522 cused2[j]++;
a64f9843 1523 ncl++;
1524
1525 } // manyXmany
1526 } // loop over Pside 1Dclusters
1527
a64f9843 1528 } // use charge matching
1529
04366a57 1530
a64f9843 1531 // recover all the other crosses
04366a57 1532 //
1533 for (Int_t i=0; i<np; i++) {
8be4e1b1 1534 Float_t xbest=1000,zbest=1000,qbest=0;
d036ccd3 1535 Float_t yp=pos[i].GetY();
a64f9843 1536 if ((pos[i].GetQ()>0)&&(pos[i].GetQ()<3)) continue;
04366a57 1537 for (Int_t j=0; j<nn; j++) {
1538 // for (Int_t di = 0;di<cpositive[i];di++){
1539 // Int_t j = negativepair[10*i+di];
a64f9843 1540 if ((neg[j].GetQ()>0)&&(neg[j].GetQ()<3)) continue;
1541
1542 if ((neg[j].GetQ()==0)&&(pos[i].GetQ()==0)) continue; // reject crosses of bad strips!!
1543
04366a57 1544 if (cused2[j]||cused1[i]) continue;
308b5ea4 1545 if (fgPairs[i*nn+j]>0 &&fgPairs[i*nn+j]<100) continue;
04366a57 1546 ratio = (pos[i].GetQ()-neg[j].GetQ())/(pos[i].GetQ()+neg[j].GetQ());
d036ccd3 1547 Float_t yn=neg[j].GetY();
1548
8be4e1b1 1549 Float_t xt, zt;
1550 seg->GetPadCxz(yn, yp, xt, zt);
d036ccd3 1551
5846520e 1552 if (TMath::Abs(xt)<hwSSD)
1553 if (TMath::Abs(zt)<hlSSD) {
8be4e1b1 1554 xbest=xt; zbest=zt;
1555
04366a57 1556 qbest=0.5*(pos[i].GetQ()+neg[j].GetQ());
d036ccd3 1557
b4704be3 1558 {
8be4e1b1 1559 Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
b4704be3 1560 mT2L->MasterToLocal(loc,trk);
1561 lp[0]=trk[1];
1562 lp[1]=trk[2];
1563 }
04366a57 1564 lp[4]=qbest; //Q
1565 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1566 for (Int_t ilab=0;ilab<3;ilab++){
1567 milab[ilab] = pos[i].GetLabel(ilab);
1568 milab[ilab+3] = neg[j].GetLabel(ilab);
1569 }
1570 //
1571 CheckLabels2(milab);
1572 milab[3]=(((i<<10) + j)<<10) + idet; // pos|neg|det
1573 Int_t info[3] = {pos[i].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
d695268b 1574
7eb157d7 1575 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
1576 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
d695268b 1577 // out-of-diagonal element of covariance matrix
1578 if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1579 else if ( (info[0]>1) && (info[1]>1) ) {
7eb157d7 1580 lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2
1581 lp[3]=0.0065; // 0.08*0.08; //SigmaZ2
1582 lp[5]=-6.48e-05;
d695268b 1583 }
1584 else {
7eb157d7 1585 lp[2]=4.80e-06; // 0.00219*0.00219
1586 lp[3]=0.0093; // 0.0964*0.0964;
1587 if (info[0]==1) {
1588 lp[5]=-0.00014;
1589 }
1590 else {
1591 lp[2]=2.79e-06; // 0.0017*0.0017;
1592 lp[3]=0.00935; // 0.967*0.967;
1593 lp[5]=-4.32e-05;
1594 }
d695268b 1595 }
1596
00a7cc50 1597 AliITSRecPoint * cl2;
138df073 1598 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
308b5ea4 1599
00a7cc50 1600 cl2->SetChargeRatio(ratio);
1601 cl2->SetType(100+cpositive[j]+cnegative[i]);
a64f9843 1602
1603 if(pos[i].GetQ()==0) cl2->SetType(200+cpositive[j]+cnegative[i]);
1604 if(neg[j].GetQ()==0) cl2->SetType(300+cpositive[j]+cnegative[i]);
04366a57 1605 ncl++;
04366a57 1606 }
1607 }
1608 }
d695268b 1609
1610
d695268b 1611
7101948c 1612 if(repa->GetUseBadChannelsInClusterFinderSSD()==kTRUE) {
1613
1614 //---------------------------------------------------------
1615 // recover crosses of good 1D clusters with bad strips on the other side
1616 // Note1: at first iteration skip modules with a bad side (or almost), (would produce too many fake!)
1617 // Note2: for modules with a bad side see below
1618
1619 AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)GetResp(fModule);
1620 Int_t countPbad=0, countNbad=0;
1621 for(Int_t ib=0; ib<768; ib++) {
1622 if(cal->IsPChannelBad(ib)) countPbad++;
1623 if(cal->IsNChannelBad(ib)) countNbad++;
1624 }
1625 // AliInfo(Form("module %d has %d P- and %d N-bad strips",fModule,countPbad,countNbad));
d695268b 1626
7101948c 1627 if( (countPbad<100) && (countNbad<100) ) { // no bad side!!
d695268b 1628
7101948c 1629 for (Int_t i=0; i<np; i++) { // loop over Nside 1Dclusters with no crosses
1630 if(cnegative[i]) continue; // if intersecting Pside clusters continue;
d695268b 1631
7101948c 1632 // for(Int_t ib=0; ib<768; ib++) { // loop over all Pstrips
1633 for(Int_t ib=15; ib<753; ib++) { // loop over all Pstrips
d695268b 1634
7101948c 1635 if(cal->IsPChannelBad(ib)) { // check if strips is bad
1636 Float_t yN=pos[i].GetY();
1637 Float_t xt, zt;
1638 seg->GetPadCxz(1.*ib, yN, xt, zt);
d695268b 1639
7101948c 1640 //----------
1641 // bad Pstrip is crossing the Nside 1Dcluster -> create recpoint
1642 //
5846520e 1643 if ( (TMath::Abs(xt)<hwSSD) && (TMath::Abs(zt)<hlSSD) ) {
7101948c 1644 Double_t loc[3]={xt,0.,zt},trk[3]={0.,0.,0.};
1645 mT2L->MasterToLocal(loc,trk);
1646 lp[0]=trk[1];
1647 lp[1]=trk[2];
1648 lp[4]=pos[i].GetQ(); //Q
1649 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1650 for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = pos[i].GetLabel(ilab);
1651 CheckLabels2(milab);
1652 milab[3]=( (i<<10) << 10 ) + idet; // pos|neg|det
1653 Int_t info[3] = {pos[i].GetNd(),0,fNlayer[fModule]};
1654
7eb157d7 1655 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
1656 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
7101948c 1657 lp[5]=-0.00012; // out-of-diagonal element of covariance matrix
7eb157d7 1658 if (info[0]>1) {
1659 lp[2]=4.80e-06;
1660 lp[3]=0.0093;
1661 lp[5]=0.00014;
1662 }
1663
7101948c 1664 AliITSRecPoint * cl2;
138df073 1665 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
7101948c 1666 cl2->SetChargeRatio(1.);
d695268b 1667 cl2->SetType(50);
7101948c 1668 ncl++;
1669 } // cross is within the detector
1670 //
1671 //--------------
1672
1673 } // bad Pstrip
d695268b 1674
7101948c 1675 } // end loop over Pstrips
d695268b 1676
7101948c 1677 } // end loop over Nside 1D clusters
d695268b 1678
7101948c 1679 for (Int_t j=0; j<nn; j++) { // loop over Pside 1D clusters with no crosses
1680 if(cpositive[j]) continue;
1681
1682 // for(Int_t ib=0; ib<768; ib++) { // loop over all Nside strips
1683 for(Int_t ib=15; ib<753; ib++) { // loop over all Nside strips
1684
1685 if(cal->IsNChannelBad(ib)) { // check if strip is bad
1686 Float_t yP=neg[j].GetY();
1687 Float_t xt, zt;
1688 seg->GetPadCxz(yP, 1.*ib, xt, zt);
1689
1690 //----------
1691 // bad Nstrip is crossing the Pside 1Dcluster -> create recpoint
1692 //
5846520e 1693 if ( (TMath::Abs(xt)<hwSSD) && (TMath::Abs(zt)<hlSSD) ) {
7101948c 1694 Double_t loc[3]={xt,0.,zt},trk[3]={0.,0.,0.};
1695 mT2L->MasterToLocal(loc,trk);
1696 lp[0]=trk[1];
1697 lp[1]=trk[2];
1698 lp[4]=neg[j].GetQ(); //Q
1699 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1700 for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = neg[j].GetLabel(ilab);
1701 CheckLabels2(milab);
1702 milab[3]=( j << 10 ) + idet; // pos|neg|det
1703 Int_t info[3]={0,(Int_t)neg[j].GetNd(),fNlayer[fModule]};
1704
7eb157d7 1705 lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2
1706 lp[3]=0.012; // 0.110*0.110; //SigmaZ2
7101948c 1707 lp[5]=-0.00012; // out-of-diagonal element of covariance matrix
7eb157d7 1708 if (info[0]>1) {
1709 lp[2]=2.79e-06;
1710 lp[3]=0.00935;
1711 lp[5]=-4.32e-05;
1712 }
7101948c 1713
1714 AliITSRecPoint * cl2;
138df073 1715 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
7101948c 1716 cl2->SetChargeRatio(1.);
1717 cl2->SetType(60);
7101948c 1718 ncl++;
1719 } // cross is within the detector
1720 //
1721 //--------------
1722
1723 } // bad Nstrip
1724 } // end loop over Nstrips
1725 } // end loop over Pside 1D clusters
1726
1727 } // no bad sides
1728
1729 //---------------------------------------------------------
d695268b 1730
7101948c 1731 else if( (countPbad>700) && (countNbad<100) ) { // bad Pside!!
d695268b 1732
7101948c 1733 for (Int_t i=0; i<np; i++) { // loop over Nside 1Dclusters with no crosses
1734 if(cnegative[i]) continue; // if intersecting Pside clusters continue;
d695268b 1735
7101948c 1736 Float_t xt, zt;
1737 Float_t yN=pos[i].GetY();
1738 Float_t yP=0.;
1739 if (seg->GetLayer()==5) yP = yN + (7.6/1.9);
1740 else yP = yN - (7.6/1.9);
1741 seg->GetPadCxz(yP, yN, xt, zt);
1742
5846520e 1743 if ( (TMath::Abs(xt)<hwSSD) && (TMath::Abs(zt)<hlSSD) ) {
7101948c 1744 Double_t loc[3]={xt,0.,zt},trk[3]={0.,0.,0.};
1745 mT2L->MasterToLocal(loc,trk);
1746 lp[0]=trk[1];
1747 lp[1]=trk[2];
1748 lp[4]=pos[i].GetQ(); //Q
1749 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1750 for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = pos[i].GetLabel(ilab);
1751 CheckLabels2(milab);
1752 milab[3]=( (i<<10) << 10 ) + idet; // pos|neg|det
1753 Int_t info[3] = {(Int_t)pos[i].GetNd(),0,fNlayer[fModule]};
d695268b 1754
7eb157d7 1755 lp[2]=0.00098; // 0.031*0.031; //SigmaY2
1756 lp[3]=1.329; // 1.15*1.15; //SigmaZ2
1757 lp[5]=-0.0359;
1758 if(info[0]>1) lp[2]=0.00097;
1759
7101948c 1760 AliITSRecPoint * cl2;
138df073 1761 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
7101948c 1762 cl2->SetChargeRatio(1.);
1763 cl2->SetType(70);
7101948c 1764 ncl++;
1765 } // cross is within the detector
1766 //
1767 //--------------
1768
1769 } // end loop over Nside 1D clusters
1770
1771 } // bad Pside module
d695268b 1772
7101948c 1773 else if( (countNbad>700) && (countPbad<100) ) { // bad Nside!!
d695268b 1774
7101948c 1775 for (Int_t j=0; j<nn; j++) { // loop over Pside 1D clusters with no crosses
1776 if(cpositive[j]) continue;
1777
1778 Float_t xt, zt;
1779 Float_t yP=neg[j].GetY();
1780 Float_t yN=0.;
1781 if (seg->GetLayer()==5) yN = yP - (7.6/1.9);
1782 else yN = yP + (7.6/1.9);
1783 seg->GetPadCxz(yP, yN, xt, zt);
1784
5846520e 1785 if ( (TMath::Abs(xt)<hwSSD) && (TMath::Abs(zt)<hlSSD) ) {
7101948c 1786 Double_t loc[3]={xt,0.,zt},trk[3]={0.,0.,0.};
1787 mT2L->MasterToLocal(loc,trk);
1788 lp[0]=trk[1];
1789 lp[1]=trk[2];
1790 lp[4]=neg[j].GetQ(); //Q
1791 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1792 for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = neg[j].GetLabel(ilab);
1793 CheckLabels2(milab);
1794 milab[3]=( j << 10 ) + idet; // pos|neg|det
1795 Int_t info[3] = {0,(Int_t)neg[j].GetNd(),fNlayer[fModule]};
d695268b 1796
7eb157d7 1797 lp[2]=7.27e-05; // 0.0085*0.0085; //SigmaY2
1798 lp[3]=1.33; // 1.15*1.15; //SigmaZ2
1799 lp[5]=0.00931;
1800 if(info[1]>1) lp[2]=6.91e-05;
7101948c 1801
1802 AliITSRecPoint * cl2;
138df073 1803 cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info);
7101948c 1804 cl2->SetChargeRatio(1.);
1805 cl2->SetType(80);
7101948c 1806 ncl++;
1807 } // cross is within the detector
1808 //
1809 //--------------
1810
1811 } // end loop over Pside 1D clusters
d695268b 1812
7101948c 1813 } // bad Nside module
d695268b 1814
7101948c 1815 //---------------------------------------------------------
d695268b 1816
7101948c 1817 } // use bad channels
d695268b 1818
d695268b 1819 //cout<<ncl<<" clusters for this module"<<endl;
1820
8be4e1b1 1821 delete [] cnegative;
1822 delete [] cused1;
1823 delete [] negativepair;
1824 delete [] cpositive;
1825 delete [] cused2;
1826 delete [] positivepair;
04366a57 1827
0a56760a 1828}