X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=ITS%2FAliITSClusterFinderV2SSD.cxx;h=0d622498c82ee27fd2bd4cb998f80b563e2f8f37;hb=387f1a9ed96bc434a1315e06cfe3ab7f2c0bd040;hp=5dd545210fb54e73f50e30176e8b4d83a33b5d15;hpb=6490e1c59c3b02c8918d9c6442721d9639dbff55;p=u%2Fmrichter%2FAliRoot.git diff --git a/ITS/AliITSClusterFinderV2SSD.cxx b/ITS/AliITSClusterFinderV2SSD.cxx index 5dd545210fb..0d622498c82 100644 --- a/ITS/AliITSClusterFinderV2SSD.cxx +++ b/ITS/AliITSClusterFinderV2SSD.cxx @@ -1,5 +1,5 @@ /************************************************************************** - * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. * + * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. * * * * Author: The ALICE Off-line Project. * * Contributors are mentioned in the code where appropriate. * @@ -12,39 +12,104 @@ * about the suitability of this software for any purpose. It is * * provided "as is" without express or implied warranty. * **************************************************************************/ + +/* $Id$ */ + //////////////////////////////////////////////////////////////////////////// // Implementation of the ITS clusterer V2 class // // // // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch // +// Last revision: 13-05-09 Enrico Fragiacomo // +// enrico.fragiacomo@ts.infn.it // // // /////////////////////////////////////////////////////////////////////////// - #include "AliITSClusterFinderV2SSD.h" + +#include +#include + +#include "AliLog.h" +#include "AliMagF.h" #include "AliITSRecPoint.h" +#include "AliITSRecPointContainer.h" +#include "AliITSgeomTGeo.h" #include "AliITSDetTypeRec.h" #include "AliRawReader.h" #include "AliITSRawStreamSSD.h" #include +#include #include "AliITSdigitSSD.h" +#include "AliITSReconstructor.h" +#include "AliITSCalibrationSSD.h" +#include "AliITSsegmentationSSD.h" -ClassImp(AliITSClusterFinderV2SSD) +Short_t *AliITSClusterFinderV2SSD::fgPairs = 0x0; +Int_t AliITSClusterFinderV2SSD::fgPairsSize = 0; +const Float_t AliITSClusterFinderV2SSD::fgkThreshold = 5.; +const Float_t AliITSClusterFinderV2SSD::fgkCosmic2008StripShifts[16][9] = + {{-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35}, // DDL 512 + {-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35}, // DDL 513 + {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 514 + {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 515 + { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00}, // DDL 516 + { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00}, // DDL 517 + {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 518 + {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 519 + {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.25,-0.15}, // DDL 520 + {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15}, // DDL 521 + {-0.10,-0.10,-0.10,-0.40,-0.40,-0.40,-0.10,-0.10,-0.45}, // DDL 522 + {-0.10,-0.10,-0.10,-0.35,-0.35,-0.35,-0.10,-0.35,-0.50}, // DDL 523 + { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00}, // DDL 524 + { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00}, // DDL 525 + { 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35}, // DDL 526 + { 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45}}; // DDL 527 -AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(AliITSDetTypeRec* dettyp):AliITSClusterFinderV2(dettyp), -fLastSSD1(0), -fYpitchSSD(0.0095), -fHwSSD(3.65), -fHlSSD(2.00), -fTanP(0.0275), -fTanN(0.0075){ +ClassImp(AliITSClusterFinderV2SSD) - //Default constructor - fLastSSD1=fDetTypeRec->GetITSgeom()->GetModuleIndex(6,1,1)-1; + AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(AliITSDetTypeRec* dettyp):AliITSClusterFinder(dettyp),fLastSSD1(AliITSgeomTGeo::GetModuleIndex(6,1,1)-1), fLorentzShiftP(0), fLorentzShiftN(0) +{ +//Default constructor + static AliITSRecoParam *repa = NULL; + if(!repa){ + repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam(); + if(!repa){ + repa = AliITSRecoParam::GetHighFluxParam(); + AliWarning("Using default AliITSRecoParam class"); + } + } + if (repa->GetCorrectLorentzAngleSSD()) { + AliMagF* field = dynamic_cast(TGeoGlobalMagField::Instance()->GetField()); + if (field == 0) + AliError("Cannot get magnetic field from TGeoGlobalMagField"); + Float_t Bfield = field->SolenoidField(); + // NB: spatial shift has opposite sign for lay 5 and 6, but strip numbering also changes direction, so no sign-change + // Shift due to ExB on drift N-side, units: strip width + fLorentzShiftP = -repa->GetTanLorentzAngleElectronsSSD() * 150.e-4/95.e-4 * Bfield / 5.0; + // Shift due to ExB on drift P-side, units: strip width + fLorentzShiftN = -repa->GetTanLorentzAngleHolesSSD() * 150.e-4/95.e-4 * Bfield / 5.0; + AliDebug(1,Form("Bfield %f Lorentz Shift P-side %f N-side %f",Bfield,fLorentzShiftN,fLorentzShiftP)); + } } +//______________________________________________________________________ +AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(const AliITSClusterFinderV2SSD &cf) : AliITSClusterFinder(cf), fLastSSD1(cf.fLastSSD1), fLorentzShiftP(cf.fLorentzShiftP), fLorentzShiftN(cf.fLorentzShiftN) +{ + // Copy constructor +} + +//______________________________________________________________________ +AliITSClusterFinderV2SSD& AliITSClusterFinderV2SSD::operator=(const AliITSClusterFinderV2SSD& cf ){ + // Assignment operator + + this->~AliITSClusterFinderV2SSD(); + new(this) AliITSClusterFinderV2SSD(cf); + return *this; +} + void AliITSClusterFinderV2SSD::FindRawClusters(Int_t mod){ @@ -60,32 +125,105 @@ void AliITSClusterFinderV2SSD::FindClustersSSD(TClonesArray *alldigits) { //------------------------------------------------------------ Int_t smaxall=alldigits->GetEntriesFast(); if (smaxall==0) return; - TObjArray *digits = new TObjArray; + + + //--------------------------------------- + // load recoparam and calibration + // + static AliITSRecoParam *repa = NULL; + if(!repa){ + repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam(); + if(!repa){ + repa = AliITSRecoParam::GetHighFluxParam(); + AliWarning("Using default AliITSRecoParam class"); + } + } + + AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)GetResp(fModule); + Float_t gain=0; + Float_t noise=0; + //--------------------------------------- + + + //------------------------------------ + // fill the digits array with zero-suppression condition + // Signal is converted in KeV + // + TObjArray digits; for (Int_t i=0;iUncheckedAt(i); - Float_t q=d->GetSignal()/4.29;// temp. fix (for PID purposed - normalis. to be checked) + + if(d->IsSideP()) noise = cal->GetNoiseP(d->GetStripNumber()); + else noise = cal->GetNoiseN(d->GetStripNumber()); + if (d->GetSignal()<3.*noise) continue; + + if(d->IsSideP()) gain = cal->GetGainP(d->GetStripNumber()); + else gain = cal->GetGainN(d->GetStripNumber()); + + Float_t q=gain*d->GetSignal(); // + q=cal->ADCToKeV(q); // converts the charge in KeV from ADC units d->SetSignal(Int_t(q)); - if (d->GetSignal()<3) continue; - digits->AddLast(d); + + digits.AddLast(d); } - Int_t smax = digits->GetEntriesFast(); + Int_t smax = digits.GetEntriesFast(); if (smax==0) return; + //------------------------------------ + const Int_t kMax=1000; Int_t np=0, nn=0; Ali1Dcluster pos[kMax], neg[kMax]; Float_t y=0., q=0., qmax=0.; - Int_t lab[4]={-2,-2,-2,-2}; + Int_t lab[4]={-2,-2,-2,-2}; + Bool_t flag5 = 0; - AliITSdigitSSD *d=(AliITSdigitSSD*)digits->UncheckedAt(0); + /* + cout<<"-----------------------------"<fLastSSD1) + layer = 5; + + //-------------------------------------------------------- + // start 1D-clustering from the first digit in the digits array + // + AliITSdigitSSD *d=(AliITSdigitSSD*)digits.UncheckedAt(0); q += d->GetSignal(); y += d->GetCoord2()*d->GetSignal(); qmax=d->GetSignal(); lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2); + + if(d->IsSideP()) { + noise = cal->GetNoiseP(d->GetStripNumber()); + gain = cal->GetGainP(d->GetStripNumber()); + } + else { + noise = cal->GetNoiseN(d->GetStripNumber()); + gain = cal->GetGainN(d->GetStripNumber()); + } + noise*=gain; + noise=cal->ADCToKeV(noise); // converts noise in KeV from ADC units + + if(qmax>fgkThreshold*noise) flag5=1; // seed for the cluster + + /* + cout<GetSignal()<<" "<GetCoord1()<<" "<GetCoord2()<GetCoord2(); Int_t flag=d->GetCoord1(); + + // Note: the first side which will be processed is supposed to be the + // P-side which is neg Int_t *n=&nn; Ali1Dcluster *c=neg; + if(flag) {n=&np; c=pos;} // in case we have only Nstrips (P was bad!) + Int_t nd=1; Int_t milab[10]; for (Int_t ilab=0;ilab<10;ilab++){ @@ -93,73 +231,101 @@ void AliITSClusterFinderV2SSD::FindClustersSSD(TClonesArray *alldigits) { } milab[0]=d->GetTrack(0); milab[1]=d->GetTrack(1); milab[2]=d->GetTrack(2); + + //---------------------------------------------------------- + // search for neighboring digits + // for (Int_t s=1; sUncheckedAt(s); + d=(AliITSdigitSSD*)digits.UncheckedAt(s); Int_t strip=d->GetCoord2(); - if ((strip-curr) > 1 || flag!=d->GetCoord1()) { - c[*n].SetY(y/q); + + // if digits is not a neighbour or side did not change + // and at least one of the previous digits met the seed condition + // then creates a new 1D cluster + if ( ( ((strip-curr) > 1) || (flag!=d->GetCoord1()) ) ) { + + if(flag5) { + //cout<<"here1"<10&&nd<16){ - c[*n].SetY(y/q-0.3*nd); - c[*n].SetQ(0.5*q); - (*n)++; - if (*n==MAX) { - Error("FindClustersSSD","Too many 1D clusters !"); - return; - } - c[*n].SetY(y/q-0.0*nd); - c[*n].SetQ(0.5*q); - c[*n].SetNd(nd); - (*n)++; - if (*n==MAX) { - Error("FindClustersSSD","Too many 1D clusters !"); - return; - } - // - c[*n].SetY(y/q+0.3*nd); - c[*n].SetQ(0.5*q); - c[*n].SetNd(nd); - c[*n].SetLabels(milab); - } - else{ - */ - if (nd>4&&nd<25) { - c[*n].SetY(y/q-0.25*nd); - c[*n].SetQ(0.5*q); - (*n)++; - if (*n==kMax) { - Error("FindClustersSSD","Too many 1D clusters !"); - return; - } - c[*n].SetY(y/q+0.25*nd); - c[*n].SetQ(0.5*q); - c[*n].SetNd(nd); - c[*n].SetLabels(milab); - } + + if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) { + // Note: fUseUnfoldingInClusterFinderSSD=kFALSE by default in RecoParam + + //Split suspiciously big cluster + if (nd>4&&nd<25) { + c[*n].SetY(y/q-0.25*nd+dLorentz); + c[*n].SetQ(0.5*q); + (*n)++; + if (*n==kMax) { + Error("FindClustersSSD","Too many 1D clusters !"); + return; + } + c[*n].SetY(y/q+0.25*nd+dLorentz); + c[*n].SetQ(0.5*q); + c[*n].SetNd(nd); + c[*n].SetLabels(milab); + } + + } // unfolding is on + (*n)++; if (*n==kMax) { Error("FindClustersSSD","Too many 1D clusters !"); return; } + + } // flag5 set + + // reset everything y=q=qmax=0.; nd=0; + flag5=0; lab[0]=lab[1]=lab[2]=-2; - // - for (Int_t ilab=0;ilab<10;ilab++){ - milab[ilab]=-2; - } - // - if (flag!=d->GetCoord1()) { n=&np; c=pos; } - } + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + + // if side changed from P to N, switch to pos 1D clusters + // (if for some reason the side changed from N to P then do the opposite) + if (flag!=d->GetCoord1()) + { if(!flag) {n=&np; c=pos;} else {n=&nn; c=neg;} } + + } // end create new 1D cluster from previous neighboring digits + + // continues adding digits to the previous cluster + // or start a new one flag=d->GetCoord1(); q += d->GetSignal(); y += d->GetCoord2()*d->GetSignal(); nd++; + + if(d->IsSideP()) { + noise = cal->GetNoiseP(d->GetStripNumber()); + gain = cal->GetGainP(d->GetStripNumber()); + } + else { + noise = cal->GetNoiseN(d->GetStripNumber()); + gain = cal->GetGainN(d->GetStripNumber()); + } + noise*=gain; + noise=cal->ADCToKeV(noise); // converts the charge in KeV from ADC units + + if(d->GetSignal()>fgkThreshold*noise) flag5=1; + + /* + cout<GetSignal()<<" "<GetCoord1()<<" "<GetCoord2()<GetSignal()>qmax) { qmax=d->GetSignal(); lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2); @@ -168,135 +334,373 @@ void AliITSClusterFinderV2SSD::FindClustersSSD(TClonesArray *alldigits) { if (d->GetTrack(ilab)>=0) AddLabel(milab, (d->GetTrack(ilab))); } curr=strip; - } - c[*n].SetY(y/q); - c[*n].SetQ(q); - c[*n].SetNd(nd); - c[*n].SetLabels(lab); - //Split suspiciously big cluster - if (nd>4 && nd<25) { - c[*n].SetY(y/q-0.25*nd); - c[*n].SetQ(0.5*q); - (*n)++; - if (*n==kMax) { - Error("FindClustersSSD","Too many 1D clusters !"); - return; - } - c[*n].SetY(y/q+0.25*nd); - c[*n].SetQ(0.5*q); - c[*n].SetNd(nd); - c[*n].SetLabels(lab); - } - (*n)++; - if (*n==kMax) { - Error("FindClustersSSD","Too many 1D clusters !"); - return; - } - FindClustersSSD(neg, nn, pos, np); + + } // loop over digits, no more digits in the digits array + + + // add the last 1D cluster + if(flag5) { + + // cout<<"here2"<GetUseUnfoldingInClusterFinderSSD()==kTRUE) { + + //Split suspiciously big cluster + if (nd>4 && nd<25) { + c[*n].SetY(y/q-0.25*nd + dLorentz); + c[*n].SetQ(0.5*q); + (*n)++; + if (*n==kMax) { + Error("FindClustersSSD","Too many 1D clusters !"); + return; + } + c[*n].SetY(y/q+0.25*nd + dLorentz); + c[*n].SetQ(0.5*q); + c[*n].SetNd(nd); + c[*n].SetLabels(lab); + } + } // unfolding is on + + (*n)++; + if (*n==kMax) { + Error("FindClustersSSD","Too many 1D clusters !"); + return; + } + + } // if flag5 last 1D cluster added + + + //------------------------------------------------------ + // call FindClustersSSD to pair neg and pos 1D clusters + // and create recpoints from the crosses + // Note1: neg are Pside and pos are Nside!! + // Note2: if there are no Pside digits nn=0 (bad strips!!) (same for Nside) + // + // cout< 0) { + TClonesArray* clusters = rpc->UncheckedGetClusters(fModule); + clusters->Clear(); + FindClustersSSD(neg, nn, pos, np, clusters); + TIter itr(clusters); + AliITSRecPoint *irp; + while ((irp = (AliITSRecPoint*)itr.Next())) fDetTypeRec->AddRecPoint(*irp); + } + //----------------------------------------------------- } -void AliITSClusterFinderV2SSD::RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters){ +void AliITSClusterFinderV2SSD::RawdataToClusters(AliRawReader* rawReader){ - //------------------------------------------------------------ + //------------------------------------------------------------ // This function creates ITS clusters from raw data //------------------------------------------------------------ rawReader->Reset(); AliITSRawStreamSSD inputSSD(rawReader); - FindClustersSSD(&inputSSD,clusters); + FindClustersSSD(&inputSSD); } -void AliITSClusterFinderV2SSD::FindClustersSSD(AliITSRawStream* input, - TClonesArray** clusters) + +void AliITSClusterFinderV2SSD::FindClustersSSD(AliITSRawStreamSSD* input) { //------------------------------------------------------------ // Actual SSD cluster finder for raw data //------------------------------------------------------------ + + AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance(); + static AliITSRecoParam *repa = NULL; + if(!repa){ + repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam(); + if(!repa){ + repa = AliITSRecoParam::GetHighFluxParam(); + AliWarning("Using default AliITSRecoParam class"); + } + } Int_t nClustersSSD = 0; - const Int_t kMax = 1000; - Ali1Dcluster clusters1D[2][kMax]; - Int_t nClusters[2] = {0, 0}; - Int_t lab[3]={-2,-2,-2}; - Float_t q = 0.; - Float_t y = 0.; - Int_t nDigits = 0; - Int_t prevStrip = -1; - Int_t prevFlag = -1; - Int_t prevModule = -1; - - // read raw data input stream + const Int_t kNADC = 12; + const Int_t kMaxADCClusters = 1000; + + Int_t strips[kNADC][2][kMaxADCClusters][2]; // [ADC],[side],[istrip], [0]=istrip [1]=signal + Int_t nStrips[kNADC][2]; + + for( int i=0; iNext(); - - if(input->GetSignal()<3 && next) continue; - // check if a new cluster starts - Int_t strip = input->GetCoord2(); - Int_t flag = input->GetCoord1(); - if ((!next || (input->GetModuleID() != prevModule)|| - (strip-prevStrip > 1) || (flag != prevFlag)) && - (nDigits > 0)) { - if (nClusters[prevFlag] == kMax) { - Error("FindClustersSSD", "Too many 1D clusters !"); - return; + + bool next = input->Next(); + + //* + //* Continue if corrupted input + //* + + if( (!next)&&(input->flag) ){ + AliWarning(Form("ITSClustersFinderSSD: Corrupted data: warning from RawReader")); + continue; + } + + Int_t newDDL = input->GetDDL(); + Int_t newAD = input->GetAD(); + + if( next ){ + if( newDDL<0 || newDDL>15 ){ + AliWarning(Form("ITSClustersFinderSSD: Corrupted data: wrong DDL number (%d)",newDDL)); + continue; + } + + if( newAD<1 || newAD>9 ){ + AliWarning(Form("ITSClustersFinderSSD: Corrupted data: wrong AD number (%d)",newAD)); + continue; } - Ali1Dcluster& cluster = clusters1D[prevFlag][nClusters[prevFlag]++]; - cluster.SetY(y/q); - cluster.SetQ(q); - cluster.SetNd(nDigits); - cluster.SetLabels(lab); + } - //Split suspiciously big cluster - if (nDigits > 4&&nDigits < 25) { - cluster.SetY(y/q - 0.25*nDigits); - cluster.SetQ(0.5*q); - if (nClusters[prevFlag] == kMax) { - Error("FindClustersSSD", "Too many 1D clusters !"); - return; + bool newModule = ( !next || ddl!= newDDL || ad!=newAD ); + + if( newModule && ddl>=0 && ad>=0 ){ + + //* + //* Reconstruct the previous block of 12 modules --- actual clusterfinder + //* + //cout<fLastSSD1) + layer = 5; + + AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)fDetTypeRec->GetCalibrationModel(fModule); + if( !cal ){ + AliWarning(Form("ITSClustersFinderSSD: No calibration found for module (ddl %d ad %d adc %d)",ddl,ad,adc)); + continue; + } + + Float_t dStrip = 0; + + if( repa->GetUseCosmicRunShiftsSSD()) { // Special condition for 2007/2008 cosmic data + dStrip = fgkCosmic2008StripShifts[ddl][ad-1]; + if (TMath::Abs(dStrip) > 1.5){ + AliWarning(Form("Indexing error in Cosmic calibration: ddl = %d, dStrip %f\n",ddl,dStrip)); + dStrip = 0; + } + } + + for( int side=0; side<=1; side++ ){ + + Int_t lab[3]={-2,-2,-2}; + Float_t q = 0.; + Float_t y = 0.; + Int_t nDigits = 0; + Int_t ostrip = -2; + Bool_t snFlag = 0; + + Float_t dLorentz = 0; + if (side==0) { // P-side is neg clust + dLorentz = fLorentzShiftN; + } + else { // N-side is pos clust + dLorentz = fLorentzShiftP; + } + + Int_t n = nStrips[adc][side]; + for( int istr = 0; istrGetNoiseN(strip) :cal->GetNoiseP(strip); + gain = side ?cal->GetGainN(strip) :cal->GetGainP(strip); + stripOK = ( noise>=1. && signal>=3.0*noise + //&& !cal->IsPChannelBad(strip) + ); + } + } else stripOK = 0; // end of data + + bool newCluster = ( (abs(strip-ostrip)!=1) || !stripOK ); + + if( newCluster ){ + + //* Store the previous cluster + + if( nDigits>0 && q>0 && snFlag ){ + + if (nClusters1D[side] >= kMaxADCClusters-1 ) { + AliWarning("HLT ClustersFinderSSD: Too many 1D clusters !"); + }else { + + Ali1Dcluster &cluster = clusters1D[side][nClusters1D[side]++]; + cluster.SetY( y / q + dStrip + dLorentz); + cluster.SetQ(q); + cluster.SetNd(nDigits); + cluster.SetLabels(lab); + //cout<<"cluster 1D side "<GetModuleID() != prevModule)) { - Int_t iModule = prevModule; - - // when all data from a module was read, search for clusters - if (prevFlag >= 0) { - clusters[iModule] = new TClonesArray("AliITSRecPoint"); - fModule = iModule; - FindClustersSSD(&clusters1D[0][0], nClusters[0], - &clusters1D[1][0], nClusters[1], clusters[iModule]); - Int_t nClusters = clusters[iModule]->GetEntriesFast(); - nClustersSSD += nClusters; - } + if( adc>7 ) adc-= 2; // shift ADC numbers 8-13 to 6-11 + + Bool_t side = input->GetSideFlag(); + Int_t strip = input->GetStrip(); + Int_t signal = input->GetSignal(); + + + //cout<<"SSD: "<767 ){ + AliWarning(Form("HLT ClustersFinderSSD: Corrupted data: wrong strip number (ddl %d ad %d adc %d side %d, strip %d", + ddl, ad, adc, side,strip) ); + continue; } + if (strip < 0) continue; + + int &n = nStrips[adc][side]; + if( n >0 ){ + Int_t oldStrip = strips[adc][side][n-1][0]; - // add digit to current cluster - q += input->GetSignal(); - y += strip * input->GetSignal(); - nDigits++; - prevStrip = strip; - prevFlag = flag; - prevModule = input->GetModuleID(); + if( strip==oldStrip ){ + AliWarning(Form("HLT ClustersFinderSSD: Corrupted data: duplicated signal: ddl %d ad %d adc %d, side %d, strip %d", + ddl, ad, adc, side, strip )); + continue; + } + } + strips[adc][side][n][0] = strip; + strips[adc][side][n][1] = signal; + n++; - } + //cout<<"SSD: "<GetDDL()<<" "<GetAD()<<" " + //<GetADC()<<" "<GetSideFlag()<<" "<<((int)input->GetStrip())<<" "<GetSignal()<fLastSSD1) {tann=fTanP; tanp=fTanN;} + + const TGeoHMatrix *mT2L=AliITSgeomTGeo::GetTracking2LocalMatrix(fModule); + + //--------------------------------------- + // load recoparam + // + static AliITSRecoParam *repa = NULL; + if(!repa){ + repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam(); + if(!repa){ + repa = AliITSRecoParam::GetHighFluxParam(); + AliWarning("Using default AliITSRecoParam class"); + } + } + +// TClonesArray &cl=*clusters; + + AliITSsegmentationSSD *seg = dynamic_cast(fDetTypeRec->GetSegmentationModel(2)); + if (fModule>fLastSSD1) + seg->SetLayer(6); + else + seg->SetLayer(5); + + Float_t hwSSD = seg->Dx()*1e-4/2; + Float_t hlSSD = seg->Dz()*1e-4/2; + Int_t idet=fNdet[fModule]; Int_t ncl=0; + // - Int_t negativepair[30000]; - Int_t cnegative[3000]; - Int_t cused1[3000]; - Int_t positivepair[30000]; - Int_t cpositive[3000]; - Int_t cused2[3000]; - for (Int_t i=0;i<3000;i++) {cnegative[i]=0; cused1[i]=0;} - for (Int_t i=0;i<3000;i++) {cpositive[i]=0; cused2[i]=0;} - for (Int_t i=0;i<30000;i++) {negativepair[i]=0; positivepair[i]=0;} - static Short_t pairs[1000][1000]; - memset(pairs,0,sizeof(Short_t)*1000000); -// Short_t ** pairs = new Short_t*[1000]; -// for (Int_t i=0; i<1000; i++) { -// pairs[i] = new Short_t[1000]; -// memset(pairs[i],0,sizeof(Short_t)*1000); -// } + Int_t *cnegative = new Int_t[np]; + Int_t *cused1 = new Int_t[np]; + Int_t *negativepair = new Int_t[10*np]; + Int_t *cpositive = new Int_t[nn]; + Int_t *cused2 = new Int_t[nn]; + Int_t *positivepair = new Int_t[10*nn]; + for (Int_t i=0;i fgPairsSize) { + + if (fgPairs) delete [] fgPairs; + fgPairsSize = 4*np*nn; + fgPairs = new Short_t[fgPairsSize]; + } + memset(fgPairs,0,sizeof(Short_t)*np*nn); + // // find available pairs // + Int_t ncross = 0; for (Int_t i=0; i0) && (pos[i].GetQ()<3) ) continue; for (Int_t j=0; j0) && (neg[j].GetQ()<3) ) continue; + Float_t yn=neg[j].GetY(); + + Float_t xt, zt; + seg->GetPadCxz(yn, yp, xt, zt); + //cout<0) && (pos[i].GetQ()<3) ) continue; for (Int_t j=0; j0) && (neg[j].GetQ()<3) ) continue; + // if both 1Dclusters have an other cross continue if (cpositive[j]&&cnegative[i]) continue; - Float_t yn=neg[j].GetY()*fYpitchSSD; - Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp); - Float_t yt=yn + tann*zt; - zt-=fHlSSD; yt-=fHwSSD; - if (TMath::Abs(yt)GetPadCxz(yn, yp, xt, zt); + + if (TMath::Abs(xt)GetUseChargeMatchingInClusterFinderSSD()==kTRUE) { + + // - // select gold clusters - if ( (cnegative[ip]==1) && cpositive[negativepair[10*ip]]==1){ - Float_t yp=pos[ip].GetY()*fYpitchSSD; - Int_t j = negativepair[10*ip]; - ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ()); - // - Float_t yn=neg[j].GetY()*fYpitchSSD; - Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp); - Float_t yt=yn + tann*zt; - zt-=fHlSSD; yt-=fHwSSD; - ybest=yt; zbest=zt; - qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ()); - lp[0]=-(-ybest+fYshift[fModule]); - lp[1]= -zbest+fZshift[fModule]; - lp[2]=0.0025*0.0025; //SigmaY2 - lp[3]=0.110*0.110; //SigmaZ2 - - lp[4]=qbest; //Q - for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; - for (Int_t ilab=0;ilab<3;ilab++){ - milab[ilab] = pos[ip].GetLabel(ilab); - milab[ilab+3] = neg[j].GetLabel(ilab); - } + // sign gold tracks + // + for (Int_t ip=0;ipGetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(1); - pairs[ip][j]=1; - if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster - cl2->SetType(2); - pairs[ip][j]=2; + // select gold clusters + if ( (cnegative[ip]==1) && cpositive[negativepair[10*ip]]==1){ + Float_t yp=pos[ip].GetY(); + Int_t j = negativepair[10*ip]; + + if( (pos[ip].GetQ()==0) && (neg[j].GetQ() ==0) ) { + // both bad, hence continue; + // mark both as used (to avoid recover at the end) + cused1[ip]++; + cused2[j]++; + continue; } - cused1[ip]++; - cused2[j]++; + + ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ()); + //cout<<"ratio="< ratio=1 and the following condition is met + if (TMath::Abs(ratio)>0.2) continue; // note: 0.2=3xsigma_ratio calculated in cosmics tests + + // + Float_t yn=neg[j].GetY(); - } - else{ - cl2 = new AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(1); - pairs[ip][j]=1; - if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster - cl2->SetType(2); - pairs[ip][j]=2; + Float_t xt, zt; + seg->GetPadCxz(yn, yp, xt, zt); + + xbest=xt; zbest=zt; + + + qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ()); + 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 + + { + Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.}; + mT2L->MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; } - cused1[ip]++; - cused2[j]++; - fDetTypeRec->AddRecPoint(*cl2); + lp[4]=qbest; //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++){ + milab[ilab] = pos[ip].GetLabel(ilab); + milab[ilab+3] = neg[j].GetLabel(ilab); + } + // + CheckLabels2(milab); + milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det + Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]}; + + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + // out-of-diagonal element of covariance matrix + if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012; + else if ( (info[0]>1) && (info[1]>1) ) { + lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2 + lp[3]=0.0065; // 0.08*0.08; //SigmaZ2 + lp[5]=-6.48e-05; + } + else { + lp[2]=4.80e-06; // 0.00219*0.00219 + lp[3]=0.0093; // 0.0964*0.0964; + if (info[0]==1) { + lp[5]=-0.00014; + } + else { + lp[2]=2.79e-06; // 0.0017*0.0017; + lp[3]=0.00935; // 0.967*0.967; + lp[5]=-4.32e-05; + } + } + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + + cl2->SetChargeRatio(ratio); + cl2->SetType(1); + fgPairs[ip*nn+j]=1; + + if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster + cl2->SetType(2); + fgPairs[ip*nn+j]=2; + } + + if(pos[ip].GetQ()==0) cl2->SetType(3); + if(neg[j].GetQ()==0) cl2->SetType(4); + + cused1[ip]++; + cused2[j]++; + + ncl++; } - ncl++; } - } - for (Int_t ip=0;ipGetPadCxz(yn, yp, xt, zt); + + xbest=xt; zbest=zt; - AliITSRecPoint * cl2; - if(clusters){ - cl2 = new (cl[ncl]) AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(5); - pairs[ip][in] = 5; - if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster - cl2->SetType(6); - pairs[ip][in] = 6; - } - } - else{ - cl2 = new AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(5); - pairs[ip][in] = 5; - if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster - cl2->SetType(6); - pairs[ip][in] = 6; + qbest =pos[ip].GetQ(); + Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.}; + mT2L->MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + + lp[4]=qbest; //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++){ + milab[ilab] = pos[ip].GetLabel(ilab); + milab[ilab+3] = neg[in].GetLabel(ilab); } + // + CheckLabels2(milab); + ratio = (pos[ip].GetQ()-neg[in].GetQ())/(pos[ip].GetQ()+neg[in].GetQ()); + milab[3]=(((ip<<10) + in)<<10) + idet; // pos|neg|det + Int_t info[3] = {pos[ip].GetNd(),neg[in].GetNd(),fNlayer[fModule]}; - fDetTypeRec->AddRecPoint(*cl2); + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + // out-of-diagonal element of covariance matrix + if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012; + else if ( (info[0]>1) && (info[1]>1) ) { + lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2 + lp[3]=0.0065; // 0.08*0.08; //SigmaZ2 + lp[5]=-6.48e-05; + } + else { + lp[2]=4.80e-06; // 0.00219*0.00219 + lp[3]=0.0093; // 0.0964*0.0964; + if (info[0]==1) { + lp[5]=-0.00014; + } + else { + lp[2]=2.79e-06; // 0.0017*0.0017; + lp[3]=0.00935; // 0.967*0.967; + lp[5]=-4.32e-05; } - ncl++; } - - // - // add second pair - - // if (!(cused1[ip2] || cused2[in])){ // - if (pairs[ip2][in]==100){ - Float_t yp=pos[ip2].GetY()*fYpitchSSD; - Float_t yn=neg[in].GetY()*fYpitchSSD; - Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp); - Float_t yt=yn + tann*zt; - zt-=fHlSSD; yt-=fHwSSD; - ybest =yt; zbest=zt; - qbest =pos[ip2].GetQ(); - lp[0]=-(-ybest+fYshift[fModule]); - lp[1]= -zbest+fZshift[fModule]; - lp[2]=0.0025*0.0025; //SigmaY2 - lp[3]=0.110*0.110; //SigmaZ2 - - lp[4]=qbest; //Q - for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; - for (Int_t ilab=0;ilab<3;ilab++){ - milab[ilab] = pos[ip2].GetLabel(ilab); - milab[ilab+3] = neg[in].GetLabel(ilab); + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + cl2->SetChargeRatio(ratio); + cl2->SetType(5); + fgPairs[ip*nn+in] = 5; + if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster + cl2->SetType(6); + fgPairs[ip*nn+in] = 6; + } + ncl++; } + + // - CheckLabels2(milab); - ratio = (pos[ip2].GetQ()-neg[in].GetQ())/(pos[ip2].GetQ()+neg[in].GetQ()); - milab[3]=(((ip2<<10) + in)<<10) + idet; // pos|neg|det - Int_t info[3] = {pos[ip2].GetNd(),neg[in].GetNd(),fNlayer[fModule]}; + // add second pair - AliITSRecPoint * cl2; - if(clusters){ - cl2 = new (cl[ncl]) AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(5); - pairs[ip2][in] =5; - if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster - cl2->SetType(6); - pairs[ip2][in] =6; + // if (!(cused1[ip2] || cused2[in])){ // + if ( (fgPairs[ip2*nn+in]==100) && (pos[ip2].GetQ()) ) { + + Float_t yp=pos[ip2].GetY(); + Float_t yn=neg[in].GetY(); + + Float_t xt, zt; + seg->GetPadCxz(yn, yp, xt, zt); + + xbest=xt; zbest=zt; + + qbest =pos[ip2].GetQ(); + + Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.}; + mT2L->MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + + lp[4]=qbest; //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++){ + milab[ilab] = pos[ip2].GetLabel(ilab); + milab[ilab+3] = neg[in].GetLabel(ilab); } + // + CheckLabels2(milab); + ratio = (pos[ip2].GetQ()-neg[in].GetQ())/(pos[ip2].GetQ()+neg[in].GetQ()); + milab[3]=(((ip2<<10) + in)<<10) + idet; // pos|neg|det + Int_t info[3] = {pos[ip2].GetNd(),neg[in].GetNd(),fNlayer[fModule]}; + + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + // out-of-diagonal element of covariance matrix + if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012; + else if ( (info[0]>1) && (info[1]>1) ) { + lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2 + lp[3]=0.0065; // 0.08*0.08; //SigmaZ2 + lp[5]=-6.48e-05; + } + else { + lp[2]=4.80e-06; // 0.00219*0.00219 + lp[3]=0.0093; // 0.0964*0.0964; + if (info[0]==1) { + lp[5]=-0.00014; } - else{ - cl2 = new AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(5); - pairs[ip2][in] =5; - if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster - cl2->SetType(6); - pairs[ip2][in] =6; - } + else { + lp[2]=2.79e-06; // 0.0017*0.0017; + lp[3]=0.00935; // 0.967*0.967; + lp[5]=-4.32e-05; + } + } - fDetTypeRec->AddRecPoint(*cl2); + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + + cl2->SetChargeRatio(ratio); + cl2->SetType(5); + fgPairs[ip2*nn+in] =5; + if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster + cl2->SetType(6); + fgPairs[ip2*nn+in] =6; + } + ncl++; } - ncl++; - } - cused1[ip]++; - cused1[ip2]++; - cused2[in]++; - } - } - } + + cused1[ip]++; + cused1[ip2]++; + cused2[in]++; - // - for (Int_t jn=0;jnGetPadCxz(yn, yp, xt, zt); + + xbest=xt; zbest=zt; + + qbest =neg[jn].GetQ(); + + { + Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.}; + mT2L->MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + } lp[4]=qbest; //Q for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; @@ -606,48 +1161,62 @@ FindClustersSSD(Ali1Dcluster* neg, Int_t nn, milab[3]=(((ip<<10) + jn)<<10) + idet; // pos|neg|det Int_t info[3] = {pos[ip].GetNd(),neg[jn].GetNd(),fNlayer[fModule]}; + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + // out-of-diagonal element of covariance matrix + if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012; + else if ( (info[0]>1) && (info[1]>1) ) { + lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2 + lp[3]=0.0065; // 0.08*0.08; //SigmaZ2 + lp[5]=-6.48e-05; + } + else { + lp[2]=4.80e-06; // 0.00219*0.00219 + lp[3]=0.0093; // 0.0964*0.0964; + if (info[0]==1) { + lp[5]=-0.00014; + } + else { + lp[2]=2.79e-06; // 0.0017*0.0017; + lp[3]=0.00935; // 0.967*0.967; + lp[5]=-4.32e-05; + } + } + AliITSRecPoint * cl2; - if(clusters){ - cl2 = new (cl[ncl]) AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(7); - pairs[ip][jn] =7; - if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster - cl2->SetType(8); - pairs[ip][jn]=8; - } + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); - } - else{ - cl2 = new AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); cl2->SetChargeRatio(ratio); cl2->SetType(7); - pairs[ip][jn] =7; + fgPairs[ip*nn+jn] =7; if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster cl2->SetType(8); - pairs[ip][jn]=8; + fgPairs[ip*nn+jn]=8; } - - fDetTypeRec->AddRecPoint(*cl2); - } ncl++; } // // add second pair // if (!(cused1[ip]||cused2[jn2])){ - if (pairs[ip][jn2]==100){ - Float_t yn=neg[jn2].GetY()*fYpitchSSD; - Double_t yp=pos[ip].GetY()*fYpitchSSD; - Double_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp); - Double_t yt=yn + tann*zt; - zt-=fHlSSD; yt-=fHwSSD; - ybest =yt; zbest=zt; - qbest =neg[jn2].GetQ(); - lp[0]=-(-ybest+fYshift[fModule]); - lp[1]= -zbest+fZshift[fModule]; - lp[2]=0.0025*0.0025; //SigmaY2 - lp[3]=0.110*0.110; //SigmaZ2 + if ( (fgPairs[ip*nn+jn2]==100)&&(neg[jn2].GetQ() ) ) { // + + Float_t yn=neg[jn2].GetY(); + Double_t yp=pos[ip].GetY(); + + Float_t xt, zt; + seg->GetPadCxz(yn, yp, xt, zt); + xbest=xt; zbest=zt; + + qbest =neg[jn2].GetQ(); + + { + Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.}; + mT2L->MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + } + lp[4]=qbest; //Q for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; for (Int_t ilab=0;ilab<3;ilab++){ @@ -659,164 +1228,336 @@ FindClustersSSD(Ali1Dcluster* neg, Int_t nn, ratio = (pos[ip].GetQ()-neg[jn2].GetQ())/(pos[ip].GetQ()+neg[jn2].GetQ()); milab[3]=(((ip<<10) + jn2)<<10) + idet; // pos|neg|det Int_t info[3] = {pos[ip].GetNd(),neg[jn2].GetNd(),fNlayer[fModule]}; - AliITSRecPoint * cl2; - if(clusters){ - cl2 = new (cl[ncl]) AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - pairs[ip][jn2]=7; - cl2->SetType(7); - if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster - cl2->SetType(8); - pairs[ip][jn2]=8; - } - + + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + // out-of-diagonal element of covariance matrix + if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012; + else if ( (info[0]>1) && (info[1]>1) ) { + lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2 + lp[3]=0.0065; // 0.08*0.08; //SigmaZ2 + lp[5]=-6.48e-05; + } + else { + lp[2]=4.80e-06; // 0.00219*0.00219 + lp[3]=0.0093; // 0.0964*0.0964; + if (info[0]==1) { + lp[5]=-0.00014; } - else{ - cl2 = new AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); + else { + lp[2]=2.79e-06; // 0.0017*0.0017; + lp[3]=0.00935; // 0.967*0.967; + lp[5]=-4.32e-05; + } + } + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + + cl2->SetChargeRatio(ratio); - pairs[ip][jn2]=7; + fgPairs[ip*nn+jn2]=7; cl2->SetType(7); if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster cl2->SetType(8); - pairs[ip][jn2]=8; + fgPairs[ip*nn+jn2]=8; } - - fDetTypeRec->AddRecPoint(*cl2); - } - ncl++; } cused1[ip]++; cused2[jn]++; cused2[jn2]++; - } - } - } + + } // charge matching condition + + } // 2 Nside cross 1 Pside + } // loop over Pside clusters + - for (Int_t ip=0;ip1) continue; // more than one "proper" cluster for positive + + for (Int_t ip=0;ip1) continue; // more than one "proper" cluster for negative - - Int_t jp = 0; - - count =0; - for (Int_t dj=0;dj1) continue; - if (pairs[ip][j]<100) continue; + // 2x2 clusters // - //almost gold clusters - Float_t yp=pos[ip].GetY()*fYpitchSSD; - Float_t yn=neg[j].GetY()*fYpitchSSD; - Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp); - Float_t yt=yn + tann*zt; - zt-=fHlSSD; yt-=fHwSSD; - ybest=yt; zbest=zt; - qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ()); - lp[0]=-(-ybest+fYshift[fModule]); - lp[1]= -zbest+fZshift[fModule]; - lp[2]=0.0025*0.0025; //SigmaY2 - lp[3]=0.110*0.110; //SigmaZ2 - lp[4]=qbest; //Q - for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; - for (Int_t ilab=0;ilab<3;ilab++){ - milab[ilab] = pos[ip].GetLabel(ilab); - milab[ilab+3] = neg[j].GetLabel(ilab); - } + if ( (cnegative[ip]==2) && cpositive[negativepair[10*ip]]==2){ + Float_t minchargediff =4.; + Float_t minchargeratio =0.2; + + Int_t j=-1; + for (Int_t di=0;di1) continue; // more than one "proper" cluster for positive + // + + count =0; + for (Int_t dj=0;dj1) continue; // more than one "proper" cluster for negative + + Int_t jp = 0; + + count =0; + for (Int_t dj=0;dj1) continue; + if (fgPairs[ip*nn+j]<100) continue; + // + + + + //almost gold clusters + Float_t yp=pos[ip].GetY(); + Float_t yn=neg[j].GetY(); + Float_t xt, zt; + seg->GetPadCxz(yn, yp, xt, zt); + xbest=xt; zbest=zt; + qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ()); + { + Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.}; + mT2L->MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + } + lp[4]=qbest; //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++){ + milab[ilab] = pos[ip].GetLabel(ilab); + milab[ilab+3] = neg[j].GetLabel(ilab); + } + // + CheckLabels2(milab); + if ((neg[j].GetQ()==0)&&(pos[ip].GetQ()==0)) continue; // reject crosses of bad strips!! + ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ()); + milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det + Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]}; + + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + // out-of-diagonal element of covariance matrix + if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012; + else if ( (info[0]>1) && (info[1]>1) ) { + lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2 + lp[3]=0.0065; // 0.08*0.08; //SigmaZ2 + lp[5]=-6.48e-05; + } + else { + lp[2]=4.80e-06; // 0.00219*0.00219 + lp[3]=0.0093; // 0.0964*0.0964; + if (info[0]==1) { + lp[5]=-0.00014; + } + else { + lp[2]=2.79e-06; // 0.0017*0.0017; + lp[3]=0.00935; // 0.967*0.967; + lp[5]=-4.32e-05; + } + } + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + + cl2->SetChargeRatio(ratio); + cl2->SetType(10); + fgPairs[ip*nn+j]=10; + if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster + cl2->SetType(11); + fgPairs[ip*nn+j]=11; + } + cused1[ip]++; + cused2[j]++; + ncl++; + + } // 2X2 + } // loop over Pside 1Dclusters + + + + for (Int_t ip=0;ipGetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(10); - pairs[ip][j]=10; - if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster - cl2->SetType(11); - pairs[ip][j]=11; + // manyxmany clusters + // + if ( (cnegative[ip]<5) && cpositive[negativepair[10*ip]]<5){ + Float_t minchargediff =4.; + Int_t j=-1; + for (Int_t di=0;diGetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(10); - pairs[ip][j]=10; - if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster - cl2->SetType(11); - pairs[ip][j]=11; + if (j<0) continue; // not proper cluster + + Int_t count =0; + for (Int_t di=0;di1) continue; // more than one "proper" cluster for positive + // - fDetTypeRec->AddRecPoint(*cl2); - } - ncl++; - } + count =0; + for (Int_t dj=0;dj1) continue; // more than one "proper" cluster for negative + + Int_t jp = 0; + + count =0; + for (Int_t dj=0;dj1) continue; + if (fgPairs[ip*nn+j]<100) continue; + // + + //almost gold clusters + Float_t yp=pos[ip].GetY(); + Float_t yn=neg[j].GetY(); + - } + Float_t xt, zt; + seg->GetPadCxz(yn, yp, xt, zt); + + xbest=xt; zbest=zt; + + qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ()); + + { + Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.}; + mT2L->MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + } + lp[4]=qbest; //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++){ + milab[ilab] = pos[ip].GetLabel(ilab); + milab[ilab+3] = neg[j].GetLabel(ilab); + } + // + CheckLabels2(milab); + if ((neg[j].GetQ()==0)&&(pos[ip].GetQ()==0)) continue; // reject crosses of bad strips!! + ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ()); + milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det + Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]}; + + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + // out-of-diagonal element of covariance matrix + if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012; + else if ( (info[0]>1) && (info[1]>1) ) { + lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2 + lp[3]=0.0065; // 0.08*0.08; //SigmaZ2 + lp[5]=-6.48e-05; + } + else { + lp[2]=4.80e-06; // 0.00219*0.00219 + lp[3]=0.0093; // 0.0964*0.0964; + if (info[0]==1) { + lp[5]=-0.00014; + } + else { + lp[2]=2.79e-06; // 0.0017*0.0017; + lp[3]=0.00935; // 0.967*0.967; + lp[5]=-4.32e-05; + } + } + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + + cl2->SetChargeRatio(ratio); + cl2->SetType(12); + fgPairs[ip*nn+j]=12; + if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster + cl2->SetType(13); + fgPairs[ip*nn+j]=13; + } + cused1[ip]++; + cused2[j]++; + ncl++; + + } // manyXmany + } // loop over Pside 1Dclusters + + } // use charge matching + + // recover all the other crosses // for (Int_t i=0; i0)&&(pos[i].GetQ()<3)) continue; for (Int_t j=0; j0)&&(neg[j].GetQ()<3)) continue; + + if ((neg[j].GetQ()==0)&&(pos[i].GetQ()==0)) continue; // reject crosses of bad strips!! + if (cused2[j]||cused1[i]) continue; - if (pairs[i][j]>0 &&pairs[i][j]<100) continue; + if (fgPairs[i*nn+j]>0 &&fgPairs[i*nn+j]<100) continue; ratio = (pos[i].GetQ()-neg[j].GetQ())/(pos[i].GetQ()+neg[j].GetQ()); - Float_t yn=neg[j].GetY()*fYpitchSSD; - Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp); - Float_t yt=yn + tann*zt; - zt-=fHlSSD; yt-=fHwSSD; - if (TMath::Abs(yt)GetPadCxz(yn, yp, xt, zt); + + if (TMath::Abs(xt)MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + } lp[4]=qbest; //Q for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; for (Int_t ilab=0;ilab<3;ilab++){ @@ -827,36 +1568,258 @@ FindClustersSSD(Ali1Dcluster* neg, Int_t nn, CheckLabels2(milab); milab[3]=(((i<<10) + j)<<10) + idet; // pos|neg|det Int_t info[3] = {pos[i].GetNd(),neg[j].GetNd(),fNlayer[fModule]}; + + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + // out-of-diagonal element of covariance matrix + if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012; + else if ( (info[0]>1) && (info[1]>1) ) { + lp[2]=2.63e-06; // 0.0016*0.0016; //SigmaY2 + lp[3]=0.0065; // 0.08*0.08; //SigmaZ2 + lp[5]=-6.48e-05; + } + else { + lp[2]=4.80e-06; // 0.00219*0.00219 + lp[3]=0.0093; // 0.0964*0.0964; + if (info[0]==1) { + lp[5]=-0.00014; + } + else { + lp[2]=2.79e-06; // 0.0017*0.0017; + lp[3]=0.00935; // 0.967*0.967; + lp[5]=-4.32e-05; + } + } + AliITSRecPoint * cl2; - if(clusters){ - cl2 = new (cl[ncl]) AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + cl2->SetChargeRatio(ratio); cl2->SetType(100+cpositive[j]+cnegative[i]); - } - else{ - cl2 = new AliITSRecPoint(fModule,fDetTypeRec->GetITSgeom(),milab,lp,info); - cl2->SetChargeRatio(ratio); - cl2->SetType(100+cpositive[j]+cnegative[i]); - fDetTypeRec->AddRecPoint(*cl2); - } + + if(pos[i].GetQ()==0) cl2->SetType(200+cpositive[j]+cnegative[i]); + if(neg[j].GetQ()==0) cl2->SetType(300+cpositive[j]+cnegative[i]); ncl++; - //cl2->SetType(0); - /* - if (pairs[i][j]<100){ - printf("problem:- %d\n", pairs[i][j]); - } - if (cnegative[i]<2&&cpositive[j]<2){ - printf("problem:- %d\n", pairs[i][j]); - } - */ } } } -// for (Int_t i=0; i<1000; i++) delete [] pairs[i]; -// delete [] pairs; -} + if(repa->GetUseBadChannelsInClusterFinderSSD()==kTRUE) { + + //--------------------------------------------------------- + // recover crosses of good 1D clusters with bad strips on the other side + // Note1: at first iteration skip modules with a bad side (or almost), (would produce too many fake!) + // Note2: for modules with a bad side see below + + AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)GetResp(fModule); + Int_t countPbad=0, countNbad=0; + for(Int_t ib=0; ib<768; ib++) { + if(cal->IsPChannelBad(ib)) countPbad++; + if(cal->IsNChannelBad(ib)) countNbad++; + } + // AliInfo(Form("module %d has %d P- and %d N-bad strips",fModule,countPbad,countNbad)); + + if( (countPbad<100) && (countNbad<100) ) { // no bad side!! + + for (Int_t i=0; iIsPChannelBad(ib)) { // check if strips is bad + Float_t yN=pos[i].GetY(); + Float_t xt, zt; + seg->GetPadCxz(1.*ib, yN, xt, zt); + + //---------- + // bad Pstrip is crossing the Nside 1Dcluster -> create recpoint + // + if ( (TMath::Abs(xt)MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + lp[4]=pos[i].GetQ(); //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = pos[i].GetLabel(ilab); + CheckLabels2(milab); + milab[3]=( (i<<10) << 10 ) + idet; // pos|neg|det + Int_t info[3] = {pos[i].GetNd(),0,fNlayer[fModule]}; + + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + lp[5]=-0.00012; // out-of-diagonal element of covariance matrix + if (info[0]>1) { + lp[2]=4.80e-06; + lp[3]=0.0093; + lp[5]=0.00014; + } + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + cl2->SetChargeRatio(1.); + cl2->SetType(50); + ncl++; + } // cross is within the detector + // + //-------------- + + } // bad Pstrip + + } // end loop over Pstrips + + } // end loop over Nside 1D clusters + + for (Int_t j=0; jIsNChannelBad(ib)) { // check if strip is bad + Float_t yP=neg[j].GetY(); + Float_t xt, zt; + seg->GetPadCxz(yP, 1.*ib, xt, zt); + + //---------- + // bad Nstrip is crossing the Pside 1Dcluster -> create recpoint + // + if ( (TMath::Abs(xt)MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + lp[4]=neg[j].GetQ(); //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = neg[j].GetLabel(ilab); + CheckLabels2(milab); + milab[3]=( j << 10 ) + idet; // pos|neg|det + Int_t info[3]={0,(Int_t)neg[j].GetNd(),fNlayer[fModule]}; + + lp[2]=4.968e-06; // 0.00223*0.00223; //SigmaY2 + lp[3]=0.012; // 0.110*0.110; //SigmaZ2 + lp[5]=-0.00012; // out-of-diagonal element of covariance matrix + if (info[0]>1) { + lp[2]=2.79e-06; + lp[3]=0.00935; + lp[5]=-4.32e-05; + } + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + cl2->SetChargeRatio(1.); + cl2->SetType(60); + ncl++; + } // cross is within the detector + // + //-------------- + + } // bad Nstrip + } // end loop over Nstrips + } // end loop over Pside 1D clusters + + } // no bad sides + + //--------------------------------------------------------- + + else if( (countPbad>700) && (countNbad<100) ) { // bad Pside!! + + for (Int_t i=0; iGetLayer()==5) yP = yN + (7.6/1.9); + else yP = yN - (7.6/1.9); + seg->GetPadCxz(yP, yN, xt, zt); + + if ( (TMath::Abs(xt)MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + lp[4]=pos[i].GetQ(); //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = pos[i].GetLabel(ilab); + CheckLabels2(milab); + milab[3]=( (i<<10) << 10 ) + idet; // pos|neg|det + Int_t info[3] = {(Int_t)pos[i].GetNd(),0,fNlayer[fModule]}; + + lp[2]=0.00098; // 0.031*0.031; //SigmaY2 + lp[3]=1.329; // 1.15*1.15; //SigmaZ2 + lp[5]=-0.0359; + if(info[0]>1) lp[2]=0.00097; + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + cl2->SetChargeRatio(1.); + cl2->SetType(70); + ncl++; + } // cross is within the detector + // + //-------------- + + } // end loop over Nside 1D clusters + + } // bad Pside module + + else if( (countNbad>700) && (countPbad<100) ) { // bad Nside!! + + for (Int_t j=0; jGetLayer()==5) yN = yP - (7.6/1.9); + else yN = yP + (7.6/1.9); + seg->GetPadCxz(yP, yN, xt, zt); + + if ( (TMath::Abs(xt)MasterToLocal(loc,trk); + lp[0]=trk[1]; + lp[1]=trk[2]; + lp[4]=neg[j].GetQ(); //Q + for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2; + for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = neg[j].GetLabel(ilab); + CheckLabels2(milab); + milab[3]=( j << 10 ) + idet; // pos|neg|det + Int_t info[3] = {0,(Int_t)neg[j].GetNd(),fNlayer[fModule]}; + + lp[2]=7.27e-05; // 0.0085*0.0085; //SigmaY2 + lp[3]=1.33; // 1.15*1.15; //SigmaZ2 + lp[5]=0.00931; + if(info[1]>1) lp[2]=6.91e-05; + + AliITSRecPoint * cl2; + cl2 = new ((*clusters)[ncl]) AliITSRecPoint(milab,lp,info); + cl2->SetChargeRatio(1.); + cl2->SetType(80); + ncl++; + } // cross is within the detector + // + //-------------- + + } // end loop over Pside 1D clusters + + } // bad Nside module + + //--------------------------------------------------------- + + } // use bad channels + + //cout<