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","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 "Riostream.h"
-#include "TChain.h"
-#include "TTree.h"
+#include "TDatabasePDG.h"
+#include "TGraphErrors.h"
#include "TH1F.h"
-#include "TH2F.h"
-#include "TH3F.h"
#include "THnSparse.h"
#include "TList.h"
#include "TMath.h"
-#include "TCanvas.h"
-#include "TFile.h"
-#include "TF1.h"
+#include "TTimeStamp.h"
+#include "TTree.h"
#include "TVectorD.h"
-#include "TProfile.h"
-#include "TGraphErrors.h"
-#include "TCanvas.h"
-#include "AliTPCclusterMI.h"
-#include "AliTPCseed.h"
-#include "AliESDVertex.h"
+//#include "TChain.h"
+//#include "TFile.h"
+
+#include "AliDCSSensor.h"
+#include "AliDCSSensorArray.h"
#include "AliESDEvent.h"
-#include "AliESDfriend.h"
#include "AliESDInputHandler.h"
-#include "AliAnalysisManager.h"
-
-#include "AliTracker.h"
-#include "AliMagF.h"
-#include "AliTPCCalROC.h"
-
+#include "AliESDVertex.h"
+#include "AliESDfriend.h"
#include "AliLog.h"
-
-#include "AliTPCcalibTime.h"
#include "AliRelAlignerKalman.h"
-
-#include "TTreeStream.h"
+#include "AliTPCCalROC.h"
+#include "AliTPCParam.h"
#include "AliTPCTracklet.h"
-#include "TTimeStamp.h"
#include "AliTPCcalibDB.h"
#include "AliTPCcalibLaser.h"
-#include "AliDCSSensorArray.h"
-#include "AliDCSSensor.h"
-
-#include "TDatabasePDG.h"
+#include "AliTPCcalibTime.h"
+#include "AliTPCclusterMI.h"
+#include "AliTPCseed.h"
#include "AliTrackPointArray.h"
+#include "AliTracker.h"
ClassImp(AliTPCcalibTime)
fCutMaxDz(25), // maximal distance in rfi ditection
fCutTheta(0.03), // maximal distan theta
fCutMinDir(-0.99), // direction vector products
- fCutTracks(10),
+ fCutTracks(100),
fArrayDz(0), //NEW! Tmap of V drifts for different triggers
fAlignITSTPC(0), //alignemnt array ITS TPC match
fAlignTRDTPC(0), //alignemnt array TRD TPC match
// fXminVdrift(fTimeStart,fPtStart,fVdriftStart),
// fXmaxVdrift(fTimeEnd,fPtEnd,fVdriftEnd)
{
+ //
+ // default constructor
+ //
AliInfo("Default Constructor");
for (Int_t i=0;i<3;i++) {
fHistVdriftLaserA[i]=0;
fCutMaxDz(40), // maximal distance in rfi ditection
fCutTheta(5*0.004644),// maximal distan theta
fCutMinDir(-0.99), // direction vector products
- fCutTracks(10),
+ fCutTracks(100),
fArrayDz(0), //Tmap of V drifts for different triggers
fAlignITSTPC(0), //alignemnt array ITS TPC match
fAlignTRDTPC(0), //alignemnt array TRD TPC match
fRunStart(0),
fRunEnd(0)
{
+ //
+ // Non deafaul constructor - to be used in the Calibration setups
+ //
+
SetName(name);
SetTitle(title);
for (Int_t i=0;i<3;i++) {
fVdriftBins = 500;
fVdriftStart= -0.1;
fVdriftEnd = 0.1;
- fRunBins = 100001;
+ fRunBins = 1000001;
fRunStart = -1.5;
- fRunEnd = 99999.5;
+ fRunEnd = 999999.5;
Int_t binsVdriftLaser[4] = {fTimeBins , fPtBins , fVdriftBins*20, fRunBins };
Double_t xminVdriftLaser[4] = {fTimeStart, fPtStart, fVdriftStart , fRunStart};
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]);
AliTPCcalibTime::~AliTPCcalibTime(){
//
- // Destructor
+ // Virtual Destructor
//
for(Int_t i=0;i<3;i++){
if(fHistVdriftLaserA[i]){
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*/){
+Bool_t AliTPCcalibTime::IsLaser(const AliESDEvent *const /*event*/){
+ //
+ // Indicator is laser event not yet implemented - to be done using trigger info or event specie
+ //
return kTRUE; //More accurate creteria to be added
}
-Bool_t AliTPCcalibTime::IsCosmics(AliESDEvent */*event*/){
+Bool_t AliTPCcalibTime::IsCosmics(const AliESDEvent *const /*event*/){
+ //
+ // Indicator is cosmic event not yet implemented - to be done using trigger info or event specie
+ //
+
return kTRUE; //More accurate creteria to be added
}
-Bool_t AliTPCcalibTime::IsBeam(AliESDEvent */*event*/){
+Bool_t AliTPCcalibTime::IsBeam(const AliESDEvent *const /*event*/){
+ //
+ // Indicator is physic event not yet implemented - to be done using trigger info or event specie
+ //
+
return kTRUE; //More accurate creteria to be added
}
void AliTPCcalibTime::ResetCurrent(){
fDz=0; //Reset current dz
}
+
+
+
void AliTPCcalibTime::Process(AliESDEvent *event){
+ //
+ // main function to make calibration
+ //
if(!event) return;
if (event->GetNumberOfTracks()<2) return;
ResetCurrent();
Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
+ Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
TVectorD vecGoofie(20);
AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
if (goofieArray){
"temp0="<<temp0<<
"temp1="<<temp1<<
"vecGoofie.="<<&vecGoofie<<
+ "vdcorr="<<vdcorr<<
//laser
"rejectA="<<isReject[0]<<
"rejectC="<<isReject[1]<<
// curHist->Fill(vecDrift);
}
-void AliTPCcalibTime::ProcessCosmic(AliESDEvent *event){
+void AliTPCcalibTime::ProcessCosmic(const AliESDEvent *const event){
+ //
+ // process Cosmic event - track matching A side C side
+ //
if (!event) {
Printf("ERROR: ESD not available");
return;
// Track0 is choosen in upper TPC part
// Track1 is choosen in lower TPC part
//
+ const Int_t kMinClustersCross =30;
+ const Int_t kMinClusters =80;
Int_t ntracks=event->GetNumberOfTracks();
if (ntracks==0) return;
if (ntracks > fCutTracks) return;
- if (GetDebugLevel()>1) printf("Hallo world: Im here\n");
- AliESDfriend *
ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
+ if (GetDebugLevel()>20) printf("Hallo world: Im here\n");
+ AliESDfriend *
esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
TObjArray tpcSeeds(ntracks);
Double_t vtxx[3]={0,0,0};
//
// track loop
//
+ TArrayI clusterSideA(ntracks);
+ TArrayI clusterSideC(ntracks);
for (Int_t i=0;i<ntracks;++i) {
+ clusterSideA[i]=0;
+ clusterSideC[i]=0;
AliESDtrack *track = event->GetTrack(i);
const AliExternalTrackParam * trackIn = track->GetInnerParam();
if (!trackIn) continue;
if (!trackOut) continue;
- AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
- if (friendTrack) ProcessAlignITS(track,friendTrack);
+ AliESDfriendTrack *friendTrack = esdFriend->GetTrack(i);
+ if (!friendTrack) continue;
+ if (friendTrack) ProcessSame(track,friendTrack,event);
+ if (friendTrack) ProcessAlignITS(track,friendTrack,event,esdFriend);
+ if (friendTrack) ProcessAlignTRD(track,friendTrack);
if (friendTrack) ProcessAlignTOF(track,friendTrack);
TObject *calibObject;
AliTPCseed *seed = 0;
for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
- if (seed) tpcSeeds.AddAt(seed,i);
+ if (seed) {
+ tpcSeeds.AddAt(seed,i);
+ Int_t nA=0, nC=0;
+ for (Int_t irow=159;irow>0;irow--) {
+ AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+ if (!cl) continue;
+ if ((cl->GetDetector()%36)<18) nA++;
+ if ((cl->GetDetector()%36)>=18) nC++;
+ }
+ clusterSideA[i]=nA;
+ clusterSideC[i]=nC;
+ }
}
if (ntracks<2) return;
//
//track 1 lower part
if (!track1) continue;
if (!track1->GetOuterParam()) continue;
- if (track1->GetOuterParam()->GetAlpha()>0) continue;
+ if (track0->GetTPCNcls()+ track1->GetTPCNcls()< kMinClusters) continue;
+ Int_t nAC = TMath::Max( TMath::Min(clusterSideA[i], clusterSideC[j]),
+ TMath::Min(clusterSideC[i], clusterSideA[j]));
+ if (nAC<kMinClustersCross) continue;
+ Int_t nA0=clusterSideA[i];
+ Int_t nC0=clusterSideC[i];
+ Int_t nA1=clusterSideA[j];
+ Int_t nC1=clusterSideC[j];
+ // 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;
+ Double_t sign=(nA0>nA1)? 1:-1;
+ fDz*=sign;
TTimeStamp tstamp(fTime);
Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
"tr1.="<<track1<<
"p0.="<<¶m0<<
"p1.="<<¶m1<<
+ "nAC="<<nAC<<
+ "nA0="<<nA0<<
+ "nA1="<<nA1<<
+ "nC0="<<nC0<<
+ "nC1="<<nC1<<
"isPair="<<isPair<<
"isCross="<<isCross<<
"isSame="<<isSame<<
Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
+ Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
TVectorD vecGoofie(20);
AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
if (goofieArray){
"temp0="<<temp0<<
"temp1="<<temp1<<
"vecGoofie.=<<"<<&vecGoofie<<
+ "vdcorr="<<vdcorr<<
//
// accumulated values
//
"\n";
}
}
- printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
+ if (GetDebugLevel()>20) printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
}
-void AliTPCcalibTime::ProcessBeam(AliESDEvent */*event*/){
+void AliTPCcalibTime::ProcessBeam(const AliESDEvent *const /*event*/){
+ //
+ // Not special treatment yet - the same for cosmic and physic event
+ //
}
-void AliTPCcalibTime::Analyze(){}
+void AliTPCcalibTime::Analyze(){
+ //
+ // Special macro to analyze result of calibration and extract calibration entries
+ // Not yet ported to the Analyze function yet
+ //
+}
-THnSparse* AliTPCcalibTime::GetHistoDrift(const char* name){
+THnSparse* AliTPCcalibTime::GetHistoDrift(const char* name) const
+{
+ //
+ // Get histogram for given trigger mask
+ //
TIterator* iterator = fArrayDz->MakeIterator();
iterator->Reset();
TString newName=name;
return newHist;
}
-TObjArray* AliTPCcalibTime::GetHistoDrift(){
+TObjArray* AliTPCcalibTime::GetHistoDrift() const
+{
+ //
+ // return array of histograms
+ //
return fArrayDz;
}
TGraphErrors* AliTPCcalibTime::GetGraphDrift(const char* name){
+ //
+ // Make a drift velocity (delta Z) graph
+ //
THnSparse* histoDrift=GetHistoDrift(name);
TGraphErrors* graphDrift=NULL;
if(histoDrift){
}
TObjArray* AliTPCcalibTime::GetGraphDrift(){
+ //
+ // make a array of drift graphs
+ //
TObjArray* arrayGraphDrift=new TObjArray();
TIterator* iterator=fArrayDz->MakeIterator();
iterator->Reset();
}
AliSplineFit* AliTPCcalibTime::GetFitDrift(const char* name){
+ //
+ // Make a fit AliSplinefit of drift velocity
+ //
TGraph* graphDrift=GetGraphDrift(name);
AliSplineFit* fitDrift=NULL;
if(graphDrift && graphDrift->GetN()){
// return arrayFitDrift;
//}
-Long64_t AliTPCcalibTime::Merge(TCollection *li) {
+Long64_t AliTPCcalibTime::Merge(TCollection *const li) {
+ //
+ // Object specific merging procedure
+ //
TIterator* iter = li->MakeIterator();
AliTPCcalibTime* cal = 0;
// 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 kTRUE;
}
-Bool_t AliTPCcalibTime::IsCross(AliESDtrack *tr0, AliESDtrack *tr1){
- 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::IsCross(AliESDtrack *const tr0, AliESDtrack *const tr1){
+ //
+ // check if the cosmic pair of tracks crossed A/C side
+ //
+ Bool_t result= tr0->GetOuterParam()->GetZ()*tr1->GetOuterParam()->GetZ()<0;
+ if (result==kFALSE) return result;
+ result=kTRUE;
+ return result;
}
-Bool_t AliTPCcalibTime::IsSame(AliESDtrack */*tr0*/, AliESDtrack */*tr1*/){
- // To be implemented
- return kFALSE;
+Bool_t AliTPCcalibTime::IsSame(AliESDtrack *const tr0, AliESDtrack *const 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 *const track, AliESDfriendTrack *const friendTrack, const AliESDEvent *const 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);
+ Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,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<<
+ "vdcorr="<<vdcorr<< // drift correction applied
+ //
+ "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){
+void AliTPCcalibTime::ProcessAlignITS(AliESDtrack *const track, AliESDfriendTrack *const friendTrack, const AliESDEvent *const event, AliESDfriend *const esdFriend){
//
- // 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 = 2.; // maximal dy distance
+ const Double_t kMaxAngle= 0.015; // 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 Double_t kMinPt = 0.3; // minimal pt
+ 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 (track->GetITSclusters(dummycl)<kMinITS) return; // minimal amount of clusters
+ if (!track->IsOn(AliESDtrack::kTPCrefit)) return;
+ if (!friendTrack->GetITSOut()) return;
if (!track->GetInnerParam()) return;
if (!track->GetOuterParam()) return;
+ if (track->GetInnerParam()->Pt()<kMinPt) return;
// exclude crossing track
if (track->GetOuterParam()->GetZ()*track->GetInnerParam()->GetZ()<0) return;
if (TMath::Abs(track->GetInnerParam()->GetZ())<kMinZ) return;
+ if (track->GetInnerParam()->GetX()>90) return;
//
AliExternalTrackParam &pTPC=(AliExternalTrackParam &)(*(track->GetInnerParam()));
- AliExternalTrackParam pITS(*(friendTrack->GetITSOut()));
+ AliExternalTrackParam pITS(*(friendTrack->GetITSOut())); // ITS standalone if possible
+ AliExternalTrackParam pITS2(*(friendTrack->GetITSOut())); //TPC-ITS track
+ pITS2.Rotate(pTPC.GetAlpha());
+ pITS2.PropagateTo(pTPC.GetX(),fMagF);
+ AliESDfriendTrack *itsfriendTrack=0;
+ //
+ // try to find standalone ITS track corresponing to the TPC if possible
+ //
+ Bool_t hasAlone=kFALSE;
+ Int_t ntracks=event->GetNumberOfTracks();
+ for (Int_t i=0; i<ntracks; i++){
+ AliESDtrack *trackS = event->GetTrack(i);
+ if (trackS->GetTPCNcls()>0) continue; //continue if has TPC info
+ itsfriendTrack = esdFriend->GetTrack(i);
+ if (!itsfriendTrack) continue;
+ if (!itsfriendTrack->GetITSOut()) continue;
+ if (TMath::Abs(pITS2.GetTgl()-itsfriendTrack->GetITSOut()->GetTgl())> kMaxAngle) continue;
+ pITS=(*(itsfriendTrack->GetITSOut()));
+ //
+ pITS.Rotate(pTPC.GetAlpha());
+ pITS.PropagateTo(pTPC.GetX(),fMagF);
+ if (TMath::Abs(pITS2.GetY()-pITS.GetY())> kMaxDy) continue;
+ hasAlone=kTRUE;
+ }
+ if (!hasAlone) pITS=pITS2;
+ //
+ 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=-1;
- 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());
pITS.GetXYZ(gpITS.GetMatrixArray());
pITS.GetDirection(gdITS.GetMatrixArray());
+ Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
(*cstream)<<"itstpc"<<
"run="<<fRun<< // run number
"event="<<fEvent<< // event number
"temp0="<<temp0<<
"temp1="<<temp1<<
"vecGoofie.="<<&vecGoofie<<
- "entry="<<entry<< // current entry
+ "vdcorr="<<vdcorr<< // drift correction applied
//
- "a.="<<align<< // current alignment
- "pITS.="<<&pITS<< // track param ITS
- "pTPC.="<<&pTPC<< // track param TPC
- "gpTPC.="<<&gpTPC<<
- "gdTPC.="<<&gdTPC<<
- "gpITS.="<<&gpITS<<
- "gdITS.="<<&gdITS<<
+ "hasAlone="<<hasAlone<< // has ITS standalone ?
+ "track.="<<track<< // track info
+ "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
+ "pITS2.="<<&pITS2<< // track param ITS+TPC
+ "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::ProcessAlignTOF(AliESDtrack* track, AliESDfriendTrack *friendTrack){
+
+
+void AliTPCcalibTime::ProcessAlignTRD(AliESDtrack *const track, AliESDfriendTrack *const friendTrack){
//
- //process TOF-TPC 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
//
- 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;
+ 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 = 2.; // maximal dy distance
+ const Double_t kMaxAngle= 0.015; // 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. 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 (!track->GetInnerParam()) return;
+ if (!track->GetOuterParam()) return;
+ // exclude crossing track
+ if (track->GetOuterParam()->GetZ()*track->GetInnerParam()->GetZ()<0) 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
+ }
+ }
//
- AliExternalTrackParam *paramOut = new AliExternalTrackParam(*(track->GetOuterParam()));
- AliExternalTrackParam *param=0;
+ // 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){
+ fAlignTRDTPC->Expand(htime*2+20);
+ }
+ AliRelAlignerKalman* align = (AliRelAlignerKalman*)fAlignTRDTPC->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);
+ fAlignTRDTPC->AddAt(align,htime);
+ }
+ align->AddTrackParams(&pTRD,&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 gpTRD(3), gdTRD(3);
+ pTPC.GetXYZ(gpTPC.GetMatrixArray());
+ pTPC.GetDirection(gdTPC.GetMatrixArray());
+ pTRD.GetXYZ(gpTRD.GetMatrixArray());
+ pTRD.GetDirection(gdTRD.GetMatrixArray());
+ Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
+ (*cstream)<<"trdtpc"<<
+ "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<<
+ "vdcorr="<<vdcorr<< // drift correction applied
+ //
+ "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 *const track, AliESDfriendTrack *const friendTrack){
+ //
+ //
+ // 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.015; // 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]};
+ //
+ // -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.;
- 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);
+ ((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;
+ 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=-1;
- 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());
+ Double_t vdcorr = AliTPCcalibDB::Instance()->GetVDriftCorrectionTime(tstamp,fRun,0,1);
+ (*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<<
+ "vdcorr="<<vdcorr<< // drift correction applied
+ //
+ "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 / blobdiff