+//------------------------------------------------------------------------
+void AliITStrackerMI::BuildMaterialLUT(TString material) {
+ //--------------------------------------------------------------------
+ // Fill a look-up table with mean material
+ //--------------------------------------------------------------------
+
+ Int_t n=1000;
+ Double_t mparam[7];
+ Double_t point1[3],point2[3];
+ Double_t phi,cosphi,sinphi,z;
+ // 0-5 layers, 6 pipe, 7-8 shields
+ Double_t rmin[9]={ 3.5, 5.5,13.0,22.0,35.0,41.0, 2.0, 7.5,25.0};
+ Double_t rmax[9]={ 5.5, 7.3,17.0,26.0,41.0,47.0, 3.0,10.5,30.0};
+
+ Int_t ifirst=0,ilast=0;
+ if(material.Contains("Pipe")) {
+ ifirst=6; ilast=6;
+ } else if(material.Contains("Shields")) {
+ ifirst=7; ilast=8;
+ } else if(material.Contains("Layers")) {
+ ifirst=0; ilast=5;
+ } else {
+ Error("BuildMaterialLUT","Wrong layer name\n");
+ }
+
+ for(Int_t imat=ifirst; imat<=ilast; imat++) {
+ Double_t param[5]={0.,0.,0.,0.,0.};
+ for (Int_t i=0; i<n; i++) {
+ phi = 2.*TMath::Pi()*gRandom->Rndm();
+ cosphi = TMath::Cos(phi); sinphi = TMath::Sin(phi);
+ z = 14.*(-1.+2.*gRandom->Rndm()); // SPD barrel
+ point1[0] = rmin[imat]*cosphi;
+ point1[1] = rmin[imat]*sinphi;
+ point1[2] = z;
+ point2[0] = rmax[imat]*cosphi;
+ point2[1] = rmax[imat]*sinphi;
+ point2[2] = z;
+ AliTracker::MeanMaterialBudget(point1,point2,mparam);
+ for(Int_t j=0;j<5;j++) param[j]+=mparam[j];
+ }
+ for(Int_t j=0;j<5;j++) param[j]/=(Float_t)n;
+ if(imat<=5) {
+ fxOverX0Layer[imat] = param[1];
+ fxTimesRhoLayer[imat] = param[0]*param[4];
+ } else if(imat==6) {
+ fxOverX0Pipe = param[1];
+ fxTimesRhoPipe = param[0]*param[4];
+ } else if(imat==7) {
+ fxOverX0Shield[0] = param[1];
+ fxTimesRhoShield[0] = param[0]*param[4];
+ } else if(imat==8) {
+ fxOverX0Shield[1] = param[1];
+ fxTimesRhoShield[1] = param[0]*param[4];
+ }
+ }
+ /*
+ printf("%s\n",material.Data());
+ printf("%f %f\n",fxOverX0Pipe,fxTimesRhoPipe);
+ printf("%f %f\n",fxOverX0Shield[0],fxTimesRhoShield[0]);
+ printf("%f %f\n",fxOverX0Shield[1],fxTimesRhoShield[1]);
+ printf("%f %f\n",fxOverX0Layer[0],fxTimesRhoLayer[0]);
+ printf("%f %f\n",fxOverX0Layer[1],fxTimesRhoLayer[1]);
+ printf("%f %f\n",fxOverX0Layer[2],fxTimesRhoLayer[2]);
+ printf("%f %f\n",fxOverX0Layer[3],fxTimesRhoLayer[3]);
+ printf("%f %f\n",fxOverX0Layer[4],fxTimesRhoLayer[4]);
+ printf("%f %f\n",fxOverX0Layer[5],fxTimesRhoLayer[5]);
+ */
+ return;
+}
+//------------------------------------------------------------------------
+Int_t AliITStrackerMI::CorrectForPipeMaterial(AliITStrackMI *t,
+ TString direction) {
+ //-------------------------------------------------------------------
+ // Propagate beyond beam pipe and correct for material
+ // (material budget in different ways according to fUseTGeo value)
+ //-------------------------------------------------------------------
+
+ // Define budget mode:
+ // 0: material from AliITSRecoParam (hard coded)
+ // 1: material from TGeo (on the fly)
+ // 2: material from lut
+ // 3: material from TGeo (same for all hypotheses)
+ Int_t mode;
+ switch(fUseTGeo) {
+ case 0:
+ mode=0;
+ break;
+ case 1:
+ mode=1;
+ break;
+ case 2:
+ mode=2;
+ break;
+ case 3:
+ if(fTrackingPhase.Contains("Clusters2Tracks"))
+ { mode=3; } else { mode=1; }
+ break;
+ case 4:
+ if(fTrackingPhase.Contains("Clusters2Tracks"))
+ { mode=3; } else { mode=2; }
+ break;
+ default:
+ mode=0;
+ break;
+ }
+ if(fTrackingPhase.Contains("Default")) mode=0;
+
+ Int_t index=fCurrentEsdTrack;
+
+ Float_t dir = (direction.Contains("inward") ? 1. : -1.);
+ Double_t rToGo=(dir>0 ? AliITSRecoParam::GetrInsidePipe() : AliITSRecoParam::GetrOutsidePipe());
+ Double_t xToGo; t->GetLocalXat(rToGo,xToGo);
+
+ Double_t xOverX0,x0,lengthTimesMeanDensity;
+ Bool_t anglecorr=kTRUE;
+
+ switch(mode) {
+ case 0:
+ xOverX0 = AliITSRecoParam::GetdPipe();
+ x0 = AliITSRecoParam::GetX0Be();
+ lengthTimesMeanDensity = xOverX0*x0;
+ break;
+ case 1:
+ if (!t->PropagateToTGeo(xToGo,1)) return 0;
+ return 1;
+ break;
+ case 2:
+ if(fxOverX0Pipe<0) BuildMaterialLUT("Pipe");
+ xOverX0 = fxOverX0Pipe;
+ lengthTimesMeanDensity = fxTimesRhoPipe;
+ break;
+ case 3:
+ if(!fxOverX0PipeTrks || index<0 || index>=fNtracks) Error("CorrectForPipeMaterial","Incorrect usage of UseTGeo option!\n");
+ if(fxOverX0PipeTrks[index]<0) {
+ if (!t->PropagateToTGeo(xToGo,1,xOverX0,lengthTimesMeanDensity)) return 0;
+ Double_t angle=TMath::Sqrt((1.+t->GetTgl()*t->GetTgl())/
+ (1.-t->GetSnp()*t->GetSnp()));
+ fxOverX0PipeTrks[index] = TMath::Abs(xOverX0)/angle;
+ fxTimesRhoPipeTrks[index] = TMath::Abs(lengthTimesMeanDensity)/angle;
+ return 1;
+ }
+ xOverX0 = fxOverX0PipeTrks[index];
+ lengthTimesMeanDensity = fxTimesRhoPipeTrks[index];
+ break;
+ }
+
+ lengthTimesMeanDensity *= dir;
+
+ if (!t->AliExternalTrackParam::PropagateTo(xToGo,GetBz())) return 0;
+ if (!t->CorrectForMeanMaterial(xOverX0,lengthTimesMeanDensity,anglecorr)) return 0;
+
+ return 1;
+}
+//------------------------------------------------------------------------
+Int_t AliITStrackerMI::CorrectForShieldMaterial(AliITStrackMI *t,
+ TString shield,
+ TString direction) {
+ //-------------------------------------------------------------------
+ // Propagate beyond SPD or SDD shield and correct for material
+ // (material budget in different ways according to fUseTGeo value)
+ //-------------------------------------------------------------------
+
+ // Define budget mode:
+ // 0: material from AliITSRecoParam (hard coded)
+ // 1: material from TGeo (on the fly)
+ // 2: material from lut
+ // 3: material from TGeo (same for all hypotheses)
+ Int_t mode;
+ switch(fUseTGeo) {
+ case 0:
+ mode=0;
+ break;
+ case 1:
+ mode=1;
+ break;
+ case 2:
+ mode=2;
+ break;
+ case 3:
+ if(fTrackingPhase.Contains("Clusters2Tracks"))
+ { mode=3; } else { mode=1; }
+ break;
+ case 4:
+ if(fTrackingPhase.Contains("Clusters2Tracks"))
+ { mode=3; } else { mode=2; }
+ break;
+ default:
+ mode=0;
+ break;
+ }
+ if(fTrackingPhase.Contains("Default")) mode=0;
+
+ Float_t dir = (direction.Contains("inward") ? 1. : -1.);
+ Double_t rToGo;
+ Int_t shieldindex=0;
+ if (shield.Contains("SDD")) { // SDDouter
+ rToGo=(dir>0 ? AliITSRecoParam::GetrInsideShield(1) : AliITSRecoParam::GetrOutsideShield(1));
+ shieldindex=1;
+ } else if (shield.Contains("SPD")) { // SPDouter
+ rToGo=(dir>0 ? AliITSRecoParam::GetrInsideShield(0) : AliITSRecoParam::GetrOutsideShield(0));
+ shieldindex=0;
+ } else {
+ Error("CorrectForShieldMaterial"," Wrong shield name\n");
+ return 0;
+ }
+ Double_t xToGo; t->GetLocalXat(rToGo,xToGo);
+
+ Int_t index=2*fCurrentEsdTrack+shieldindex;
+
+ Double_t xOverX0,x0,lengthTimesMeanDensity;
+ Bool_t anglecorr=kTRUE;
+
+ switch(mode) {
+ case 0:
+ xOverX0 = AliITSRecoParam::Getdshield(shieldindex);
+ x0 = AliITSRecoParam::GetX0shield(shieldindex);
+ lengthTimesMeanDensity = xOverX0*x0;
+ break;
+ case 1:
+ if (!t->PropagateToTGeo(xToGo,1)) return 0;
+ return 1;
+ break;
+ case 2:
+ if(fxOverX0Shield[shieldindex]<0) BuildMaterialLUT("Shields");
+ xOverX0 = fxOverX0Shield[shieldindex];
+ lengthTimesMeanDensity = fxTimesRhoShield[shieldindex];
+ break;
+ case 3:
+ if(!fxOverX0ShieldTrks || index<0 || index>=2*fNtracks) Error("CorrectForShieldMaterial","Incorrect usage of UseTGeo option!\n");
+ if(fxOverX0ShieldTrks[index]<0) {
+ if (!t->PropagateToTGeo(xToGo,1,xOverX0,lengthTimesMeanDensity)) return 0;
+ Double_t angle=TMath::Sqrt((1.+t->GetTgl()*t->GetTgl())/
+ (1.-t->GetSnp()*t->GetSnp()));
+ fxOverX0ShieldTrks[index] = TMath::Abs(xOverX0)/angle;
+ fxTimesRhoShieldTrks[index] = TMath::Abs(lengthTimesMeanDensity)/angle;
+ return 1;
+ }
+ xOverX0 = fxOverX0ShieldTrks[index];
+ lengthTimesMeanDensity = fxTimesRhoShieldTrks[index];
+ break;
+ }
+
+ lengthTimesMeanDensity *= dir;
+
+ if (!t->AliExternalTrackParam::PropagateTo(xToGo,GetBz())) return 0;
+ if (!t->CorrectForMeanMaterial(xOverX0,lengthTimesMeanDensity,anglecorr)) return 0;
+
+ return 1;
+}
+//------------------------------------------------------------------------
+Int_t AliITStrackerMI::CorrectForLayerMaterial(AliITStrackMI *t,
+ Int_t layerindex,
+ Double_t oldGlobXYZ[3],
+ TString direction) {
+ //-------------------------------------------------------------------
+ // Propagate beyond layer and correct for material
+ // (material budget in different ways according to fUseTGeo value)
+ //-------------------------------------------------------------------
+
+ // Define budget mode:
+ // 0: material from AliITSRecoParam (hard coded)
+ // 1: material from TGeo (on the fly)
+ // 2: material from lut
+ // 3: material from TGeo (same for all hypotheses)
+ Int_t mode;
+ switch(fUseTGeo) {
+ case 0:
+ mode=0;
+ break;
+ case 1:
+ mode=1;
+ break;
+ case 2:
+ mode=2;
+ break;
+ case 3:
+ if(fTrackingPhase.Contains("Clusters2Tracks"))
+ { mode=3; } else { mode=1; }
+ break;
+ case 4:
+ if(fTrackingPhase.Contains("Clusters2Tracks"))
+ { mode=3; } else { mode=2; }
+ break;
+ default:
+ mode=0;
+ break;
+ }
+ if(fTrackingPhase.Contains("Default")) mode=0;
+
+ Float_t dir = (direction.Contains("inward") ? 1. : -1.);
+
+ Double_t r=fgLayers[layerindex].GetR();
+ Double_t deltar=(layerindex<2 ? 0.10*r : 0.05*r);
+
+ Double_t rToGo=TMath::Sqrt(t->GetX()*t->GetX()+t->GetY()*t->GetY())-deltar*dir;
+ Double_t xToGo; t->GetLocalXat(rToGo,xToGo);
+
+ Int_t index=6*fCurrentEsdTrack+layerindex;
+
+ // Bring the track beyond the material
+ if (!t->AliExternalTrackParam::PropagateTo(xToGo,GetBz())) return 0;
+ Double_t globXYZ[3];
+ t->GetXYZ(globXYZ);
+
+ Double_t xOverX0=0.0,x0=0.0,lengthTimesMeanDensity=0.0;
+ Double_t mparam[7];
+ Bool_t anglecorr=kTRUE;
+
+ switch(mode) {
+ case 0:
+ xOverX0 = fgLayers[layerindex].GetThickness(t->GetY(),t->GetZ(),x0);
+ lengthTimesMeanDensity = xOverX0*x0;
+ break;
+ case 1:
+ AliTracker::MeanMaterialBudget(oldGlobXYZ,globXYZ,mparam);
+ if(mparam[1]>900000) return 0;
+ xOverX0=mparam[1];
+ lengthTimesMeanDensity=mparam[0]*mparam[4];
+ anglecorr=kFALSE;
+ break;
+ case 2:
+ if(fxOverX0Layer[layerindex]<0) BuildMaterialLUT("Layers");
+ xOverX0 = fxOverX0Layer[layerindex];
+ lengthTimesMeanDensity = fxTimesRhoLayer[layerindex];
+ break;
+ case 3:
+ if(!fxOverX0LayerTrks || index<0 || index>=6*fNtracks) Error("CorrectForLayerMaterial","Incorrect usage of UseTGeo option!\n");
+ if(fxOverX0LayerTrks[index]<0) {
+ AliTracker::MeanMaterialBudget(oldGlobXYZ,globXYZ,mparam);
+ if(mparam[1]>900000) return 0;
+ Double_t angle=TMath::Sqrt((1.+t->GetTgl()*t->GetTgl())/
+ (1.-t->GetSnp()*t->GetSnp()));
+ xOverX0=mparam[1]/angle;
+ lengthTimesMeanDensity=mparam[0]*mparam[4]/angle;
+ fxOverX0LayerTrks[index] = TMath::Abs(xOverX0);
+ fxTimesRhoLayerTrks[index] = TMath::Abs(lengthTimesMeanDensity);
+ }
+ xOverX0 = fxOverX0LayerTrks[index];
+ lengthTimesMeanDensity = fxTimesRhoLayerTrks[index];
+ break;
+ }
+
+ lengthTimesMeanDensity *= dir;
+
+ if (!t->CorrectForMeanMaterial(xOverX0,lengthTimesMeanDensity,anglecorr)) return 0;
+
+ return 1;
+}
+//------------------------------------------------------------------------
+void AliITStrackerMI::MakeTrksMaterialLUT(Int_t ntracks) {
+ //-----------------------------------------------------------------
+ // Initialize LUT for storing material for each prolonged track
+ //-----------------------------------------------------------------
+ fxOverX0PipeTrks = new Float_t[ntracks];
+ fxTimesRhoPipeTrks = new Float_t[ntracks];
+ fxOverX0ShieldTrks = new Float_t[ntracks*2];
+ fxTimesRhoShieldTrks = new Float_t[ntracks*2];
+ fxOverX0LayerTrks = new Float_t[ntracks*6];
+ fxTimesRhoLayerTrks = new Float_t[ntracks*6];
+
+ for(Int_t i=0; i<ntracks; i++) {
+ fxOverX0PipeTrks[i] = -1.;
+ fxTimesRhoPipeTrks[i] = -1.;
+ }
+ for(Int_t j=0; j<ntracks*2; j++) {
+ fxOverX0ShieldTrks[j] = -1.;
+ fxTimesRhoShieldTrks[j] = -1.;
+ }
+ for(Int_t k=0; k<ntracks*6; k++) {
+ fxOverX0LayerTrks[k] = -1.;
+ fxTimesRhoLayerTrks[k] = -1.;
+ }
+
+ fNtracks = ntracks;
+
+ return;
+}
+//------------------------------------------------------------------------
+void AliITStrackerMI::DeleteTrksMaterialLUT() {
+ //-----------------------------------------------------------------
+ // Delete LUT for storing material for each prolonged track
+ //-----------------------------------------------------------------
+ if(fxOverX0PipeTrks) {
+ delete [] fxOverX0PipeTrks; fxOverX0PipeTrks = 0;
+ }
+ if(fxOverX0ShieldTrks) {
+ delete [] fxOverX0ShieldTrks; fxOverX0ShieldTrks = 0;
+ }
+
+ if(fxOverX0LayerTrks) {
+ delete [] fxOverX0LayerTrks; fxOverX0LayerTrks = 0;
+ }
+ if(fxTimesRhoPipeTrks) {
+ delete [] fxTimesRhoPipeTrks; fxTimesRhoPipeTrks = 0;
+ }
+ if(fxTimesRhoShieldTrks) {
+ delete [] fxTimesRhoShieldTrks; fxTimesRhoShieldTrks = 0;
+ }
+ if(fxTimesRhoLayerTrks) {
+ delete [] fxTimesRhoLayerTrks; fxTimesRhoLayerTrks = 0;
+ }
+ return;
+}
+//------------------------------------------------------------------------
+Int_t AliITStrackerMI::CheckSkipLayer(AliITStrackMI *track,
+ Int_t ilayer,Int_t idet) const {
+ //-----------------------------------------------------------------
+ // This method is used to decide whether to allow a prolongation
+ // without clusters, because we want to skip the layer.
+ // In this case the return value is > 0:
+ // return 1: the user requested to skip a layer
+ // return 2: track outside z acceptance of SSD/SDD and will cross both SPD
+ //-----------------------------------------------------------------
+
+ if (AliITSReconstructor::GetRecoParam()->GetLayersToSkip(ilayer)) return 1;
+
+ if (idet<0 && ilayer>1 && AliITSReconstructor::GetRecoParam()->GetExtendedEtaAcceptance()) {
+ // check if track will cross SPD outer layer
+ Double_t phiAtSPD2,zAtSPD2;
+ if (track->GetPhiZat(fgLayers[1].GetR(),phiAtSPD2,zAtSPD2)) {
+ if (TMath::Abs(zAtSPD2)<2.*AliITSRecoParam::GetSPDdetzlength()) return 2;
+ }
+ }
+
+ return 0;
+}
+//------------------------------------------------------------------------
+Int_t AliITStrackerMI::CheckDeadZone(/*AliITStrackMI *track,*/
+ Int_t ilayer,Int_t idet,
+ Double_t zmin,Double_t zmax/*,Double_t ymin,Double_t ymax*/) const {
+ //-----------------------------------------------------------------
+ // This method is used to decide whether to allow a prolongation
+ // without clusters, because there is a dead zone in the road.
+ // In this case the return value is > 0:
+ // return 1: dead zone at z=0,+-7cm in SPD
+ // return 2: dead area from the OCDB // NOT YET IMPLEMENTED
+ //-----------------------------------------------------------------
+
+ // check dead zones at z=0,+-7cm in the SPD
+ if (ilayer<2 && !AliITSReconstructor::GetRecoParam()->GetAddVirtualClustersInDeadZone()) {
+ Double_t zmindead[3]={fSPDdetzcentre[0] + 0.5*AliITSRecoParam::GetSPDdetzlength(),
+ fSPDdetzcentre[1] + 0.5*AliITSRecoParam::GetSPDdetzlength(),
+ fSPDdetzcentre[2] + 0.5*AliITSRecoParam::GetSPDdetzlength()};
+ Double_t zmaxdead[3]={fSPDdetzcentre[1] - 0.5*AliITSRecoParam::GetSPDdetzlength(),
+ fSPDdetzcentre[2] - 0.5*AliITSRecoParam::GetSPDdetzlength(),
+ fSPDdetzcentre[3] - 0.5*AliITSRecoParam::GetSPDdetzlength()};
+ for (Int_t i=0; i<3; i++)
+ if (zmin<zmaxdead[i] && zmax>zmindead[i]) return 1;
+ }
+
+ // check dead zones from OCDB
+ if (!AliITSReconstructor::GetRecoParam()->GetUseDeadZonesFromOCDB()) return 0;
+
+ if(idet<0) return 0;
+
+ // look in OCDB (only entire dead modules for the moment)
+ if (ilayer==0 || ilayer==1) { // SPD
+ AliCDBEntry* cdbEntry = AliCDBManager::Instance()->Get("ITS/Calib/SPDDead");
+ if (!cdbEntry) {
+ Error("CheckDeadZone","Cannot get CDB entry for SPD\n");
+ return 0;
+ }
+ TObjArray* spdEntry = (TObjArray*)cdbEntry->GetObject();
+ if (!spdEntry) {
+ Error("CheckDeadZone","Cannot get CDB entry for SPD\n");
+ return 0;
+ }
+ if(ilayer==1) idet += AliITSgeomTGeo::GetNLadders(1)*AliITSgeomTGeo::GetNDetectors(1);
+ //printf("SPD det: %d\n",idet);
+ AliITSCalibrationSPD *calibSPD = (AliITSCalibrationSPD*)spdEntry->At(idet);
+ if (calibSPD->IsBad()) return 2;
+ } else if (ilayer==2 || ilayer==3) { // SDD
+ AliCDBEntry* cdbEntry = AliCDBManager::Instance()->Get("ITS/Calib/CalibSDD");
+ if (!cdbEntry) {
+ Error("CheckDeadZone","Cannot get CDB entry for SDD\n");
+ return 0;
+ }
+ TObjArray* sddEntry = (TObjArray*)cdbEntry->GetObject();
+ if (!sddEntry) {
+ Error("CheckDeadZone","Cannot get CDB entry for SDD\n");
+ return 0;
+ }
+ if(ilayer==3) idet += AliITSgeomTGeo::GetNLadders(3)*AliITSgeomTGeo::GetNDetectors(3);
+ //printf("SDD det: %d\n",idet);
+ AliITSCalibrationSDD *calibSDD = (AliITSCalibrationSDD*)sddEntry->At(idet);
+ if (calibSDD->IsDead()) return 2;
+ } else if (ilayer==4 || ilayer==5) { // SSD
+ } else {
+ Error("CheckDeadZone","Wrong layer number\n");
+ if(ilayer==5) idet += AliITSgeomTGeo::GetNLadders(5)*AliITSgeomTGeo::GetNDetectors(5);
+ return 0;
+ }
+
+ return 0;
+}
+//------------------------------------------------------------------------
+Bool_t AliITStrackerMI::LocalModuleCoord(Int_t ilayer,Int_t idet,
+ AliITStrackMI *track,
+ Float_t &xloc,Float_t &zloc) const {
+ //-----------------------------------------------------------------
+ // Gives position of track in local module ref. frame
+ //-----------------------------------------------------------------
+
+ xloc=0.;
+ zloc=0.;