**************************************************************************/
/*
- Comments to be written here:
+ Comments to be written here:
1. What do we calibrate.
2. How to interpret results
3. Simple example
#include "AliLog.h"
#include "AliTPCcalibCosmic.h"
-
+#include "AliExternalComparison.h"
#include "TTreeStream.h"
#include "AliTPCTracklet.h"
fDeDxMIP(0),
fMIPvalue(1),
fCutMaxD(5), // maximal distance in rfi ditection
+ fCutMaxDz(40), // maximal distance in z ditection
fCutTheta(0.03), // maximal distan theta
fCutMinDir(-0.99) // direction vector products
{
fDeDx(0),
fDeDxMIP(0),
fMIPvalue(1),
- fCutMaxD(5), // maximal distance in rfi ditection
+ fCutMaxD(5), // maximal distance in rfi ditection
+ fCutMaxDz(40), // maximal distance in z ditection
fCutTheta(0.03), // maximal distan theta
fCutMinDir(-0.99) // direction vector products
{
SetName(name);
SetTitle(title);
- AliMagFMaps * field = new AliMagFMaps("dummy1", "dummy2",0,5,0);
- AliTracker::SetFieldMap(field, kTRUE);
fHistNTracks = new TH1F("ntracks","Number of Tracks per Event; number of tracks per event; number of tracks",501,-0.5,500.5);
fClusters = new TH1F("signal","Number of Clusters per track; number of clusters per track n_{cl}; counts",160,0,160);
return;
}
+
FindPairs(event); // nearly everything takes place in find pairs...
if (GetDebugLevel()>20) printf("Hallo world: Im here and processing an event\n");
const AliExternalTrackParam * trackOut = track->GetOuterParam();
if (!trackIn) continue;
if (!trackOut) continue;
-
+ if (ntracks>4 && TMath::Abs(trackIn->GetTgl())<0.0015) continue; // filter laser
+
+
AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
TObject *calibObject;
AliTPCseed *seed = 0;
if (!track0) continue;
if (!track0->GetOuterParam()) continue;
if (track0->GetOuterParam()->GetAlpha()<0) continue;
- Double_t d1[3];
- track0->GetDirection(d1);
+ Double_t dir0[3];
+ track0->GetDirection(dir0);
for (Int_t j=0;j<ntracks;++j) {
if (i==j) continue;
AliESDtrack *track1 = event->GetTrack(j);
if (!track1->GetOuterParam()) continue;
if (track1->GetOuterParam()->GetAlpha()>0) continue;
//
- Double_t d2[3];
- track1->GetDirection(d2);
+ Double_t dir1[3];
+ track1->GetDirection(dir1);
AliTPCseed * seed0 = (AliTPCseed*) tpcSeeds.At(i);
AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j);
Float_t dedx0O = seed0->CookdEdxNorm(0.05,0.55,0,64,159,fGainMap);
Float_t dedx1O = seed1->CookdEdxNorm(0.05,0.55,0,64,159,fGainMap);
//
- Float_t dir = (d1[0]*d2[0] + d1[1]*d2[1] + d1[2]*d2[2]);
+ Float_t dir = (dir0[0]*dir1[0] + dir0[1]*dir1[1] + dir0[2]*dir1[2]);
Float_t d0 = track0->GetLinearD(0,0);
Float_t d1 = track1->GetLinearD(0,0);
//
//
param0.GetDZ(0,0,0,bz,dvertex0);
param1.GetDZ(0,0,0,bz,dvertex1);
+ if (TMath::Abs(param0.GetZ()-param1.GetZ())>fCutMaxDz) continue;
//
Double_t xyz0[3];//,pxyz0[3];
Double_t xyz1[3];//,pxyz1[3];
//
if (fStreamLevel>0){
TTreeSRedirector * cstream = GetDebugStreamer();
- printf("My stream=%p\n",(void*)cstream);
+ //printf("My stream=%p\n",(void*)cstream);
+ AliExternalTrackParam *ip0 = (AliExternalTrackParam *)track0->GetInnerParam();
+ AliExternalTrackParam *ip1 = (AliExternalTrackParam *)track1->GetInnerParam();
+ AliExternalTrackParam *op0 = (AliExternalTrackParam *)track0->GetOuterParam();
+ AliExternalTrackParam *op1 = (AliExternalTrackParam *)track1->GetOuterParam();
+ Bool_t isCrossI = ip0->GetZ()*ip1->GetZ()<0;
+ Bool_t isCrossO = op0->GetZ()*op1->GetZ()<0;
+ Double_t alpha0 = TMath::ATan2(dir0[1],dir0[0]);
+ Double_t alpha1 = TMath::ATan2(dir1[1],dir1[0]);
if (cstream) {
(*cstream) << "Track0" <<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
"dir="<<dir<< // direction
"OK="<<isPair<< // will be accepted
"b0="<<b0<< // propagate status
"b1="<<b1<< // propagate status
+ "crossI="<<isCrossI<< // cross inner
+ "crossO="<<isCrossO<< // cross outer
+ //
"Orig0.=" << track0 << // original track 0
"Orig1.=" << track1 << // original track 1
"Tr0.="<<¶m0<< // track propagated to the DCA 0,0
"Tr1.="<<¶m1<< // track propagated to the DCA 0,0
+ "Ip0.="<<ip0<< // inner param - upper
+ "Ip1.="<<ip1<< // inner param - lower
+ "Op0.="<<op0<< // outer param - upper
+ "Op1.="<<op1<< // outer param - lower
+ //
"v00="<<dvertex0[0]<< // distance using kalman
"v01="<<dvertex0[1]<< //
"v10="<<dvertex1[0]<< //
"d0="<<d0<< // linear distance to 0,0
"d1="<<d1<< // linear distance to 0,0
//
- "x00="<<xyz0[0]<<
+ //
+ //
+ "x00="<<xyz0[0]<< // global position close to vertex
"x01="<<xyz0[1]<<
"x02="<<xyz0[2]<<
//
- "x10="<<xyz1[0]<<
+ "x10="<<xyz1[0]<< // global position close to vertex
"x11="<<xyz1[1]<<
"x12="<<xyz1[2]<<
//
+ "alpha0="<<alpha0<<
+ "alpha1="<<alpha1<<
+ "dir00="<<dir0[0]<< // direction upper
+ "dir01="<<dir0[1]<<
+ "dir02="<<dir0[2]<<
+ //
+ "dir10="<<dir1[0]<< // direction lower
+ "dir11="<<dir1[1]<<
+ "dir12="<<dir1[2]<<
+ //
+ //
"Seed0.=" << seed0 << // original seed 0
"Seed1.=" << seed1 << // original seed 1
//
while ((cal = (AliTPCcalibCosmic*)iter->Next())) {
if (!cal->InheritsFrom(AliTPCcalibCosmic::Class())) {
- Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
+ //Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
return -1;
}
const Double_t *p0 = tr0->GetParameter();
const Double_t *p1 = tr1->GetParameter();
if (TMath::Abs(p0[3]+p1[3])>fCutTheta) return kFALSE;
+ if (TMath::Abs(p0[1]-p1[1])>fCutMaxDz) return kFALSE;
if (TMath::Abs(p0[0]+p1[0])>fCutMaxD) return kFALSE;
+
Double_t d0[3], d1[3];
tr0->GetDirection(d0);
tr1->GetDirection(d1);
-void AliTPCcalibCosmic::CalculateBetheParams(TH2F *hist, Double_t * initialParam) {
-
+void AliTPCcalibCosmic::CalculateBetheParams(TH2F */*hist*/, Double_t * /*initialParam*/) {
+ //
+ // Not implemented yet
+ //
return;
}
}
+AliExternalTrackParam *AliTPCcalibCosmic::Invert(AliExternalTrackParam *input)
+{
+ //
+ // Invert paramerameter - not covariance yet
+ //
+ AliExternalTrackParam *output = new AliExternalTrackParam(*input);
+ Double_t * param = (Double_t*)output->GetParameter();
+ param[0]*=-1;
+ param[3]*=-1;
+ param[4]*=-1;
+ //
+ return output;
+}
+
+AliExternalTrackParam *AliTPCcalibCosmic::MakeTrack(const AliExternalTrackParam *track0, const AliExternalTrackParam *track1){
+ //
+ //
+ //
+ AliExternalTrackParam *par1R= new AliExternalTrackParam(*track1);
+ par1R->Rotate(track0->GetAlpha());
+ //
+ //
+ Double_t * param = (Double_t*)par1R->GetParameter();
+ Double_t * covar = (Double_t*)par1R->GetCovariance();
+ param[0]*=1; //OK
+ param[1]*=1; //OK
+ param[2]*=1; //?
+ param[3]*=-1; //OK
+ param[4]*=-1; //OK
+ //
+ covar[6] *=-1.; covar[7] *=-1.; covar[8] *=-1.;
+ //covar[10]*=-1.; covar[11]*=-1.; covar[12]*=-1.;
+ covar[13]*=-1.;
+ par1R->PropagateTo(track0->GetX(),0); // bz shold be set -
+ //if (1){
+ // printf("Print param\n");
+ // track1->Print();
+ // par1R->Print();
+ //}
+ return par1R;
+}
+
+void AliTPCcalibCosmic::UpdateTrack(AliExternalTrackParam &track1, const AliExternalTrackParam &track2){
+ //
+ // Update track 1 with track 2
+ //
+ //
+ //
+ TMatrixD vecXk(5,1); // X vector
+ TMatrixD covXk(5,5); // X covariance
+ TMatrixD matHk(5,5); // vector to mesurement
+ TMatrixD measR(5,5); // measurement error
+ TMatrixD vecZk(5,1); // measurement
+ //
+ TMatrixD vecYk(5,1); // Innovation or measurement residual
+ TMatrixD matHkT(5,5);
+ TMatrixD matSk(5,5); // Innovation (or residual) covariance
+ TMatrixD matKk(5,5); // Optimal Kalman gain
+ TMatrixD mat1(5,5); // update covariance matrix
+ TMatrixD covXk2(5,5); //
+ TMatrixD covOut(5,5);
+ //
+ Double_t *param1=(Double_t*) track1.GetParameter();
+ Double_t *covar1=(Double_t*) track1.GetCovariance();
+ Double_t *param2=(Double_t*) track2.GetParameter();
+ Double_t *covar2=(Double_t*) track2.GetCovariance();
+ //
+ // copy data to the matrix
+ for (Int_t ipar=0; ipar<5; ipar++){
+ vecXk(ipar,0) = param1[ipar];
+ vecZk(ipar,0) = param2[ipar];
+ for (Int_t jpar=0; jpar<5; jpar++){
+ covXk(ipar,jpar) = covar1[track1.GetIndex(ipar, jpar)];
+ measR(ipar,jpar) = covar2[track2.GetIndex(ipar, jpar)];
+ }
+ }
+ //
+ //
+ //
+ //
+ matHk(0,0)=1; matHk(1,1)= 1; matHk(2,2)= 1;
+ matHk(3,3)= 1; matHk(4,4)= 1; // vector to measurement
+ //
+ vecYk = vecZk-matHk*vecXk; // Innovation or measurement residual
+ matHkT=matHk.T(); matHk.T();
+ matSk = (matHk*(covXk*matHkT))+measR; // Innovation (or residual) covariance
+ matSk.Invert();
+ matKk = (covXk*matHkT)*matSk; // Optimal Kalman gain
+ vecXk += matKk*vecYk; // updated vector
+ mat1(0,0)=1; mat1(1,1)=1; mat1(2,2)=1; mat1(3,3)=1; mat1(4,4)=1;
+ covXk2 = (mat1-(matKk*matHk));
+ covOut = covXk2*covXk;
+ //
+ //
+ //
+ // copy from matrix to parameters
+ if (0) {
+ vecXk.Print();
+ vecZk.Print();
+ //
+ measR.Print();
+ covXk.Print();
+ covOut.Print();
+ //
+ track1.Print();
+ track2.Print();
+ }
+
+ for (Int_t ipar=0; ipar<5; ipar++){
+ param1[ipar]= vecXk(ipar,0) ;
+ for (Int_t jpar=0; jpar<5; jpar++){
+ covar1[track1.GetIndex(ipar, jpar)]=covOut(ipar,jpar);
+ }
+ }
+}
+
+void AliTPCcalibCosmic::ProcessTree(TTree * chainTracklet, AliExternalComparison *comp){
+ //
+ // Process the debug streamer tree
+ // Possible to modify selection criteria
+ //
+ TTreeSRedirector * cstream = new TTreeSRedirector("cosmicdump.root");
+ //AliTPCcalibCosmic *cosmic = this;
+ //
+ AliExternalTrackParam * tr0 = 0;
+ AliExternalTrackParam * tr1 = 0;
+ Int_t npoints =0;
+ {
+ Int_t entries=chainTracklet->GetEntries();
+ for (Int_t i=0; i< entries; i++){
+ chainTracklet->GetBranch("Tr0.")->SetAddress(&tr0);
+ chainTracklet->GetBranch("Tr1.")->SetAddress(&tr1);
+ chainTracklet->GetEntry(i);
+ if (!tr0) continue;
+ if (!tr1) continue;
+ if (tr0->GetY()==0) continue;
+ if (tr1->GetY()==0) continue;
+ // make a local copy
+ AliExternalTrackParam par0(*tr0);
+ AliExternalTrackParam par1(*tr1);
+ AliExternalTrackParam par1R(*tr1);
+ par1R.Rotate(par1.GetAlpha()+TMath::Pi());
+ AliExternalTrackParam *par1T = MakeTrack(tr0,tr1);
+ if (0) {
+ printf("%d\t%d\t\n",i, npoints);
+ par1R.Print();
+ par1T->Print();
+ }
+ AliExternalTrackParam par0U=par0;
+ AliExternalTrackParam par1U=*par1T;
+ //
+ UpdateTrack(par0U,*par1T);
+ UpdateTrack(par1U,par0);
+ //
+ //
+ if (i%100==0) printf("%d\t%d\tt\n",i, npoints);
+ Bool_t accept = comp->AcceptPair(&par0,par1T);
+
+ if (1||fStreamLevel>0){
+ (*cstream)<<"Tracklet"<<
+ "accept="<<accept<<
+ "tr0.="<<&par0<< //original track up
+ "tr1.="<<&par1<< //original track down
+ "tr1R.="<<&par1R<< //track1 rotated to 0 frame
+ "tr1T.="<<par1T<< //track1 transformed to the track 0 frame
+ //
+ "tr0U.="<<&par0U<< //track 0 updated with track 1
+ "tr1U.="<<&par1U<< //track 1 updated with track 0
+ "\n";
+ }
+ //
+ if (accept) {
+ npoints++;
+ if (comp) comp->Process(&par0,par1T);
+ }
+ delete par1T;
+ }
+ }
+ delete cstream;
+}
+
+
+
+
+
+
+
+/*
+
+void Init(){
+
+.x ~/UliStyle.C
+.x ~/rootlogon.C
+gSystem->Load("libSTAT.so");
+gSystem->Load("libANALYSIS");
+gSystem->Load("libTPCcalib");
+gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
+
+gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
+AliXRDPROOFtoolkit tool;
+TChain * chainCosmic = tool.MakeChain("cosmic.txt","Track0",0,1000000);
+chainCosmic->Lookup();
+
+TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.01"); // OK
+TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<2"); // OK
+TCut cutP1("cutP1","abs(Tr0.fP[1]-Tr1.fP[1])<20"); // OK
+TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
+TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>50");
+TCut cutM("cutM","abs(mag)>0.01");
+TCut cutA=cutT+cutD+cutPt+cutN+cutP1+"trigger!=16";
+
+TCut cuthpt("abs(Tr0.fP[4])+abs(Tr1.fP[4])<0.2");
+TCut cutS("cutS","Orig0.fIp.fP[1]*Orig1.fIp.fP[1]>0");
+
+//
+chainCosmic->Draw(">>listELP",cutA,"entryList");
+//TEntryList *elist = (TEntryList*)gDirectory->Get("listEL");
+//TEntryList *elist = (TEntryList*)gProof->GetOutputList()->At(1);
+chainCosmic->SetEntryList(elist);
+//
+chainCosmic->Draw(">>listV40Z100","abs(d0)<40&&abs(v01)<100","entryList");
+TEntryList *elistV40Z100 = (TEntryList*)gDirectory->Get("listV40Z100");
+chainCosmic->SetEntryList(elistV40Z100);
+
+//
+// Aliases
+//
+chainCosmic->SetAlias("side","(-1+(Tr0.fP[1]>0)*2)");
+chainCosmic->SetAlias("hpt","abs(Tr0.fP[4])<0.2");
+chainCosmic->SetAlias("signy","(-1+(Tr0.fP[0]>0)*2)");
+
+chainCosmic->SetAlias("dy","Tr0.fP[0]+Tr1.fP[0]");
+chainCosmic->SetAlias("dz","Tr0.fP[1]-Tr1.fP[1]");
+chainCosmic->SetAlias("d1pt","Tr0.fP[4]+Tr1.fP[4]");
+chainCosmic->SetAlias("dtheta","(Tr0.fP[3]+Tr1.fP[3])");
+chainCosmic->SetAlias("dphi","(Tr0.fAlpha-Tr1.fAlpha-pi)");
+
+chainCosmic->SetAlias("mtheta","(Tr0.fP[3]-Tr1.fP[3])*0.5")
+chainCosmic->SetAlias("sa","(sin(Tr0.fAlpha+0.))");
+chainCosmic->SetAlias("ca","(cos(Tr0.fAlpha+0.))");
+
+
+
+chainCosmic->Draw("dy:sqrt(abs(Tr0.fP[4]))>>hisdyA(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0");
+hisdyA->FitSlicesY();
+hisdyA_2->SetXTitle("#sqrt{1/p_{t}}");
+hisdyA_2->SetYTitle("#sigma_{y}(cm)");
+hisdyA_2->SetTitle("Cosmic - Y matching");
+hisdyA_2->SetMaximum(0.5);
+
+
+chainCosmic->Draw("dy:sqrt(abs(Tr0.fP[4]))>>hisdyC(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0");
+hisdyC->FitSlicesY();
+hisdyC_2->SetXTitle("#sqrt{1/p_{t}}");
+hisdyC_2->SetYTitle("#sigma_{y}(cm)");
+hisdyC_2->SetTitle("Cosmic - Y matching");
+hisdyC_2->SetMaximum(1);
+hisdyC_2->SetMinimum(0);
+hisdyC_2->SetMarkerStyle(22);
+hisdyA_2->SetMarkerStyle(21);
+hisdyC_2->SetMarkerSize(1.5);
+hisdyA_2->SetMarkerSize(1.5);
+hisdyC_2->Draw();
+hisdyA_2->Draw("same");
+gPad->SaveAs("~/Calibration/Cosmic/pic/ymatching.gif")
+
+chainCosmic->Draw("dz:sqrt(abs(Tr0.fP[4]))>>hisdzA(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0");
+hisdzA->FitSlicesY();
+hisdzA_2->SetXTitle("#sqrt{1/p_{t}}");
+hisdzA_2->SetYTitle("#sigma_{z}(cm)");
+hisdzA_2->SetTitle("Cosmic - Z matching - A side ");
+hisdzA_2->SetMaximum(0.5);
+
+chainCosmic->Draw("dz:sqrt(abs(Tr0.fP[4]))>>hisdzC(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0");
+hisdzC->FitSlicesY();
+hisdzC_2->SetXTitle("#sqrt{1/p_{t}}");
+hisdzC_2->SetYTitle("#sigma_{z}(cm)");
+hisdzC_2->SetTitle("Cosmic - Z matching");
+hisdzC_2->SetMaximum(0.5);
+hisdzC_2->SetMarkerStyle(22);
+hisdzA_2->SetMarkerStyle(21);
+hisdzC_2->SetMarkerSize(1.5);
+hisdzA_2->SetMarkerSize(1.5);
+
+hisdzC_2->Draw();
+hisdzA_2->Draw("same");
+
+
+//
+// PICTURE 1/pt
+//
+chainCosmic->Draw("d1pt:sqrt(abs(Tr0.fP[4]))>>hisd1ptA(5,0,1,30,-0.1,0.1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0"+cutM);
+hisd1ptA->FitSlicesY();
+hisd1ptA_2->SetXTitle("#sqrt{1/p_{t}}");
+hisd1ptA_2->SetYTitle("#sigma_{z}(cm)");
+hisd1ptA_2->SetTitle("Cosmic - Z matching - A side ");
+hisd1ptA_2->SetMaximum(0.5);
+
+chainCosmic->Draw("d1pt:sqrt(abs(Tr0.fP[4]))>>hisd1ptC(5,0,1,30,-0.1,0.1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0"+cutM);
+hisd1ptC->FitSlicesY();
+hisd1ptC_2->SetXTitle("#sqrt{1/p_{t}}");
+hisd1ptC_2->SetYTitle("#sigma_{1/pt}(1/GeV)");
+hisd1ptC_2->SetTitle("Cosmic - 1/pt matching");
+hisd1ptC_2->SetMaximum(0.05);
+hisd1ptC_2->SetMarkerStyle(22);
+hisd1ptA_2->SetMarkerStyle(21);
+hisd1ptC_2->SetMarkerSize(1.5);
+hisd1ptA_2->SetMarkerSize(1.5);
+
+hisd1ptC_2->Draw();
+hisd1ptA_2->Draw("same");
+gPad->SaveAs("~/Calibration/Cosmic/pic/1ptmatching.gif")
+
+//
+// Theta
+//
+chainCosmic->Draw("dtheta:sqrt(abs(Tr0.fP[4]))>>hisdthetaA(5,0,1,30,-0.01,0.01)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0");
+hisdthetaA->FitSlicesY();
+hisdthetaA_2->SetXTitle("#sqrt{1/p_{t}}");
+hisdthetaA_2->SetYTitle("#sigma_{#theta}(cm)");
+hisdthetaA_2->SetTitle("Cosmic - Z matching - A side ");
+hisdthetaA_2->SetMaximum(0.5);
+
+chainCosmic->Draw("dtheta:sqrt(abs(Tr0.fP[4]))>>hisdthetaC(5,0,1,30,-0.01,0.01)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0");
+hisdthetaC->FitSlicesY();
+hisdthetaC_2->SetXTitle("#sqrt{1/p_{t}}");
+hisdthetaC_2->SetYTitle("#sigma_{#theta}(rad)");
+hisdthetaC_2->SetTitle("Cosmic - Theta matching");
+hisdthetaC_2->SetMaximum(0.01);
+hisdthetaC_2->SetMinimum(0.0);
+hisdthetaC_2->SetMarkerStyle(22);
+hisdthetaA_2->SetMarkerStyle(21);
+hisdthetaC_2->SetMarkerSize(1.5);
+hisdthetaA_2->SetMarkerSize(1.5);
+
+hisdthetaC_2->Draw();
+hisdthetaA_2->Draw("same");
+gPad->SaveAs("~/Calibration/Cosmic/pic/thetamatching.gif")
+//
+// Phi
+//
+chainCosmic->Draw("dphi:sqrt(abs(Tr0.fP[4]))>>hisdphiA(5,0,1,30,-0.01,0.01)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0");
+hisdphiA->FitSlicesY();
+hisdphiA_2->SetXTitle("#sqrt{1/p_{t}}");
+hisdphiA_2->SetYTitle("#sigma_{#phi}(rad)");
+hisdphiA_2->SetTitle("Cosmic - Z matching - A side ");
+hisdphiA_2->SetMaximum(0.5);
+
+chainCosmic->Draw("dphi:sqrt(abs(Tr0.fP[4]))>>hisdphiC(5,0,1,30,-0.01,0.01)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]<0");
+hisdphiC->FitSlicesY();
+hisdphiC_2->SetXTitle("#sqrt{1/p_{t}}");
+hisdphiC_2->SetYTitle("#sigma_{#phi}(rad)");
+hisdphiC_2->SetTitle("Cosmic - Phi matching");
+hisdphiC_2->SetMaximum(0.01);
+hisdphiC_2->SetMinimum(0.0);
+hisdphiC_2->SetMarkerStyle(22);
+hisdphiA_2->SetMarkerStyle(21);
+hisdphiC_2->SetMarkerSize(1.5);
+hisdphiA_2->SetMarkerSize(1.5);
+
+hisdphiC_2->Draw();
+hisdphiA_2->Draw("same");
+gPad->SaveAs("~/Calibration/Cosmic/pic/phimatching.gif")
+
+
+
+}
+
+
+*/
+
+
+/*
+void MatchTheta(){
+
+TStatToolkit toolkit;
+Double_t chi2=0;
+Int_t npoints=0;
+TVectorD fitParamA0;
+TVectorD fitParamA1;
+TVectorD fitParamC0;
+TVectorD fitParamC1;
+TMatrixD covMatrix;
+
+
+TString fstring="";
+//
+fstring+="mtheta++";
+fstring+="ca++";
+fstring+="sa++";
+fstring+="ca*mtheta++";
+fstring+="sa*mtheta++";
+//
+fstring+="side++";
+fstring+="side*mtheta++";
+fstring+="side*ca++";
+fstring+="side*sa++";
+fstring+="side*ca*mtheta++";
+fstring+="side*sa*mtheta++";
+
+
+TString *strTheta0 = toolkit.FitPlane(chain,"dtheta",fstring->Data(), "hpt&&!crossI&&!crossO", chi2,npoints,fitParamA0,covMatrix,0.8);
+chainCosmic->SetAlias("dtheta0",strTheta0.Data());
+strTheta0->Tokenize("+")->Print();
+
+
+//fstring+="mtheta++";
+//fstring+="mtheta^2++";
+//fstring+="ca*mtheta^2++";
+//fstring+="sa*mtheta^2++";
+
+
+
+}
+
+*/
+
+
+
+
+/*
+ void PosCorrection()
+
+
+
+
+ TStatToolkit toolkit;
+ Double_t chi2=0;
+ Int_t npoints=0;
+ TVectorD fitParam;
+ TMatrixD covMatrix;
+
+ //Theta
+chainCosmic->SetAlias("dthe","(Tr0.fP[3]+Tr1.fP[3])");
+chainCosmic->SetAlias("sign","(-1+(Tr0.fP[1]>0)*2)");
+chainCosmic->SetAlias("di","(1-abs(Tr0.fP[1])/250)");
+//
+//
+TString strFit="";
+//
+strFit+="sign++"; //1
+strFit+="Tr0.fP[3]++"; //2
+//
+strFit+="sin(Tr0.fAlpha)*(Tr0.fP[1]>0)++"; //3
+strFit+="sin(Tr0.fAlpha)*(Tr0.fP[1]<0)++"; //4
+//
+strFit+="cos(Tr0.fAlpha)*(Tr0.fP[1]>0)++"; //5
+strFit+="cos(Tr0.fAlpha)*(Tr0.fP[1]<0)++"; //6
+//
+strFit+="sin(Tr0.fAlpha)*Tr0.fP[3]++"; //7
+strFit+="cos(Tr0.fAlpha)*Tr0.fP[3]++"; //8
+
+
+ //
+ TString * thetaParam = toolkit.FitPlane(chain,"dthe", strFit.Data(),"1", chi2,npoints,fitParam,covMatrix,0.8,0,10000)
+
+ chainCosmic->SetAlias("corTheta",thetaParam->Data());
+ chainCosmic->Draw("dthe:Tr0.fP[1]","","",50000);
+
+
+
+*/
/*
void AliTPCcalibCosmic::dEdxCorrection(){
- TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03"); // OK
- TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5"); // OK
- TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.2&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
- TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>70");
- TCut cutA=cutT+cutD+cutPt+cutN;
+ TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.01"); // OK
+ TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<2"); // OK
+ TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
+ TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>100");
+ TCut cutS("cutS","Orig0.fIp.fP[1]*Orig1.fIp.fP[1]>0");
+ TCut cutA=cutT+cutD+cutPt+cutN+cutS;
.x ~/UliStyle.C
gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
AliXRDPROOFtoolkit tool;
- TChain * chain = tool.MakeChain("cosmic.txt","Track0",0,1000000);
- chain->Lookup();
+ TChain * chainCosmic = tool.MakeChain("cosmic.txt","Track0",0,1000000);
+ chainCosmic->Lookup();
+
+ chainCosmic->Draw(">>listEL",cutA,"entryList");
+ TEntryList *elist = (TEntryList*)gDirectory->Get("listEL");
+ chainCosmic->SetEntryList(elist);
.x ~/rootlogon.C
gSystem->Load("libSTAT.so");
-
+ TStatToolkit toolkit;
Double_t chi2=0;
Int_t npoints=0;
TVectorD fitParam;
TMatrixD covMatrix;
- chain->Draw("Tr0.fP[4]+Tr1.fP[4]","OK"+cutA);
+ chainCosmic->Draw("Tr0.fP[4]+Tr1.fP[4]","OK"+cutA);
TString strFit;
strFit+="(Tr0.fP[1]/250)++";
strFit+="(Tr0.fP[3])++";
strFit+="(Tr0.fP[3])^2++";
- TString * ptParam = TStatToolkit::FitPlane(chain,"Tr0.fP[4]+Tr1.fP[4]", strFit.Data(),cutA, chi2,npoints,fitParam,covMatrix)
+ TString * ptParam = TStatToolkit::FitPlane(chain,"Tr0.fP[4]+Tr1.fP[4]", strFit.Data(),"1", chi2,npoints,fitParam,covMatrix)
+
-strFit+="(Tr0.fP[1]/250)++";
-strFit+="(Tr0.fP[1]/250)^2++";
-strFit+="(Tr0.fP[3])++";
-strFit+="(Tr0.fP[3])^2++";
-strFit+="(Tr0.fP[1]/250)^2*Tr0.fP[3]++";
-strFit+="(Tr0.fP[1]/250)^2*Tr0.fP[3]^2++";
+
+*/
+
+
+/*
+gSystem->Load("libANALYSIS");
+gSystem->Load("libSTAT");
+gSystem->Load("libTPCcalib");
+
+TStatToolkit toolkit;
+Double_t chi2;
+TVectorD fitParam;
+TMatrixD covMatrix;
+Int_t npoints;
//
+TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03"); // OK
+TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5"); // OK
+TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.2&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
+TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>110");
+TCut cutA=cutT+cutD+cutPt+cutN;
-strFit+="sign(Tr0.fP[1])++"
-strFit+="sign(Tr0.fP[1])*(1-abs(Tr0.fP[1]/250))"
-
-TString * thetaParam = TStatToolkit::FitPlane(chain,"Tr0.fP[3]+Tr1.fP[3]", strFit.Data(),cutA, chi2,npoints,fitParam,covMatrix)
+TTree * chainCosmic = Track0;
+
+
+chainCosmic->SetAlias("norm","signalTot0.fElements[3]/signalTot1.fElements[3]");
+//
+chainCosmic->SetAlias("dr1","(signalTot0.fElements[1]/signalTot0.fElements[3])");
+chainCosmic->SetAlias("dr2","(signalTot0.fElements[2]/signalTot0.fElements[3])");
+chainCosmic->SetAlias("dr4","(signalTot0.fElements[4]/signalTot0.fElements[3])");
+chainCosmic->SetAlias("dr5","(signalTot0.fElements[5]/signalTot0.fElements[3])");
+
+TString fstring="";
+fstring+="dr1++";
+fstring+="dr2++";
+fstring+="dr4++";
+fstring+="dr5++";
+//
+fstring+="dr1*dr2++";
+fstring+="dr1*dr4++";
+fstring+="dr1*dr5++";
+fstring+="dr2*dr4++";
+fstring+="dr2*dr5++";
+fstring+="dr4*dr5++";
+
-(-0.009263+(Tr0.fP[1]/250)*(-0.009693)+(Tr0.fP[1]/250)^2*(0.009062)+(Tr0.fP[3])*(0.002256)+(Tr0.fP[3])^2*(-0.052775)+(Tr0.fP[1]/250)^2*Tr0.fP[3]*(-0.020627)+(Tr0.fP[1]/250)^2*Tr0.fP[3]^2*(0.049030))
+TString *strqdedx = toolkit.FitPlane(chain,"norm",fstring->Data(), cutA, chi2,npoints,fitParam,covMatrix,-1,0,200000);
+
+chainCosmic->SetAlias("corQT",strqdedx->Data());
+
+*/
+
+
+/*
+ chainCosmic->SetProof(kTRUE);
+ chainCosmic->Draw("Seed0.CookdEdxNorm(0,0.6,1,0,159,0,kTRUE,kTRUE):Seed0.CookdEdxNorm(0,0.6,1,0,159,0,kFALSE,kTRUE)",""+cutA,"",100000);
+
+
+chainCosmic->Draw("Seed0.CookdEdxNorm(0,0.6,1,0,159,0,kTRUE,kTRUE)/Seed1.CookdEdxNorm(0,0.6,1,0,159,0,kTRUE,kTRUE)>>his(100,0.5,1.5)","min(Orig0.fTPCncls,Orig1.fTPCncls)>130"+cutA,"",50000);
+
+*/
+
+
+/*
+chainCosmic->Draw("Tr0.fP[1]-Tr1.fP[1]:sqrt(abs(Tr0.fP[4]))>>hisdzA(5,0,1,50,-1,1)","!crossO&&!crossI&&abs(d0)<40&&Tr0.fP[1]>0&&abs(mag)>0.1"+cutA);
+
+TGraph *grdzA = (TGraph*)gProof->GetOutputList()->At(1)->Clone();
+
+
+
+
*/
+
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff