gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
AliXRDPROOFtoolkit tool;
- AliXRDPROOFtoolkit::FilterList("timeitstpc.txt","* itstpc",1)
- AliXRDPROOFtoolkit::FilterList("timetoftpc.txt","* pointMatch",1)
- AliXRDPROOFtoolkit::FilterList("time.txt","* trackInfo",1)
- AliXRDPROOFtoolkit::FilterList("timelaser.txt","* laserInfo",1)
+ TChain * chainTime = tool.MakeChainRandom("time.txt","trackInfo",0,10000);
+ AliXRDPROOFtoolkit::FilterList("timetpctpc.txt","* tpctpc",1)
+ AliXRDPROOFtoolkit::FilterList("timetoftpc.txt","* toftpc",1)
+ AliXRDPROOFtoolkit::FilterList("timeitstpc.txt","* itstpc",1)
+ AliXRDPROOFtoolkit::FilterList("timelaser.txt","* laserInfo",1)
+ TChain * chainTPCTPC = tool.MakeChainRandom("timetpctpc.txt.Good","tpctpc",0,10000);
TChain * chainTPCITS = tool.MakeChainRandom("timeitstpc.txt.Good","itstpc",0,10000);
- TChain * chainTPCTOF = tool.MakeChainRandom("timetoftpc.txt.Good","pointMatch",0,500);
- TChain * chainTime = tool.MakeChainRandom("time.txt.Good","trackInfo",0,10000);
+ TChain * chainTPCTOF = tool.MakeChainRandom("timetoftpc.txt.Good","toftpc",0,10000);
TChain * chainLaser = tool.MakeChainRandom("timelaser.txt.Good","laserInfo",0,10000);
chainTime->Lookup();
chainLaser->Lookup();
#include "AliTracker.h"
#include "AliMagF.h"
#include "AliTPCCalROC.h"
+#include "AliTPCParam.h"
#include "AliLog.h"
fAlignITSTPC = new TObjArray; //alignemnt array ITS TPC match
fAlignTRDTPC = new TObjArray; //alignemnt array ITS TPC match
fAlignTOFTPC = new TObjArray; //alignemnt array ITS TPC match
-
+ fAlignITSTPC->SetOwner(kTRUE);
+ fAlignTRDTPC->SetOwner(kTRUE);
+ fAlignTOFTPC->SetOwner(kTRUE);
+
// fArrayDz->AddLast(fHistVdriftLaserA[0]);
// fArrayDz->AddLast(fHistVdriftLaserA[1]);
// fArrayDz->AddLast(fHistVdriftLaserA[2]);
fCosmiMatchingHisto[i]=NULL;
}
}
+ fAlignITSTPC->SetOwner(kTRUE);
+ fAlignTRDTPC->SetOwner(kTRUE);
+ fAlignTOFTPC->SetOwner(kTRUE);
+
fAlignITSTPC->Delete();
fAlignTRDTPC->Delete();
fAlignTOFTPC->Delete();
+ delete fAlignITSTPC;
+ delete fAlignTRDTPC;
+ delete fAlignTOFTPC;
}
Bool_t AliTPCcalibTime::IsLaser(AliESDEvent */*event*/){
if (!trackOut) continue;
AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
+ if (friendTrack) ProcessSame(track,friendTrack,event);
if (friendTrack) ProcessAlignITS(track,friendTrack);
if (friendTrack) ProcessAlignTRD(track,friendTrack);
if (friendTrack) ProcessAlignTOF(track,friendTrack);
//track 1 lower part
if (!track1) continue;
if (!track1->GetOuterParam()) continue;
- if (track1->GetOuterParam()->GetAlpha()>0) continue;
+ // if (track1->GetOuterParam()->GetAlpha()>0) continue;
//
Double_t d2[3];
track1->GetDirection(d2);
//
// conservative cuts - convergence to be guarantied
// applying before track propagation
- if (TMath::Abs(dist0+dist1)>fCutMaxD) continue; // distance to the 0,0
- if (dir>fCutMinDir) continue; // direction vector product
+ if (TMath::Abs(TMath::Abs(dist0)-TMath::Abs(dist1))>fCutMaxD) continue; // distance to the 0,0
+ if (TMath::Abs(dir)<TMath::Abs(fCutMinDir)) continue; // direction vector product
Float_t bz = AliTracker::GetBz();
Float_t dvertex0[2]; //distance to 0,0
Float_t dvertex1[2]; //distance to 0,0
delete hist;
hist=NULL;
- if(isSame || (isCross && isPair)){
- if (track0->GetTPCNcls() > 80) {
+ if((isSame) || (isCross && isPair)){
+ if (track0->GetTPCNcls()+ track1->GetTPCNcls()> 80) {
fDz = param0.GetZ() - param1.GetZ();
if(track0->GetOuterParam()->GetZ()<0) fDz=-fDz;
TTimeStamp tstamp(fTime);
// localHist->Add(addHist);
// }
for(Int_t i=0;i<10;i++) if (cal->GetCosmiMatchingHisto(i)) fCosmiMatchingHisto[i]->Add(cal->GetCosmiMatchingHisto(i));
+ //
+ // Merge alignment
+ //
+ for (Int_t itype=0; itype<3; itype++){
+ //
+ //
+ TObjArray *arr0= 0;
+ TObjArray *arr1= 0;
+ if (itype==0) {arr0=fAlignITSTPC; arr1=cal->fAlignITSTPC;}
+ if (itype==1) {arr0=fAlignTRDTPC; arr1=cal->fAlignTRDTPC;}
+ if (itype==2) {arr0=fAlignTOFTPC; arr1=cal->fAlignTOFTPC;}
+ if (!arr1) continue;
+ if (!arr0) arr0=new TObjArray(arr1->GetEntriesFast());
+ if (arr1->GetEntriesFast()>arr0->GetEntriesFast()){
+ arr0->Expand(arr1->GetEntriesFast());
+ }
+ for (Int_t i=0;i<arr1->GetEntriesFast(); i++){
+ AliRelAlignerKalman *kalman1 = (AliRelAlignerKalman *)arr1->UncheckedAt(i);
+ AliRelAlignerKalman *kalman0 = (AliRelAlignerKalman *)arr0->UncheckedAt(i);
+ if (!kalman1) continue;
+ if (!kalman0) {arr0->AddAt(new AliRelAlignerKalman(*kalman1),i); continue;}
+ kalman0->SetRejectOutliers(kFALSE);
+ kalman0->Merge(kalman1);
+ }
+ }
+
}
return 0;
}
return tr0->GetOuterParam()->GetZ()*tr1->GetOuterParam()->GetZ()<0 && tr0->GetInnerParam()->GetZ()*tr1->GetInnerParam()->GetZ()<0 && tr0->GetOuterParam()->GetZ()*tr0->GetInnerParam()->GetZ()>0 && tr1->GetOuterParam()->GetZ()*tr1->GetInnerParam()->GetZ()>0;
}
-Bool_t AliTPCcalibTime::IsSame(AliESDtrack */*tr0*/, AliESDtrack */*tr1*/){
- // To be implemented
- return kFALSE;
+Bool_t AliTPCcalibTime::IsSame(AliESDtrack *tr0, AliESDtrack *tr1){
+ //
+ // track crossing the CE
+ // 0. minimal number of clusters
+ // 1. Same sector +-1
+ // 2. Inner and outer track param on opposite side
+ // 3. Outer and inner track parameter close each to other
+ // 3.
+ Bool_t result=kTRUE;
+ //
+ // inner and outer on opposite sides in z
+ //
+ const Int_t knclCut0 = 30;
+ const Double_t kalphaCut = 0.4;
+ //
+ // 0. minimal number of clusters
+ //
+ if (tr0->GetTPCNcls()<knclCut0) return kFALSE;
+ if (tr1->GetTPCNcls()<knclCut0) return kFALSE;
+ //
+ // 1. alpha cut - sector+-1
+ //
+ if (TMath::Abs(tr0->GetOuterParam()->GetAlpha()-tr1->GetOuterParam()->GetAlpha())>kalphaCut) return kFALSE;
+ //
+ // 2. Z crossing
+ //
+ if (tr0->GetOuterParam()->GetZ()*tr0->GetInnerParam()->GetZ()>0) result&=kFALSE;
+ if (tr1->GetOuterParam()->GetZ()*tr1->GetInnerParam()->GetZ()>0) result&=kFALSE;
+ if (result==kFALSE){
+ return result;
+ }
+ //
+ //
+ const Double_t *p0I = tr0->GetInnerParam()->GetParameter();
+ const Double_t *p1I = tr1->GetInnerParam()->GetParameter();
+ const Double_t *p0O = tr0->GetOuterParam()->GetParameter();
+ const Double_t *p1O = tr1->GetOuterParam()->GetParameter();
+ //
+ if (TMath::Abs(p0I[0]-p1I[0])>fCutMaxD) result&=kFALSE;
+ if (TMath::Abs(p0I[1]-p1I[1])>fCutMaxDz) result&=kFALSE;
+ if (TMath::Abs(p0I[2]-p1I[2])>fCutTheta) result&=kFALSE;
+ if (TMath::Abs(p0I[3]-p1I[3])>fCutTheta) result&=kFALSE;
+ if (TMath::Abs(p0O[0]-p1O[0])>fCutMaxD) result&=kFALSE;
+ if (TMath::Abs(p0O[1]-p1O[1])>fCutMaxDz) result&=kFALSE;
+ if (TMath::Abs(p0O[2]-p1O[2])>fCutTheta) result&=kFALSE;
+ if (TMath::Abs(p0O[3]-p1O[3])>fCutTheta) result&=kFALSE;
+ if (result==kTRUE){
+ result=kTRUE; // just to put break point here
+ }
+ return result;
}
-/*
-chainDrift->Draw("p0.fP[0]+p1.fP[0]","isPair");
- mean ~-0.02 ~-0.03913
- RMS ~ 0.5 ~ 0.5356 --> 3 (fCutMaxD)
-
-chainDrift->Draw("p0.fP[1]-p1.fP[1]","isPair");
- mean ~ 0.1855
- RMS ~ 4.541 -->25 (fCutMaxDz)
-
-chainDrift->Draw("p1.fAlpha-p0.fAlpha+pi","isPair");
-//chainDrift->Draw("p1.fAlpha+p0.fAlpha","isPair");
-//chainDrift->Draw("p1.fP[2]-p0.fP[2]+pi","isPair");
-//chainDrift->Draw("p1.fP[2]+p0.fP[2]","isPair");
- mean ~ 0 ~ 0.001898
- RMS ~ 0.009 ~ 0.01134 --> 0.06
-
-chainDrift->Draw("p0.fP[3]+p1.fP[3]","isPair");
- mean ~ 0.0013 ~ 0.001539
- RMS ~ 0.003 ~ 0.004644 --> 0.03 (fCutTheta)
-
-chainDrift->Draw("p1.fP[4]+p0.fP[4]>>his(100,-0.2,0.2)","isPair")
- mean ~ 0.012 ~-0.0009729
- RMS ~ 0.036 ~ 0.03773 --> 0.2
-*/
+void AliTPCcalibTime::ProcessSame(AliESDtrack* track, AliESDfriendTrack *friendTrack,AliESDEvent *event){
+ //
+ // Process TPC tracks crossing CE
+ //
+ // 0. Select only track crossing the CE
+ // 1. Cut on the track length
+ // 2. Refit the terack on A and C side separatelly
+ // 3. Fill time histograms
+ const Int_t kMinNcl=100;
+ const Int_t kMinNclS=25; // minimul number of clusters on the sides
+ if (!friendTrack->GetTPCOut()) return;
+ //
+ // 0. Select only track crossing the CE
+ //
+ if (track->GetInnerParam()->GetZ()*friendTrack->GetTPCOut()->GetZ()>0) return;
+ //
+ // 1. cut on track length
+ //
+ if (track->GetTPCNcls()<kMinNcl) return;
+ //
+ // 2. Refit track sepparatel on A and C side
+ //
+ TObject *calibObject;
+ AliTPCseed *seed = 0;
+ for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
+ if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
+ }
+ if (!seed) return;
+ //
+ AliExternalTrackParam trackIn(*track->GetInnerParam());
+ AliExternalTrackParam trackOut(*track->GetOuterParam());
+ Double_t cov[3]={0.01,0.,0.01}; //use the same errors
+ Double_t xyz[3]={0,0.,0.0};
+ Double_t bz =0;
+ Int_t nclIn=0,nclOut=0;
+ trackIn.ResetCovariance(30.);
+ trackOut.ResetCovariance(30.);
+ //
+ //2.a Refit inner
+ //
+ for (Int_t irow=0;irow<159;irow++) {
+ AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+ if (!cl) continue;
+ if (cl->GetX()<80) continue;
+ if (track->GetInnerParam()->GetZ()<0 &&(cl->GetDetector()%36)<18) break;
+ if (track->GetInnerParam()->GetZ()>0 &&(cl->GetDetector()%36)>=18) break;
+ Int_t sector = cl->GetDetector();
+ Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha();
+ if (TMath::Abs(dalpha)>0.01){
+ if (!trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break;
+ }
+ Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
+ trackIn.GetXYZ(xyz);
+ bz = AliTracker::GetBz(xyz);
+ if (!trackIn.PropagateTo(r[0],bz)) break;
+ nclIn++;
+ trackIn.Update(&r[1],cov);
+ }
+ //
+ //2.b Refit outer
+ //
+ for (Int_t irow=159;irow>0;irow--) {
+ AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+ if (!cl) continue;
+ if (cl->GetX()<80) continue;
+ if (cl->GetZ()*track->GetOuterParam()->GetZ()<0) break;
+ if (friendTrack->GetTPCOut()->GetZ()<0 &&(cl->GetDetector()%36)<18) break;
+ if (friendTrack->GetTPCOut()->GetZ()>0 &&(cl->GetDetector()%36)>=18) break;
+ Int_t sector = cl->GetDetector();
+ Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha();
+ if (TMath::Abs(dalpha)>0.01){
+ if (!trackOut.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break;
+ }
+ Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
+ trackOut.GetXYZ(xyz);
+ bz = AliTracker::GetBz(xyz);
+ if (!trackOut.PropagateTo(r[0],bz)) break;
+ nclOut++;
+ trackOut.Update(&r[1],cov);
+ }
+ trackOut.Rotate(trackIn.GetAlpha());
+ Double_t meanX = (trackIn.GetX()+trackOut.GetX())*0.5;
+ trackIn.PropagateTo(meanX,bz);
+ trackOut.PropagateTo(meanX,bz);
+ TTreeSRedirector *cstream = GetDebugStreamer();
+ if (cstream){
+ TVectorD gxyz(3);
+ trackIn.GetXYZ(gxyz.GetMatrixArray());
+ TTimeStamp tstamp(fTime);
+ Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
+ Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
+ (*cstream)<<"tpctpc"<<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
+ "ptrel0.="<<ptrelative0<<
+ "ptrel1.="<<ptrelative1<<
+ //
+ "xyz.="<<&gxyz<< // global position
+ "tIn.="<<&trackIn<< // refitterd track in
+ "tOut.="<<&trackOut<< // refitter track out
+ "nclIn="<<nclIn<< //
+ "nclOut="<<nclOut<< //
+ "\n";
+ }
+ //
+ // 3. Fill time histograms
+ // Debug stremaer expression
+ // chainTPCTPC->Draw("(tIn.fP[1]-tOut.fP[1])*sign(-tIn.fP[3]):tIn.fP[3]","min(nclIn,nclOut)>30","")
+ if (TMath::Min(nclIn,nclOut)>kMinNclS){
+ fDz = trackOut.GetZ()-trackIn.GetZ();
+ if (trackOut.GetTgl()<0) fDz*=-1.;
+ TTimeStamp tstamp(fTime);
+ Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
+ Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
+ Double_t vecDrift[4]={fTime,(ptrelative0+ptrelative1)/2.0,fDz/500.0,event->GetRunNumber()};
+ //
+ // fill histograms per trigger class and itegrated
+ //
+ THnSparse* curHist=NULL;
+ for (Int_t itype=0; itype<2; itype++){
+ TString name="MEAN_VDRIFT_CROSS_";
+ if (itype==0){
+ name+=event->GetFiredTriggerClasses();
+ name.ToUpper();
+ }else{
+ name+="ALL";
+ }
+ curHist=(THnSparseF*)fArrayDz->FindObject(name);
+ if(!curHist){
+ curHist=new THnSparseF(name,"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
+ fArrayDz->AddLast(curHist);
+ }
+ curHist->Fill(vecDrift);
+ }
+ }
+
+}
void AliTPCcalibTime::ProcessAlignITS(AliESDtrack* track, AliESDfriendTrack *friendTrack){
//
- // Process track
- // Update TPC-ITS alignment
+ // Process track - Update TPC-ITS alignment
+ // Updates:
+ // 0. Apply standartd cuts
+ // 1. Recalucluate the current statistic median/RMS
+ // 2. Apply median+-rms cut
+ // 3. Update kalman filter
+ //
+ const Int_t kMinTPC = 80; // minimal number of TPC cluster
+ const Int_t kMinITS = 3; // minimal number of ITS cluster
+ const Double_t kMinZ = 10; // maximal dz distance
+ const Double_t kMaxDy = 1.; // maximal dy distance
+ const Double_t kMaxAngle= 0.01; // maximal angular distance
+ const Double_t kSigmaCut= 5; // maximal sigma distance to median
+ const Double_t kVdErr = 0.1; // initial uncertainty of the vd correction
+ const Double_t kVdYErr = 0.05; // initial uncertainty of the vd correction
+ const Double_t kOutCut = 1.0; // outlyer cut in AliRelAlgnmentKalman
+ const Int_t kN=500; // deepnes of history
+ static Int_t kglast=0;
+ static Double_t* kgdP[4]={new Double_t[kN], new Double_t[kN], new Double_t[kN], new Double_t[kN]};
+ /*
+ 0. Standrd cuts:
+ TCut cut="abs(pTPC.fP[2]-pITS.fP[2])<0.01&&abs(pTPC.fP[3]-pITS.fP[3])<0.01&&abs(pTPC.fP[2]-pITS.fP[2])<1";
+ */
//
- const Int_t kMinTPC = 80;
- const Int_t kMinITS = 3;
- const Double_t kMinZ = 10;
- const Double_t kMaxDy = 2;
- const Double_t kMaxAngle= 0.02;
+ // 0. Apply standard cuts
//
Int_t dummycl[1000];
if (track->GetITSclusters(dummycl)<kMinITS) return; // minimal amount of clusters
if (track->GetTPCNcls()<kMinTPC) return; // minimal amount of clusters cut
- //
- if (!friendTrack->GetITSOut()) return;
+ if (!friendTrack->GetITSOut()) return;
if (!track->GetInnerParam()) return;
if (!track->GetOuterParam()) return;
// exclude crossing track
//
AliExternalTrackParam &pTPC=(AliExternalTrackParam &)(*(track->GetInnerParam()));
AliExternalTrackParam pITS(*(friendTrack->GetITSOut()));
+ pITS.Rotate(pTPC.GetAlpha());
+ pITS.PropagateTo(pTPC.GetX(),fMagF);
+ if (TMath::Abs(pITS.GetY()-pTPC.GetY()) >kMaxDy) return;
+ if (TMath::Abs(pITS.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
+ if (TMath::Abs(pITS.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
+ //
+ // 1. Update median and RMS info
+ //
+ TVectorD vecDelta(4),vecMedian(4), vecRMS(4);
+ TVectorD vecDeltaN(5);
+ Double_t sign=(pITS.GetParameter()[1]>0)? 1.:-1.;
+ vecDelta[4]=0;
+ for (Int_t i=0;i<4;i++){
+ vecDelta[i]=(pITS.GetParameter()[i]-pTPC.GetParameter()[i])*sign;
+ kgdP[i][kglast%kN]=vecDelta[i];
+ }
+ kglast=(kglast+1);
+ Int_t entries=(kglast<kN)?kglast:kN;
+ for (Int_t i=0;i<4;i++){
+ vecMedian[i] = TMath::Median(entries,kgdP[i]);
+ vecRMS[i] = TMath::RMS(entries,kgdP[i]);
+ vecDeltaN[i] = 0;
+ if (vecRMS[i]>0.){
+ vecDeltaN[i] = (vecDelta[i]-vecMedian[i])/vecRMS[i];
+ vecDeltaN[4]+= TMath::Abs(vecDeltaN[i]); //sum of abs residuals
+ }
+ }
//
+ // 2. Apply median+-rms cut
//
+ if (kglast<3) return; //median and RMS to be defined
+ if ( vecDeltaN[4]/4.>kSigmaCut) return;
+ //
+ // 3. Update alignment
//
Int_t htime = fTime/3600; //time in hours
if (fAlignITSTPC->GetEntries()<htime){
}
AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignITSTPC->At(htime);
if (!align){
+ // make Alignment object if doesn't exist
align=new AliRelAlignerKalman();
- align->SetOutRejSigma(2.);
- //align->SetRejectOutliers(kFALSE);
+ align->SetRunNumber(fRun);
+ (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
+ (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
+ align->SetOutRejSigma(kOutCut+kOutCut*kN);
+ align->SetRejectOutliers(kFALSE);
+
+ align->SetTPCvd(AliTPCcalibDB::Instance()->GetParameters()->GetDriftV()/1000000.);
align->SetMagField(fMagF);
fAlignITSTPC->AddAt(align,htime);
}
- pITS.Rotate(pTPC.GetAlpha());
- pITS.PropagateTo(pTPC.GetX(),fMagF);
- if (TMath::Abs(pITS.GetY()-pTPC.GetY())>kMaxDy) return;
- if (TMath::Abs(pITS.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
- if (TMath::Abs(pITS.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
align->AddTrackParams(&pITS,&pTPC);
align->SetTimeStamp(fTime);
- // align->SetRunNumber(fRun );
- static Int_t entry=0;
- entry++;
- // Int_t nupdates=align->GetNUpdates();
- Int_t nupdates=entry;
- align->SetOutRejSigma(1.+1./Double_t(nupdates));
+ align->SetRunNumber(fRun );
+ //
+ Int_t nupdates=align->GetNUpdates();
+ align->SetOutRejSigma(kOutCut+kOutCut*kN/Double_t(nupdates));
+ align->SetRejectOutliers(kFALSE);
TTreeSRedirector *cstream = GetDebugStreamer();
if (cstream && align->GetState() && align->GetState()->GetNrows()>2 ){
TTimeStamp tstamp(fTime);
"temp0="<<temp0<<
"temp1="<<temp1<<
"vecGoofie.="<<&vecGoofie<<
- "entry="<<entry<< // current entry
//
- "a.="<<align<< // current alignment
- "pITS.="<<&pITS<< // track param ITS
- "pTPC.="<<&pTPC<< // track param TPC
- "gpTPC.="<<&gpTPC<<
- "gdTPC.="<<&gdTPC<<
- "gpITS.="<<&gpITS<<
- "gdITS.="<<&gdITS<<
+ "nmed="<<kglast<< // number of entries to define median and RMS
+ "vMed.="<<&vecMedian<< // median of deltas
+ "vRMS.="<<&vecRMS<< // rms of deltas
+ "vDelta.="<<&vecDelta<< // delta in respect to median
+ "vDeltaN.="<<&vecDeltaN<< // normalized delta in respect to median
+ "t.="<<track<< // ful track - find proper cuts
+ "a.="<<align<< // current alignment
+ "pITS.="<<&pITS<< // track param ITS
+ "pTPC.="<<&pTPC<< // track param TPC
+ "gpTPC.="<<&gpTPC<< // global position TPC
+ "gdTPC.="<<&gdTPC<< // global direction TPC
+ "gpITS.="<<&gpITS<< // global position ITS
+ "gdITS.="<<&gdITS<< // global position ITS
"\n";
}
-
}
+
+
+
+
void AliTPCcalibTime::ProcessAlignTRD(AliESDtrack* track, AliESDfriendTrack *friendTrack){
//
- // Process track
- // Update TPC-TRD alignment
+ // Process track - Update TPC-TRD alignment
+ // Updates:
+ // 0. Apply standartd cuts
+ // 1. Recalucluate the current statistic median/RMS
+ // 2. Apply median+-rms cut
+ // 3. Update kalman filter
+ //
+ const Int_t kMinTPC = 80; // minimal number of TPC cluster
+ const Int_t kMinTRD = 50; // minimal number of TRD cluster
+ const Double_t kMinZ = 20; // maximal dz distance
+ const Double_t kMaxDy = 1.; // maximal dy distance
+ const Double_t kMaxAngle= 0.01; // maximal angular distance
+ const Double_t kSigmaCut= 5; // maximal sigma distance to median
+ const Double_t kVdErr = 0.1; // initial uncertainty of the vd correction
+ const Double_t kVdYErr = 0.05; // initial uncertainty of the vd correction
+ const Double_t kOutCut = 1.0; // outlyer cut in AliRelAlgnmentKalman
+ const Int_t kN=500; // deepnes of history
+ static Int_t kglast=0;
+ static Double_t* kgdP[4]={new Double_t[kN], new Double_t[kN], new Double_t[kN], new Double_t[kN]};
//
- const Int_t kMinTPC = 80;
- const Int_t kMinTRD = 60;
- const Double_t kMinZ = 10;
- const Double_t kMaxDy = 2;
- const Double_t kMaxAngle= 0.02;
+ // 0. Apply standard cuts
//
Int_t dummycl[1000];
if (track->GetTRDclusters(dummycl)<kMinTRD) return; // minimal amount of clusters
if (track->GetTPCNcls()<kMinTPC) return; // minimal amount of clusters cut
- //
- if (!friendTrack->GetTRDIn()) return;
+ if (!friendTrack->GetTRDIn()) return;
if (!track->GetInnerParam()) return;
if (!track->GetOuterParam()) return;
// exclude crossing track
if (track->GetOuterParam()->GetZ()*track->GetInnerParam()->GetZ()<0) return;
- if (TMath::Abs(track->GetOuterParam()->GetZ())<kMinZ) return;
+ if (TMath::Abs(track->GetInnerParam()->GetZ())<kMinZ) return;
//
AliExternalTrackParam &pTPC=(AliExternalTrackParam &)(*(track->GetOuterParam()));
AliExternalTrackParam pTRD(*(friendTrack->GetTRDIn()));
+ pTRD.Rotate(pTPC.GetAlpha());
+ pTRD.PropagateTo(pTPC.GetX(),fMagF);
+ ((Double_t*)pTRD.GetCovariance())[2]+=3.*3.; // increas sys errors
+ ((Double_t*)pTRD.GetCovariance())[9]+=0.1*0.1; // increse sys errors
+
+ if (TMath::Abs(pTRD.GetY()-pTPC.GetY()) >kMaxDy) return;
+ if (TMath::Abs(pTRD.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
+ if (TMath::Abs(pTRD.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
//
+ // 1. Update median and RMS info
//
+ TVectorD vecDelta(4),vecMedian(4), vecRMS(4);
+ TVectorD vecDeltaN(5);
+ Double_t sign=(pTRD.GetParameter()[1]>0)? 1.:-1.;
+ vecDelta[4]=0;
+ for (Int_t i=0;i<4;i++){
+ vecDelta[i]=(pTRD.GetParameter()[i]-pTPC.GetParameter()[i])*sign;
+ kgdP[i][kglast%kN]=vecDelta[i];
+ }
+ kglast=(kglast+1);
+ Int_t entries=(kglast<kN)?kglast:kN;
+ for (Int_t i=0;i<4;i++){
+ vecMedian[i] = TMath::Median(entries,kgdP[i]);
+ vecRMS[i] = TMath::RMS(entries,kgdP[i]);
+ vecDeltaN[i] = 0;
+ if (vecRMS[i]>0.){
+ vecDeltaN[i] = (vecDelta[i]-vecMedian[i])/vecRMS[i];
+ vecDeltaN[4]+= TMath::Abs(vecDeltaN[i]); //sum of abs residuals
+ }
+ }
+ //
+ // 2. Apply median+-rms cut
+ //
+ if (kglast<3) return; //median and RMS to be defined
+ if ( vecDeltaN[4]/4.>kSigmaCut) return;
+ //
+ // 3. Update alignment
//
Int_t htime = fTime/3600; //time in hours
if (fAlignTRDTPC->GetEntries()<htime){
}
AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignTRDTPC->At(htime);
if (!align){
+ // make Alignment object if doesn't exist
align=new AliRelAlignerKalman();
- align->SetOutRejSigma(2.);
- //align->SetRejectOutliers(kFALSE);
+ align->SetRunNumber(fRun);
+ (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
+ (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
+ align->SetOutRejSigma(kOutCut+kOutCut*kN);
+ align->SetRejectOutliers(kFALSE);
+ align->SetTPCvd(AliTPCcalibDB::Instance()->GetParameters()->GetDriftV()/1000000.);
align->SetMagField(fMagF);
fAlignTRDTPC->AddAt(align,htime);
}
- pTRD.Rotate(pTPC.GetAlpha());
- pTRD.PropagateTo(pTPC.GetX(),fMagF);
- ((Double_t*)pTRD.GetCovariance())[2]+=3.*3.;
- ((Double_t*)pTRD.GetCovariance())[9]+=0.1*0.1;
-
- if (TMath::Abs(pTRD.GetY()-pTPC.GetY())>kMaxDy) return;
- if (TMath::Abs(pTRD.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
- if (TMath::Abs(pTRD.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
align->AddTrackParams(&pTRD,&pTPC);
align->SetTimeStamp(fTime);
- // align->SetRunNumber(fRun );
- static Int_t entry=0;
- entry++;
- // Int_t nupdates=align->GetNUpdates();
- Int_t nupdates=entry;
- align->SetOutRejSigma(1.+1./Double_t(nupdates));
+ align->SetRunNumber(fRun );
+ //
+ Int_t nupdates=align->GetNUpdates();
+ align->SetOutRejSigma(kOutCut+kOutCut*kN/Double_t(nupdates));
+ align->SetRejectOutliers(kFALSE);
TTreeSRedirector *cstream = GetDebugStreamer();
if (cstream && align->GetState() && align->GetState()->GetNrows()>2 ){
TTimeStamp tstamp(fTime);
pTPC.GetDirection(gdTPC.GetMatrixArray());
pTRD.GetXYZ(gpTRD.GetMatrixArray());
pTRD.GetDirection(gdTRD.GetMatrixArray());
- (*cstream)<<"itstpc"<<
+ (*cstream)<<"trdtpc"<<
"run="<<fRun<< // run number
"event="<<fEvent<< // event number
"time="<<fTime<< // time stamp of event
"temp0="<<temp0<<
"temp1="<<temp1<<
"vecGoofie.="<<&vecGoofie<<
- "entry="<<entry<< // current entry
//
- "a.="<<align<< // current alignment
- "pTRD.="<<&pTRD<< // track param TRD
- "pTPC.="<<&pTPC<< // track param TPC
- "gpTPC.="<<&gpTPC<<
- "gdTPC.="<<&gdTPC<<
- "gpTRD.="<<&gpTRD<<
- "gdTRD.="<<&gdTRD<<
+ "nmed="<<kglast<< // number of entries to define median and RMS
+ "vMed.="<<&vecMedian<< // median of deltas
+ "vRMS.="<<&vecRMS<< // rms of deltas
+ "vDelta.="<<&vecDelta<< // delta in respect to median
+ "vDeltaN.="<<&vecDeltaN<< // normalized delta in respect to median
+ "t.="<<track<< // ful track - find proper cuts
+ "a.="<<align<< // current alignment
+ "pTRD.="<<&pTRD<< // track param TRD
+ "pTPC.="<<&pTPC<< // track param TPC
+ "gpTPC.="<<&gpTPC<< // global position TPC
+ "gdTPC.="<<&gdTPC<< // global direction TPC
+ "gpTRD.="<<&gpTRD<< // global position TRD
+ "gdTRD.="<<&gdTRD<< // global position TRD
"\n";
}
-
}
void AliTPCcalibTime::ProcessAlignTOF(AliESDtrack* track, AliESDfriendTrack *friendTrack){
//
- //process TOF-TPC alignment
//
- Int_t kminNcl=80;
- Float_t kMaxDy=6;
- Float_t kMaxDz=10;
- if (track->GetTPCNcls()<kminNcl) return;
- if (track->GetOuterParam()==0) return;
- if (track->GetInnerParam()==0) return;
- if (track->GetTOFsignal()<=0) return;
+ // Process track - Update TPC-TOF alignment
+ // Updates:
+ // -1. Make a TOF "track"
+ // 0. Apply standartd cuts
+ // 1. Recalucluate the current statistic median/RMS
+ // 2. Apply median+-rms cut
+ // 3. Update kalman filter
+ //
+ const Int_t kMinTPC = 80; // minimal number of TPC cluster
+ const Double_t kMinZ = 10; // maximal dz distance
+ const Double_t kMaxDy = 5.; // maximal dy distance
+ const Double_t kMaxAngle= 0.01; // maximal angular distance
+ const Double_t kSigmaCut= 5; // maximal sigma distance to median
+ const Double_t kVdErr = 0.1; // initial uncertainty of the vd correction
+ const Double_t kVdYErr = 0.05; // initial uncertainty of the vd correction
+
+ const Double_t kOutCut = 1.0; // outlyer cut in AliRelAlgnmentKalman
+ const Int_t kN=1000; // deepnes of history
+ static Int_t kglast=0;
+ static Double_t* kgdP[4]={new Double_t[kN], new Double_t[kN], new Double_t[kN], new Double_t[kN]};
//
- AliExternalTrackParam *paramOut = new AliExternalTrackParam(*(track->GetOuterParam()));
- AliExternalTrackParam *param=0;
+ // -1. Make a TOF track-
+ // Clusters are not in friends - use alingment points
+ //
+ if (track->GetTOFsignal()<=0) return;
+ if (!friendTrack->GetTPCOut()) return;
+ if (!track->GetInnerParam()) return;
+ if (!track->GetOuterParam()) return;
const AliTrackPointArray *points=friendTrack->GetTrackPointArray();
if (!points) return;
+ AliExternalTrackParam pTPC(*(track->GetOuterParam()));
+ AliExternalTrackParam pTOF(pTPC);
+ Double_t mass = TDatabasePDG::Instance()->GetParticle("mu+")->Mass();
Int_t npoints = points->GetNPoints();
AliTrackPoint point;
- //Double_t alpha=
- Double_t mass = TDatabasePDG::Instance()->GetParticle("mu+")->Mass();
- TTreeSRedirector * cstream = GetDebugStreamer();
- //
- //
+ Int_t naccept=0;
//
for (Int_t ipoint=0;ipoint<npoints;ipoint++){
- //
points->GetPoint(point,ipoint);
Float_t xyz[3];
point.GetXYZ(xyz);
Double_t r=TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
- if (r<300) continue;
+ if (r<350) continue;
if (r>400) continue;
- param=paramOut;
- if (!param) continue;
- AliTracker::PropagateTrackToBxByBz(param,r,mass,2.,kTRUE);
- AliTracker::PropagateTrackToBxByBz(param,r,mass,0.1,kTRUE);
- AliTrackPoint lpoint = point.Rotate(param->GetAlpha());
- param->PropagateTo(lpoint.GetX(),fMagF);
- //
- //
- // this is ugly - we need AliCluster constructor
- //
- AliExternalTrackParam &pTPC=*param;
- AliExternalTrackParam pTOF(*param);
+ AliTracker::PropagateTrackToBxByBz(&pTPC,r,mass,2.,kTRUE);
+ AliTracker::PropagateTrackToBxByBz(&pTPC,r,mass,0.1,kTRUE);
+ AliTrackPoint lpoint = point.Rotate(pTPC.GetAlpha());
+ pTPC.PropagateTo(lpoint.GetX(),fMagF);
+ pTOF=pTPC;
((Double_t*)pTOF.GetParameter())[0] =lpoint.GetY();
((Double_t*)pTOF.GetParameter())[1] =lpoint.GetZ();
- pTOF.ResetCovariance(20);
- ((Double_t*)pTOF.GetCovariance())[0]+=3.*3.;
- ((Double_t*)pTOF.GetCovariance())[2]+=3.*3.;
+ ((Double_t*)pTOF.GetCovariance())[0]+=3.*3./12.;
+ ((Double_t*)pTOF.GetCovariance())[2]+=3.*3./12.;
((Double_t*)pTOF.GetCovariance())[5]+=0.1*0.1;
((Double_t*)pTOF.GetCovariance())[9]+=0.1*0.1;
- if (TMath::Abs(pTOF.GetY()-pTPC.GetY())>kMaxDy) continue;
- if (TMath::Abs(pTOF.GetZ()-pTPC.GetZ())>kMaxDz) continue;
- //
- Int_t htime = fTime/3600; //time in hours
-
- if (fAlignTOFTPC->GetEntries()<htime){
- fAlignTOFTPC->Expand(htime*2+20);
- }
- AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignTOFTPC->At(htime);
- if (!align){
- align=new AliRelAlignerKalman();
- align->SetOutRejSigma(2.);
- //align->SetRejectOutliers(kFALSE);
- align->SetMagField(fMagF);
- fAlignTOFTPC->AddAt(align,htime);
+ naccept++;
+ }
+ if (naccept==0) return; // no tof match clusters
+ //
+ // 0. Apply standard cuts
+ //
+ if (track->GetTPCNcls()<kMinTPC) return; // minimal amount of clusters cut
+ // exclude crossing track
+ if (track->GetOuterParam()->GetZ()*track->GetInnerParam()->GetZ()<0) return;
+ //
+ if (TMath::Abs(pTOF.GetY()-pTPC.GetY()) >kMaxDy) return;
+ if (TMath::Abs(pTOF.GetSnp()-pTPC.GetSnp())>kMaxAngle) return;
+ if (TMath::Abs(pTOF.GetTgl()-pTPC.GetTgl())>kMaxAngle) return;
+ //
+ // 1. Update median and RMS info
+ //
+ TVectorD vecDelta(4),vecMedian(4), vecRMS(4);
+ TVectorD vecDeltaN(5);
+ Double_t sign=(pTOF.GetParameter()[1]>0)? 1.:-1.;
+ vecDelta[4]=0;
+ for (Int_t i=0;i<4;i++){
+ vecDelta[i]=(pTOF.GetParameter()[i]-pTPC.GetParameter()[i])*sign;
+ kgdP[i][kglast%kN]=vecDelta[i];
+ }
+ kglast=(kglast+1);
+ Int_t entries=(kglast<kN)?kglast:kN;
+ Bool_t isOK=kTRUE;
+ for (Int_t i=0;i<4;i++){
+ vecMedian[i] = TMath::Median(entries,kgdP[i]);
+ vecRMS[i] = TMath::RMS(entries,kgdP[i]);
+ vecDeltaN[i] = 0;
+ if (vecRMS[i]>0.){
+ vecDeltaN[i] = (vecDelta[i]-vecMedian[i])/(vecRMS[i]+1.);
+ vecDeltaN[4]+= TMath::Abs(vecDeltaN[i]); //sum of abs residuals
+ if (TMath::Abs(vecDeltaN[i])>kSigmaCut) isOK=kFALSE;
}
- pTOF.Rotate(pTPC.GetAlpha());
- pTOF.PropagateTo(pTPC.GetX(),fMagF);
- align->AddTrackParams(&pTOF,&pTPC);
- align->SetTimeStamp(fTime);
- static Int_t entry=0;
- entry++;
- // Int_t nupdates=align->GetNUpdates();
- Int_t nupdates=entry;
- align->SetOutRejSigma(1.+1./Double_t(nupdates));
-
- //
- //
- if (cstream) {
- (*cstream) << "pointMatch" <<
- "run="<<fRun<< // run number
- "event="<<fEvent<< // event number
- "time="<<fTime<< // time stamp of event
- "trigger="<<fTrigger<< // trigger
- "mag="<<fMagF<< // magnetic field
- //
- "a.="<<align<< // current alignment
- "p.="<<&point<<
- "lp.="<<&lpoint<<
- "pTPC.="<<&pTPC<<
- "pTOF.="<<&pTOF<<
- "\n";
+ }
+ //
+ // 2. Apply median+-rms cut
+ //
+ if (kglast<10) return; //median and RMS to be defined
+ if (!isOK) return;
+ //
+ // 3. Update alignment
+ //
+ Int_t htime = fTime/3600; //time in hours
+ if (fAlignTOFTPC->GetEntries()<htime){
+ fAlignTOFTPC->Expand(htime*2+20);
+ }
+ AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignTOFTPC->At(htime);
+ if (!align){
+ // make Alignment object if doesn't exist
+ align=new AliRelAlignerKalman();
+ align->SetRunNumber(fRun);
+ (*align->GetStateCov())(6,6)=kVdErr*kVdErr;
+ (*align->GetStateCov())(8,8)=kVdYErr*kVdYErr;
+ align->SetOutRejSigma(kOutCut+kOutCut*kN);
+ align->SetRejectOutliers(kFALSE);
+ align->SetTPCvd(AliTPCcalibDB::Instance()->GetParameters()->GetDriftV()/1000000.);
+ align->SetMagField(fMagF);
+ fAlignTOFTPC->AddAt(align,htime);
+ }
+ align->AddTrackParams(&pTOF,&pTPC);
+ align->SetTimeStamp(fTime);
+ align->SetRunNumber(fRun );
+ //
+ Int_t nupdates=align->GetNUpdates();
+ align->SetOutRejSigma(kOutCut+kOutCut*kN/Double_t(nupdates));
+ align->SetRejectOutliers(kFALSE);
+ TTreeSRedirector *cstream = GetDebugStreamer();
+ if (cstream && align->GetState() && align->GetState()->GetNrows()>2 ){
+ TTimeStamp tstamp(fTime);
+ Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
+ Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
+ Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
+ Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
+ Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
+ Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
+ TVectorD vecGoofie(20);
+ AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
+ if (goofieArray){
+ for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
+ AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
+ if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
+ }
}
-
-
-
+ TVectorD gpTPC(3), gdTPC(3);
+ TVectorD gpTOF(3), gdTOF(3);
+ pTPC.GetXYZ(gpTPC.GetMatrixArray());
+ pTPC.GetDirection(gdTPC.GetMatrixArray());
+ pTOF.GetXYZ(gpTOF.GetMatrixArray());
+ pTOF.GetDirection(gdTOF.GetMatrixArray());
+ (*cstream)<<"toftpc"<<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
+ // Environment values
+ "press0="<<valuePressure0<<
+ "press1="<<valuePressure1<<
+ "pt0="<<ptrelative0<<
+ "pt1="<<ptrelative1<<
+ "temp0="<<temp0<<
+ "temp1="<<temp1<<
+ "vecGoofie.="<<&vecGoofie<<
+ //
+ "nmed="<<kglast<< // number of entries to define median and RMS
+ "vMed.="<<&vecMedian<< // median of deltas
+ "vRMS.="<<&vecRMS<< // rms of deltas
+ "vDelta.="<<&vecDelta<< // delta in respect to median
+ "vDeltaN.="<<&vecDeltaN<< // normalized delta in respect to median
+ "t.="<<track<< // ful track - find proper cuts
+ "a.="<<align<< // current alignment
+ "pTOF.="<<&pTOF<< // track param TOF
+ "pTPC.="<<&pTPC<< // track param TPC
+ "gpTPC.="<<&gpTPC<< // global position TPC
+ "gdTPC.="<<&gdTPC<< // global direction TPC
+ "gpTOF.="<<&gpTOF<< // global position TOF
+ "gdTOF.="<<&gdTOF<< // global position TOF
+ "\n";
}
- delete paramOut;
}
+
+
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff