#include "TVectorD.h"
#include "AliLog.h"
#include "AliTPCROC.h"
+#include "AliTPCClusterParam.h"
#include "AliTPCPointCorrection.h"
AliTPCPointCorrection* AliTPCPointCorrection::fgInstance = 0;
fParamOutRVersion(0),
fErrorsOutR(38),
fErrorsOutZ(38),
- fParamOutZVersion(0)
+ fParamOutZVersion(0),
+ //
+ fXIO(0),
+ fXmiddle(0),
+ fXquadrant(0),
+ fArraySectorIntParam(36), // array of sector alignment parameters
+ fArraySectorIntCovar(36), // array of sector alignment covariances
+ //
+ // Kalman filter for global alignment
+ //
+ fSectorParam(0), // Kalman parameter
+ fSectorCovar(0) // Kalman covariance
+
{
//
// Default constructor
//
+ AliTPCROC * roc = AliTPCROC::Instance();
+ fXquadrant = roc->GetPadRowRadii(36,53);
+ fXmiddle = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
+ fXIO = ( roc->GetPadRowRadii(0,roc->GetNRows(0)-1)+roc->GetPadRowRadii(36,0))*0.5;
}
AliTPCPointCorrection::AliTPCPointCorrection(const Text_t *name, const Text_t *title):
fParamOutRVersion(0),
fErrorsOutR(38),
fErrorsOutZ(38),
- fParamOutZVersion(0)
+ fParamOutZVersion(0),
+ //
+ //
+ //
+ fXIO(0),
+ fXmiddle(0),
+ fXquadrant(0),
+ fArraySectorIntParam(36), // array of sector alignment parameters
+ fArraySectorIntCovar(36), // array of sector alignment covariances
+ //
+ // Kalman filter for global alignment
+ //
+ fSectorParam(0), // Kalman parameter for A side
+ fSectorCovar(0) // Kalman covariance for A side
+
{
//
// Non default constructor
//
+ AliTPCROC * roc = AliTPCROC::Instance();
+ fXquadrant = roc->GetPadRowRadii(36,53);
+ fXmiddle = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
+ fXIO = ( roc->GetPadRowRadii(0,roc->GetNRows(0)-1)+roc->GetPadRowRadii(36,0))*0.5;
+
}
AliTPCPointCorrection::~AliTPCPointCorrection(){
Double_t radius = TMath::Sqrt(cx*cx+cy*cy);
AliTPCROC * roc = AliTPCROC::Instance();
Double_t xout = roc->GetPadRowRadiiUp(roc->GetNRows(37)-1);
- Double_t dout = xout-radius;
+ Double_t dout = xout-radius;
+ if (dout<0) return 0;
//drift
Double_t dr = 0.5 - TMath::Abs(cz/250.);
//
result+=(*vec)[5]*eoutp*dr*dr;
result+=(*vec)[6]*eoutm*dr*dr;
return result;
+}
+Double_t AliTPCPointCorrection::RPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky,Float_t qMax, Float_t threshold){
+ //
+ // Calculates COG corection in RPHI direction
+ // cluster and track position y is supposed to be corrected before for other effects
+ // (e.g ExB and alignemnt)
+ // Rphi distortion dependeds on the distance to the edge-pad, distance to the wire edge and
+ // relative distance to the center of the pad. Therefore the y position is trnsfromed to the
+ // pad coordiante frame (correction offset (ExB alignemnt) substracted).
+ //
+ // Input parameters:
+ //
+ // sector - sector number - 0-71 - cl.GetDetector()
+ // padrow - padrow number -0-63 -IROC 0-95 OROC - cl->GetRow()
+ // pad - mean pad number - cl->GetPad()
+ // cy - cluster y - cl->GetY()
+ // y - track -or cluster y
+ // qMax - cluster max charge - cl-.GetMax()
+ // threshold - clusterer threshold
+ //
+ AliTPCClusterParam * clparam = AliTPCClusterParam::Instance();
+ Int_t padtype=0;
+ if (sector>=36) padtype = (padrow>64)?2:1;
+ Double_t padwidth=(padtype==0)? 0.4:0.6;
+
+ Double_t sigma = clparam->GetRMS0(0,padtype,250-TMath::Abs(z),ky)/padwidth;
+ //
+ Int_t npads = AliTPCROC::Instance()->GetNPads(sector,padrow);
+ Float_t yshift = TMath::Abs(cy)-TMath::Abs((-npads*0.5+pad)*padwidth); // the clusters can be corrected before
+ // shift to undo correction
+
+ y -= yshift*((y>0)?1.:-1.); // y in the sector coordinate
+ Double_t y0 = npads*0.5*padwidth;
+ Double_t dy = (TMath::Abs(y0)-TMath::Abs(y))/padwidth-0.5; // distance to the center of pad0
+ //
+ Double_t padcenter = TMath::Nint(dy);
+ Double_t sumw=0;
+ Double_t sumwi=0;
+ for (Double_t ip=padcenter-2.5; ip<=padcenter+2.5;ip++){
+ Double_t weightGaus = qMax*TMath::Exp(-(dy-ip)*(dy-ip)/(2*sigma*sigma));
+ Double_t ypad = (ip-npads*0.5)*padwidth;
+ Double_t weightGain = GetEdgeQ0(sector,padrow,ypad);
+ //
+ Double_t weight = TMath::Nint(weightGaus*weightGain);
+ if (ip<0 &&weight> threshold) weight = threshold; // this is done in cl finder
+ if (weight<0) continue;
+ sumw+=weight;
+ sumwi+=weight*(ip);
+ }
+ Double_t result =0;
+ if (sumw>0) result = sumwi/sumw;
+ result = (result-dy)*padwidth;
+ return result;
}
+
+
+
Double_t AliTPCPointCorrection::CorrectionOutZ0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector){
//
// return dR correction - for correction version 0
AliTPCROC * roc = AliTPCROC::Instance();
Double_t xout = roc->GetPadRowRadiiUp(roc->GetNRows(37)-1);
Double_t dout = xout-radius;
+ if (dout<0) return 0;
//drift
Double_t dr = 0.5 - TMath::Abs(cz/250.);
//
}
+Double_t AliTPCPointCorrection::GetEdgeQ0(Int_t sector, Int_t padrow, Float_t y){
+ //
+ /* TF1 fexp("fexp","1-exp(-[0]*(x-[1]))",0,20)
+ | param [0] | param [1]
+ IROC | 4.71413 | 1.39558
+ OROC1| 2.11437 | 1.52643
+ OROC2| 2.15082 | 1.53537
+ */
+
+ Double_t params[2]={100,0};
+ if (sector<36){
+ params[0]=4.71413; params[1]=1.39558;
+ }
+ if (sector>=36){
+ if (padrow<64) { params[0]=2.114; params[1]=1.526;}
+ if (padrow>=64){ params[0]=2.15; params[1]=1.535;}
+ }
+ Double_t result= 1;
+ Double_t xrow = AliTPCROC::Instance()->GetPadRowRadii(sector,padrow);
+ Double_t ymax = TMath::Tan(TMath::Pi()/18.)*xrow;
+ Double_t dedge = ymax-TMath::Abs(y);
+ if (dedge-params[2]<0) return 0;
+ if (dedge>10.*params[1]) return 1;
+ result = 1.-TMath::Exp(-params[0]*(dedge-params[2]));
+ return result;
+}
+
+Double_t AliTPCPointCorrection::SRPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky,Float_t qMax, Float_t threshold){
+ return fgInstance->RPhiCOGCorrection(sector, padrow, pad, cy, y,z, ky, qMax, threshold);
+}
+
+Double_t AliTPCPointCorrection::SGetEdgeQ0(Int_t sector, Int_t padrow, Float_t y){
+ //
+ //
+ return fgInstance->GetEdgeQ0(sector, padrow, y);
+}
+
+void AliTPCPointCorrection::AddCorrectionSector(TObjArray & sideAPar, TObjArray &sideCPar, TObjArray & sideACov, TObjArray &sideCCov, Bool_t reset){
+ //
+ //
+ //
+ for (Int_t isec=0;isec<36;isec++){
+ TMatrixD * mat1 = (TMatrixD*)(sideAPar.At(isec));
+ TMatrixD * mat1Covar = (TMatrixD*)(sideCPar.At(isec));
+ if (!mat1) continue;
+ TMatrixD * mat0 = (TMatrixD*)(fArraySectorIntParam.At(isec));
+ TMatrixD * mat0Covar = (TMatrixD*)(fArraySectorIntCovar.At(isec));
+ if (!mat0) {
+ fArraySectorIntParam.AddAt(mat1->Clone(),isec);
+ fArraySectorIntCovar.AddAt(mat1Covar->Clone(),isec);
+ continue;
+ }
+ (*mat0Covar)=(*mat1Covar);
+ if (reset) (*mat0)=(*mat1);
+ if (!reset) (*mat0)+=(*mat1);
+ }
+
+ for (Int_t isec=0;isec<36;isec++){
+ TMatrixD * mat1 = (TMatrixD*)(sideCPar.At(isec));
+ TMatrixD * mat1Covar = (TMatrixD*)(sideCCov.At(isec));
+ if (!mat1) continue;
+ TMatrixD * mat0 = (TMatrixD*)(fArraySectorIntParam.At(isec));
+ TMatrixD * mat0Covar = (TMatrixD*)(fArraySectorIntCovar.At(isec));
+ if (!mat0) {
+ fArraySectorIntParam.AddAt(mat1->Clone(),isec);
+ fArraySectorIntCovar.AddAt(mat1Covar->Clone(),isec);
+ continue;
+ }
+ (*mat0Covar)=(*mat1Covar);
+ if (reset) (*mat0)=(*mat1);
+ if (!reset) (*mat0)+=(*mat1);
+ }
+}
+
+
+Double_t AliTPCPointCorrection::GetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant){
+ //
+ // Get position correction for given sector
+ //
+
+ TMatrixD * param = (TMatrixD*)fArraySectorIntParam.At(sector%36);
+ if (!param) return 0;
+ if (quadrant<0){ //recalc quadrant if not specified
+ if (lx<fXIO) quadrant=0;
+ if(lx>fXIO){
+ if (lx<fXquadrant) {
+ if (ly<0) quadrant=1;
+ if (ly>0) quadrant=2;
+ }
+ if (lx>=fXquadrant) {
+ if (ly<0) quadrant=3;
+ if (ly>0) quadrant=4;
+ }
+ }
+ }
+ Double_t a10 = (*param)(quadrant*6+0,0);
+ Double_t a20 = (*param)(quadrant*6+1,0);
+ Double_t a21 = (*param)(quadrant*6+2,0);
+ Double_t dx = (*param)(quadrant*6+3,0);
+ Double_t dy = (*param)(quadrant*6+4,0);
+ Double_t dz = (*param)(quadrant*6+5,0);
+ Double_t deltaX = lx-fXmiddle;
+ if (coord==0) return dx;
+ if (coord==1) return dy+deltaX*a10;
+ if (coord==2) return dz+deltaX*a20+ly*a21;
+ if (coord==3) return a10;
+ if (coord==4) return a20;
+ if (coord==5) return a21;
+ //
+ return 0;
+}
+
+Double_t AliTPCPointCorrection::SGetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant){
+ //
+ //
+ //
+ if (!Instance()) return 0;
+ return Instance()->GetCorrectionSector(coord,sector,lx,ly,lz, quadrant);
+}
+
+
+
+Double_t AliTPCPointCorrection::GetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz){
+ //
+ // Get position correction for given sector
+ //
+ if (!fSectorParam) return 0;
+
+ Double_t a10 = (*fSectorParam)(sector*6+0,0);
+ Double_t a20 = (*fSectorParam)(sector*6+1,0);
+ Double_t a21 = (*fSectorParam)(sector*6+2,0);
+ Double_t dx = (*fSectorParam)(sector*6+3,0);
+ Double_t dy = (*fSectorParam)(sector*6+4,0);
+ Double_t dz = (*fSectorParam)(sector*6+5,0);
+ Double_t deltaX = lx-fXIO;
+ //
+ if (coord==0) return dx;
+ if (coord==1) return dy+deltaX*a10;
+ if (coord==2) return dz+deltaX*a20+ly*a21;
+ if (coord==3) return a10;
+ if (coord==4) return a20;
+ if (coord==5) return a21;
+ return 0;
+}
+
+Double_t AliTPCPointCorrection::SGetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz){
+ //
+ //
+ //
+ if (!Instance()) return 0;
+ return Instance()->GetCorrection(coord,sector,lx,ly,lz);
+}
+
+
+
+
+
#include "TObjArray.h"
#include "TVectorD.h"
-
class AliTPCPointCorrection:public TNamed {
public:
Double_t CorrectionOutR0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector);
Double_t CorrectionOutZ0(Bool_t isGlobal, Bool_t type, Double_t cx, Double_t cy, Double_t cz, Int_t sector);
+ Double_t RPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky, Float_t qMax, Float_t threhsold);
+ Double_t SRPhiCOGCorrection(Int_t sector, Int_t padrow, Float_t pad, Float_t cy, Float_t y, Float_t z, Float_t ky, Float_t qMax, Float_t threhsold);
+ //
+ Double_t GetEdgeQ0(Int_t sector, Int_t padrow, Float_t y);
+ static Double_t SGetEdgeQ0(Int_t sector, Int_t padrow, Float_t y);
+ //
+ // IROC -OROC+Quadrant alignment
+ //
+ void AddCorrectionSector(TObjArray & sideAPar, TObjArray &sideCPar, TObjArray & sideACov, TObjArray &sideCCov, Bool_t reset);
+ Double_t GetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant =-1);
+ static Double_t SGetCorrectionSector(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz, Int_t quadrant=-1);
+ //
+ // Global alignment
+ //
+ Double_t GetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz);
+ static Double_t SGetCorrection(Int_t coord, Int_t sector, Double_t lx, Double_t ly, Double_t lz);
public:
//
// Correction out
TObjArray fErrorsOutZ; // Parameters for z distortion - outer field cage
Int_t fParamOutZVersion; // version of the parameterization
//
+ // Edge rfi
+ //
+ //
+ // Alignment part
+ //
+ Double_t fXIO; // OROC-IROC boundary
+ Double_t fXmiddle; // center of the TPC sector local X
+ Double_t fXquadrant; // x quadrant
+ //
+ // IROC OROC alignment
+ //
+ TObjArray fArraySectorIntParam; // array of sector alignment parameters
+ TObjArray fArraySectorIntCovar; // array of sector alignment covariances
+ //
+ // Kalman filter for global alignment
+ //
+ TMatrixD *fSectorParam; // Kalman parameter
+ TMatrixD *fSectorCovar; // Kalman covariance
+ //
+ //
+ //
private:
+
AliTPCPointCorrection(const AliTPCPointCorrection&);
AliTPCPointCorrection& operator=(const AliTPCPointCorrection&);
static AliTPCPointCorrection* fgInstance; //! Instance of this class (singleton implementation)
- ClassDef(AliTPCPointCorrection, 1);
+ ClassDef(AliTPCPointCorrection, 3);
};
#endif
fTriggerMaskAccept(-1), //trigger mask - accept trigger
fHasLaser(kFALSE), //flag the laser is overlayed with given event
fRejectLaser(kTRUE), //flag- reject laser
+ fTriggerClass(),
fDebugLevel(0)
{
//
fTriggerMaskAccept(-1), //trigger mask - accept trigger
fHasLaser(kFALSE), //flag the laser is overlayed with given event
fRejectLaser(kTRUE), //flag- reject laser
+ fTriggerClass(),
fDebugLevel(0)
{
//
fTime = event->GetTimeStamp();
fTrigger = event->GetTriggerMask();
fMagF = event->GetMagneticField();
+ fTriggerClass = event->GetFiredTriggerClasses().Data();
fHasLaser = HasLaser(event);
+
}
////
#include "TNamed.h"
+#include "TObjString.h"
class AliTPCseed;
class AliESDEvent;
class AliESDtrack;
Int_t fTriggerMaskAccept; //trigger mask - accept
Bool_t fHasLaser; //flag the laser is overlayed with given event
Bool_t fRejectLaser; //flag- reject laser
+ TObjString fTriggerClass; // trigger class
private:
Int_t fDebugLevel; // debug level
#include "AliTPCTransform.h"
#include "AliTPCclusterMI.h"
#include "AliTPCseed.h"
+#include "AliTPCPointCorrection.h"
ClassImp(AliTPCcalibCalib)
// Refit track
//
+ //
+ // 0 - Setup transform object
+ //
+ AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
+ transform->SetCurrentRun(fRun);
+ transform->SetCurrentTimeStamp((UInt_t)fTime);
//
// First apply calibration
//
-
- AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform();
+ AliTPCPointCorrection * corr = AliTPCPointCorrection::Instance();
for (Int_t irow=0;irow<159;irow++) {
AliTPCclusterMI *cluster=seed->GetClusterPointer(irow);
if (!cluster) continue;
Float_t dz =AliTPCClusterParam::SPosCorrection(1,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
//x[1]-=dy;
//x[2]-=dz;
+ //
+ // Apply sector alignment
+ //
+ Double_t dxq = AliTPCPointCorrection::SGetCorrectionSector(0,cluster->GetDetector()%36,cluster->GetX(),
+ cluster->GetY(),cluster->GetZ());
+ Double_t dyq = AliTPCPointCorrection::SGetCorrectionSector(1,cluster->GetDetector()%36,cluster->GetX(),
+ cluster->GetY(),cluster->GetZ());
+ Double_t dzq = AliTPCPointCorrection::SGetCorrectionSector(2,cluster->GetDetector()%36,cluster->GetX(),
+ cluster->GetY(),cluster->GetZ());
+ if (kTRUE){
+ x[0]-=dxq;
+ x[1]-=dyq;
+ x[2]-=dzq;
+ }
+ //
+ // Apply r-phi correction - To be done on track level- knowing the track angle !!!
+ //
+ Double_t corrclY =
+ corr->RPhiCOGCorrection(cluster->GetDetector(),cluster->GetRow(), cluster->GetPad(),
+ cluster->GetY(),cluster->GetY(), cluster->GetZ(), 0., cluster->GetMax(),2.5);
+ // R correction
+ Double_t corrR = corr->CorrectionOutR0(kFALSE,kFALSE,cluster->GetX(),cluster->GetY(),cluster->GetZ(),cluster->GetDetector());
+
+ if (kTRUE){
+ if (cluster->GetY()>0) x[1]+=corrclY; // rphi correction
+ if (cluster->GetY()<0) x[1]-=corrclY; // rphi correction
+ x[0]+=corrR; // radial correction
+ }
+ //
+ //
//
cluster->SetX(x[0]);
cluster->SetY(x[1]);
"event="<<fEvent<< // event number
"time="<<fTime<< // time stamp of event
"trigger="<<fTrigger<< // trigger
+ "triggerClass="<<&fTriggerClass<< // trigger
"mag="<<fMagF<< // magnetic field
"cl0.="<<&cl0<<
"cl.="<<cluster<<
"cy="<<dy<<
"cz="<<dz<<
+ "cR="<<corrR<<
+ "dxq="<<dxq<<
+ "dyq="<<dyq<<
+ "dzq="<<dzq<<
"\n";
}
}
}
+ //
+ //
+ //
Int_t ncl = seed->GetNumberOfClusters();
Double_t covar[15];
for (Int_t i=0;i<15;i++) covar[i]=0;
AliExternalTrackParam * trackOutOld = (AliExternalTrackParam*)track->GetOuterParam();
AliExternalTrackParam trackIn = *trackOutOld;
- AliExternalTrackParam trackOut = *trackInOld;
trackIn.ResetCovariance(200.);
- trackOut.ResetCovariance(200.);
trackIn.AddCovariance(covar);
- trackOut.AddCovariance(covar);
Double_t xyz[3];
Int_t nclIn=0,nclOut=0;
//
+ // Refit in
+ //
+
+ for (Int_t irow=159; irow>0; irow--){
+ AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+ if (!cl) continue;
+ if (cl->GetX()<80) continue;
+ Int_t sector = cl->GetDetector();
+ Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha();
+ if (TMath::Abs(dalpha)>0.01)
+ trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.));
+ Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
+ Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
+ AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]);
+ cov[0]*=cov[0];
+ cov[2]*=cov[2];
+ trackIn.GetXYZ(xyz);
+ Double_t bz = AliTracker::GetBz(xyz);
+
+ if (!trackIn.PropagateTo(r[0],bz)) continue;
+ if (RejectCluster(cl,&trackIn)) continue;
+ nclIn++;
+ trackIn.Update(&r[1],cov);
+ }
+ //
+ AliExternalTrackParam trackOut = trackIn;
+ trackOut.ResetCovariance(200.);
+ trackOut.AddCovariance(covar);
+ //
// Refit out
//
+ //Bool_t lastEdge=kFALSE;
for (Int_t irow=0; irow<160; irow++){
AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
if (!cl) continue;
cov[2]*=cov[2];
trackOut.GetXYZ(xyz);
Double_t bz = AliTracker::GetBz(xyz);
- if (trackOut.PropagateTo(r[0],bz)) nclOut++;
+ if (!trackOut.PropagateTo(r[0],bz)) continue;
if (RejectCluster(cl,&trackOut)) continue;
- trackOut.Update(&r[1],cov);
+ nclOut++;
+ trackOut.Update(&r[1],cov);
+ //if (cl->GetType()<0) lastEdge=kTRUE;
+ //if (cl->GetType()>=0) lastEdge=kFALSE;
}
//
- // Refit in
//
-
+ //
+ nclIn=0;
+ trackIn = trackOut;
+ trackIn.ResetCovariance(10.);
+ //
+ // Refit in one more time
+ //
for (Int_t irow=159; irow>0; irow--){
AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
if (!cl) continue;
trackIn.GetXYZ(xyz);
Double_t bz = AliTracker::GetBz(xyz);
- if (trackIn.PropagateTo(r[0],bz)) nclIn++;
+ if (!trackIn.PropagateTo(r[0],bz)) continue;
if (RejectCluster(cl,&trackIn)) continue;
+ nclIn++;
trackIn.Update(&r[1],cov);
}
+
+
trackIn.Rotate(trackInOld->GetAlpha());
trackOut.Rotate(trackOutOld->GetAlpha());
//
"event="<<fEvent<< // event number
"time="<<fTime<< // time stamp of event
"trigger="<<fTrigger<< // trigger
+ "triggerClass="<<&fTriggerClass<< // trigger
"mag="<<fMagF<< // magnetic field
"nclIn="<<nclIn<<
"nclOut="<<nclOut<<
AliExternalTrackParam *t = &trackIn;
track->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
seed->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
+ seed->SetNumberOfClusters((nclIn+nclOut)*0.5);
return kTRUE;
}
// check the acceptance of cluster
// Cut on edge effects
//
+ Float_t kEdgeCut=2.5;
+ Float_t kSigmaCut=6;
+
Bool_t isReject = kFALSE;
Float_t edgeY = cl->GetX()*TMath::Tan(TMath::Pi()/18);
Float_t dist = edgeY - TMath::Abs(cl->GetY());
if (param) dist = TMath::Abs(edgeY - TMath::Abs(param->GetY()));
- if (dist<3) isReject=kTRUE;
- if (cl->GetType()<0) isReject=kTRUE;
+ if (dist<kEdgeCut) isReject=kTRUE;
+
+ Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
+ AliTPCseed::GetError(cl, param,cov[0],cov[2]);
+ Double_t py = (cl->GetY()-param->GetY())/TMath::Sqrt(cov[0]*cov[0]+param->GetSigmaY2());
+ Double_t pz = (cl->GetZ()-param->GetZ())/TMath::Sqrt(cov[2]*cov[2]+param->GetSigmaZ2());
+ //
+ if ((py*py+pz*pz)>kSigmaCut*kSigmaCut) isReject=kTRUE;
+
return isReject;
}
#include "AliESDfriend.h"
#include "AliESDInputHandler.h"
#include "AliAnalysisManager.h"
+#include "AliExternalComparison.h"
#include "AliTracker.h"
#include "AliMagF.h"
#include "AliLog.h"
#include "AliTPCcalibCosmic.h"
-#include "AliExternalComparison.h"
#include "TTreeStream.h"
#include "AliTPCTracklet.h"
+#include "AliESDcosmic.h"
+
ClassImp(AliTPCcalibCosmic)
AliTPCcalibCosmic::AliTPCcalibCosmic()
:AliTPCcalibBase(),
- fGainMap(0),
+ fComp(0),
fHistNTracks(0),
fClusters(0),
fModules(0),
AliTPCcalibCosmic::AliTPCcalibCosmic(const Text_t *name, const Text_t *title)
:AliTPCcalibBase(),
- fGainMap(0),
+ fComp(0),
fHistNTracks(0),
fClusters(0),
fModules(0),
//
//
//
+ delete fComp;
}
Int_t ntracks=event->GetNumberOfTracks();
fHistNTracks->Fill(ntracks);
if (ntracks==0) return;
-
+ AliESDcosmic cosmicESD;
+ TTreeSRedirector * cstream = GetDebugStreamer();
+ cosmicESD.SetDebugStreamer(cstream);
+ cosmicESD.ProcessEvent(event);
+ if (cstream) cosmicESD.DumpToTree();
+
+
}
Double_t meanP = 0.5*(trackIn->GetP() + trackOut->GetP());
if (seed && track->GetTPCNcls() > 80 + 60/(1+TMath::Exp(-meanP+5))) {
- fDeDx->Fill(meanP, seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap));
+ fDeDx->Fill(meanP, seed->CookdEdxNorm(0.0,0.45,0,0,159));
//
- if (meanP > 0.4 && meanP < 0.45) fDeDxMIP->Fill(seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap));
+ if (meanP > 0.4 && meanP < 0.45) fDeDxMIP->Fill(seed->CookdEdxNorm(0.0,0.45,0,0,159));
//
- if (GetDebugLevel()>0&&meanP>0.2&&seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap)>300) {
+ if (GetDebugLevel()>0&&meanP>0.2&&seed->CookdEdxNorm(0.0,0.45,0,0,159)>300) {
TFile *curfile = AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile();
if (curfile) printf(">>> p+ in file: %s \t event: %i \t Number of ESD tracks: %i \n", curfile->GetName(), (int)event->GetEventNumberInFile(), (int)ntracks);
if (track->GetOuterParam()->GetAlpha()<0) cout << " Polartiy: " << track->GetSign() << endl;
AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j);
if (! seed0) continue;
if (! seed1) continue;
- Float_t dedx0 = seed0->CookdEdxNorm(0.05,0.55,0,0,159,fGainMap);
- Float_t dedx1 = seed1->CookdEdxNorm(0.05,0.55,0,0,159,fGainMap);
+ Float_t dedx0 = seed0->CookdEdxNorm(0.05,0.55,0,0,159);
+ Float_t dedx1 = seed1->CookdEdxNorm(0.05,0.55,0,0,159);
//
- Float_t dedx0I = seed0->CookdEdxNorm(0.05,0.55,0,0,63,fGainMap);
- Float_t dedx1I = seed1->CookdEdxNorm(0.05,0.55,0,0,63,fGainMap);
+ Float_t dedx0I = seed0->CookdEdxNorm(0.05,0.55,0,0,63);
+ Float_t dedx1I = seed1->CookdEdxNorm(0.05,0.55,0,0,63);
//
- 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 dedx0O = seed0->CookdEdxNorm(0.05,0.55,0,64,159);
+ Float_t dedx1O = seed1->CookdEdxNorm(0.05,0.55,0,64,159);
//
Float_t dir = (dir0[0]*dir1[0] + dir0[1]*dir1[1] + dir0[2]*dir1[2]);
Float_t d0 = track0->GetLinearD(0,0);
Bool_t isPair = IsPair(¶m0,¶m1);
//
if (isPair) FillAcordeHist(track0);
+ if (fComp){
+ Bool_t acceptComp = fComp->AcceptPair(¶m0,¶m1);
+ if (acceptComp) fComp->Process(¶m0,¶m1);
+ }
+ //
+ // combined track params
+ //
+ AliExternalTrackParam *par0U=MakeCombinedTrack(¶m0,¶m1);
+ AliExternalTrackParam *par1U=MakeCombinedTrack(¶m1,¶m0);
+
+
//
if (fStreamLevel>0){
TTreeSRedirector * cstream = GetDebugStreamer();
"event="<<fEvent<< // event number
"time="<<fTime<< // time stamp of event
"trigger="<<fTrigger<< // trigger
+ "triggerClass="<<&fTriggerClass<< // trigger
"mag="<<fMagF<< // magnetic field
"dir="<<dir<< // direction
"OK="<<isPair<< // will be accepted
"Ip1.="<<ip1<< // inner param - lower
"Op0.="<<op0<< // outer param - upper
"Op1.="<<op1<< // outer param - lower
+ "Up0.="<<par0U<< // combined track 0
+ "Up1.="<<par1U<< // combined track 1
//
"v00="<<dvertex0[0]<< // distance using kalman
"v01="<<dvertex0[1]<< //
"dedx1O="<<dedx1O<< // dedx1 - outer ROC
"\n";
}
- }
+ }
+ delete par0U;
+ delete par1U;
}
}
}
fDeDxMIP->Add(cal->GetHistMIP());
}
-
+ //if (fComp && cal->fComp) fComp->Add(cal->fComp);
return 0;
}
}
-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());
+ par1R->PropagateTo(track0->GetX(),AliTracker::GetBz());
//
//
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; //?
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;
}
+AliExternalTrackParam *AliTPCcalibCosmic::MakeCombinedTrack(const AliExternalTrackParam *track0, const AliExternalTrackParam *track1){
+ //
+ // Make combined track
+ //
+ //
+ AliExternalTrackParam * par1T = MakeTrack(track0,track1);
+ AliExternalTrackParam * par0U = new AliExternalTrackParam(*track0);
+ //
+ UpdateTrack(*par0U,*par1T);
+ delete par1T;
+ return par0U;
+}
+
+
void AliTPCcalibCosmic::UpdateTrack(AliExternalTrackParam &track1, const AliExternalTrackParam &track2){
//
// Update track 1 with track 2
//
// 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(ipar,jpar)=0;
+ mat1(ipar,jpar)=0;
}
+ vecXk(ipar,0) = param1[ipar];
+ vecZk(ipar,0) = param2[ipar];
+ matHk(ipar,ipar)=1;
+ mat1(ipar,ipar)=0;
}
//
//
//
//
- 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.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;
//
-
-
-
-/*
-
-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.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->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();
-
- 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;
-
- chainCosmic->Draw("Tr0.fP[4]+Tr1.fP[4]","OK"+cutA);
-
- TString strFit;
- strFit+="(Tr0.fP[1]/250)++";
- strFit+="(Tr0.fP[1]/250)^2++";
- strFit+="(Tr0.fP[3])++";
- strFit+="(Tr0.fP[3])^2++";
-
- TString * ptParam = TStatToolkit::FitPlane(chain,"Tr0.fP[4]+Tr1.fP[4]", strFit.Data(),"1", chi2,npoints,fitParam,covMatrix)
-
-
-
-*/
-
-
-/*
-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;
-
-
-
-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++";
-
-
-
-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();
-
-
-
-
-*/
-
-
-
-
//
void FindPairs(AliESDEvent *event);
Bool_t IsPair(AliExternalTrackParam *tr0, AliExternalTrackParam *tr1);
- void SetGainMap(AliTPCCalPad *GainMap){fGainMap = GainMap;};
static void CalculateBetheParams(TH2F *hist, Double_t * initialParam);
static Double_t CalculateMIPvalue(TH1F * hist);
- AliExternalTrackParam *Invert(AliExternalTrackParam *input);
AliExternalTrackParam *MakeTrack(const AliExternalTrackParam *track0, const AliExternalTrackParam *track1);
+ AliExternalTrackParam *MakeCombinedTrack(const AliExternalTrackParam *track0, const AliExternalTrackParam *track1);
+
void UpdateTrack(AliExternalTrackParam &track0, const AliExternalTrackParam &track1);
+ void SetComparison(AliExternalComparison * comp) { fComp=comp;}
+ AliExternalComparison *GetComparison() { return fComp;}
+
//
TH1F * GetHistNTracks(){return fHistNTracks;};
TH1F * GetHistClusters(){return fClusters;};
void FillAcordeHist(AliESDtrack *upperTrack);
- AliTPCCalPad *fGainMap; // gain map from Krypton calibration
+ AliExternalComparison * fComp; // comparison histogram
TH1F *fHistNTracks; // histogram showing number of ESD tracks per event
TH1F *fClusters; // histogram showing the number of clusters per track
TH2F *fModules; // 2d histogram of tracks which are propagated to the ACORDE scintillator array
// amplitude
sprintf(chname,"Amp_Sector%d",i);
- his1 = new TH1F(chname,chname,250,0,500); // valgrind
+ his1 = new TH1F(chname,chname,100,0,32); // valgrind
his1->SetXTitle("Max Amplitude (ADC)");
fArrayAmp->AddAt(his1,i);
sprintf(chname,"Amp_Sector%d",i+36);
- his1 = new TH1F(chname,chname,200,0,600); // valgrind 3 13,408,208 bytes in 229 blocks are still reachable
+ his1 = new TH1F(chname,chname,100,0,32); // valgrind 3 13,408,208 bytes in 229 blocks are still reachable
his1->SetXTitle("Max Amplitude (ADC)");
fArrayAmp->AddAt(his1,i+36);
// driftlength
sprintf(chname, "driftlengt vs. charge, ROC %i", i);
- prof1 = new TProfile(chname, chname, 500, 0, 250);
+ prof1 = new TProfile(chname, chname, 25, 0, 250);
prof1->SetYTitle("Charge");
prof1->SetXTitle("Driftlength");
fArrayChargeVsDriftlength->AddAt(prof1,i);
sprintf(chname, "driftlengt vs. charge, ROC %i", i+36);
- prof1 = new TProfile(chname, chname, 500, 0, 250);
+ prof1 = new TProfile(chname, chname, 25, 0, 250);
prof1->SetYTitle("Charge");
prof1->SetXTitle("Driftlength");
fArrayChargeVsDriftlength->AddAt(prof1,i+36);
TProfile *profAmpRow = (TProfile*)fArrayAmpRow->At(sector);
if ( irow >= kFirstLargePad) max /= kLargePadSize;
- profAmpRow->Fill( (Double_t)cluster0->GetRow(), max );
+ Double_t smax = TMath::Sqrt(max);
+ profAmpRow->Fill( (Double_t)cluster0->GetRow(), smax );
TH1F *hisAmp = (TH1F*)fArrayAmp->At(sector);
- hisAmp->Fill(max);
+ hisAmp->Fill(smax);
// remove the following two lines one day:
TProfile *profDriftLength = (TProfile*)fArrayChargeVsDriftlength->At(sector);
- profDriftLength->Fill( 250.-(Double_t)TMath::Abs(cluster0->GetZ()), max );
+ profDriftLength->Fill( 250.-(Double_t)TMath::Abs(cluster0->GetZ()), smax );
TProfile *profDriftLengthTmp = (TProfile*)(fcalPadRegionChargeVsDriftlength->GetObject(isegment, padSize));
- profDriftLengthTmp->Fill( 250.-(Double_t)TMath::Abs(cluster0->GetZ()), max );
+ profDriftLengthTmp->Fill( 250.-(Double_t)TMath::Abs(cluster0->GetZ()), smax );
fHclusterPerPadrow->Fill(irow); // fill histogram showing clusters per padrow
Int_t ipad = TMath::Nint(cluster0->GetPad());