// Original: AliHLTConfMapFit.cxx,v 1.14 2005/06/14 10:55:21 cvetan
/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * This file is property of and copyright by the ALICE HLT Project *
+ * ALICE Experiment at CERN, All rights reserved. *
* *
- * Authors: Anders Vestbo *
- * Matthias Richter <Matthias.Richter@ift.uib.no> *
- * for The ALICE Off-line Project. *
+ * Primary Authors: Anders Vestbo, maintained by *
+ * Matthias Richter <Matthias.Richter@ift.uib.no> *
+ * for The ALICE HLT Project. *
* *
* Permission to use, copy, modify and distribute this software and its *
* documentation strictly for non-commercial purposes is hereby granted *
@brief Fit class for conformal mapping tracking.
*/
+// see header file for class documentation //
+// or //
+// refer to README to build package //
+// or //
+// visit http://web.ift.uib.no/~kjeks/doc/alice-hlt //
+
#include "AliHLTTPCRootTypes.h"
#include "AliHLTTPCLogging.h"
#include "AliHLTTPCVertex.h"
using namespace std;
#endif
-ClassImp(AliHLTTPCConfMapFit)
-
+ClassImp(AliHLTTPCConfMapFit);
AliHLTTPCConfMapFit::AliHLTTPCConfMapFit()
:
//constructor
}
-AliHLTTPCConfMapFit::AliHLTTPCConfMapFit(const AliHLTTPCConfMapFit&)
- :
- fTrack(NULL),
- fVertex(NULL)
-{
- // dummy copy constructor
- //HLTFatal("copy constructor untested");
-}
-
-AliHLTTPCConfMapFit& AliHLTTPCConfMapFit::operator=(const AliHLTTPCConfMapFit&)
-{
- // dummy assignment operator
- //HLTFatal("assignment operator untested");
- return *this;
-}
-
AliHLTTPCConfMapFit::~AliHLTTPCConfMapFit()
{
// destructor
return 0;
}
-// #### -B0-CHANGE-START == JMT
-
Int_t AliHLTTPCConfMapFit::FitStraightLine() {
//fit the straight line
if(FitLineXY()) {
return 0;
}
-// #### -B0-CHANGE-END == JMT
-
Int_t AliHLTTPCConfMapFit::FitCircle()
{
//-----------------------------------------------------------------
Double_t xav = 0.0 ;
Double_t yav = 0.0 ;
- Int_t num_of_hits = fTrack->GetNumberOfPoints();
+ Int_t numOfHits = fTrack->GetNumberOfPoints();
//
// Loop over hits calculating average
Int_t co=0;
xav += cHit->GetXYWeight() * cHit->GetX() ;
yav += cHit->GetXYWeight() * cHit->GetY() ;
}
- if(co!=num_of_hits)
+ if(co!=numOfHits)
LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapFit::FitCircle","TrackFit")<<AliHLTTPCLog::kDec<<
- "Mismatch of hits. Counter: "<<co<<" nHits: "<<num_of_hits<<ENDLOG;
+ "Mismatch of hits. Counter: "<<co<<" nHits: "<<numOfHits<<ENDLOG;
if (fTrack->ComesFromMainVertex() == true)
{
wsum += fVertex->GetXYWeight() ;
Double_t xixi, yiyi, riri, wiriri, xold, yold ;
- //for (hit_counter=0; hit_counter<num_of_hits; hit_counter++)
+ //for (hit_counter=0; hit_counter<numOfHits; hit_counter++)
for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit())
{
//AliHLTTPCConfMapPoint *cHit = (AliHLTTPCConfMapPoint*)hits->At(hit_counter);
rrrrav = rrrrav / wsum ;
xyav = xyav / wsum ;
- Int_t const ntry = 5 ;
+ const Int_t ntry = 5 ;
//
//--> USE THESE TO GET THE COEFFICIENTS OF THE 4-TH ORDER POLYNIMIAL
//--> DON'T PANIC - THE THIRD ORDER TERM IS ZERO !
fTrack->SetCenterX(acent);
fTrack->SetCenterY(bcent);
+ fTrack->SetY0err(0.03);
+ //set error for pT and Y. psi, Z and Tgl are set.
+
//
// Get errors from fast fit
//
Double_t radius = (Double_t)(fTrack->GetPt() / AliHLTTPCTransform::GetBFieldValue() ) ;
//TObjArray *hits = fTrack->GetHits();
- //Int_t num_of_hits = fTrack->GetNumberOfPoints();
+ //Int_t numOfHits = fTrack->GetNumberOfPoints();
if (0)// fTrack->ComesFromMainVertex() == true )
{
}
Double_t localPsi = 0.5F * sqrt ( dx*dx + dy*dy ) / radius ;
- Double_t total_s ;
+ Double_t totalS ;
if ( fabs(localPsi) < 1. )
{
- total_s = 2.0 * radius * asin ( localPsi ) ;
+ totalS = 2.0 * radius * asin ( localPsi ) ;
}
else
{
- total_s = 2.0 * radius * AliHLTTPCTransform::Pi() ;
+ totalS = 2.0 * radius * AliHLTTPCTransform::Pi() ;
}
AliHLTTPCConfMapPoint *previousHit = NULL;
//for ( startLoop() ; done() ; nextHit() ) {
Double_t dpsi,s;
- // for(hit_counter=0; hit_counter<num_of_hits; hit_counter++)
+ // for(hit_counter=0; hit_counter<numOfHits; hit_counter++)
for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit())
{
// AliHLTTPCConfMapPoint *cHit = (AliHLTTPCConfMapPoint*)hits->At(hit_counter);
cHit->SetS(s);
}
else
- cHit->SetS(total_s);
- // cHit->s = total_s ;
+ cHit->SetS(totalS);
+ // cHit->s = totalS ;
sum += cHit->GetZWeight() ;
ss += cHit->GetZWeight() * cHit->GetS() ;
chi2 = 0.;
Double_t r1 ;
- //for(hit_counter=0; hit_counter<num_of_hits; hit_counter++)
+ //for(hit_counter=0; hit_counter<numOfHits; hit_counter++)
for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit())
{
//AliHLTTPCConfMapPoint *cHit = (AliHLTTPCConfMapPoint*)hits->At(hit_counter);
fTrack->SetTglerr(dtanl);
fTrack->SetZ0err(dz0);
+
+ //The calculation of pT comes from "Perticle Physics Booklet":
+ //24.11 Measurement of particle momenta in a uniform magnetic field.
+ //(dk)^2=(dk_res)^2 + (dk_ms)^2
+ //for now k_ms is 0. Need to find length of track in 3D.
+
+ Double_t lengthXY = fTrack->GetLengthXY();
+ //Double_t lengthTot = fTrack->GetLengthTot();
+ //Double_t beta = fTrack->GetP()/TMath::Sqrt((fTrack->GetP()*fTrack->GetP())+(0.13957*0.13957));
+ //Double_t lambda = TMath::ATan(fTrack->GetTgl());
+ Double_t lengthXY2 = lengthXY*lengthXY;
+ Int_t nCluster4 = fTrack->GetNHits()+4;
+
+ Double_t Kres = 0.03/lengthXY2;
+ Kres = Kres * TMath::Sqrt(720/nCluster4);
+
+ //Double_t Kres = (0.03/(lengthXY*lengthXY))*TMath::Sqrt(720/(fTrack->GetNHits()+4));
+
+ //Double_t d = lengthTot*fTrack->GetP()*beta*TMath::Cos(lambda)*TMath::Cos(lambda);
+ //Double_t Kms = (0.016/d)*TMath::Sqrt(lengthTot/24.0);
+ Double_t Kms = 0.0;
+
+ Double_t KTot = TMath::Sqrt((Kres * Kres) + (Kms * Kms));
+ Double_t Pterr = (1/(0.3*AliHLTTPCTransform::GetBField()))*KTot;
+
+ fTrack->SetPterr(Pterr);
+
return 0 ;
}
-// #### -B0-CHANGE-START == JMT
-
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// Straight Line Fit in x-y plane
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// with y = a' + bx' , x' = x - <x>
// -----------------------------------------------------------------------------
- Double_t S = 0.;
- Double_t Sx = 0.;
+ Double_t s = 0.;
+ Double_t sx = 0.;
- Double_t SPrime = 0.;
- Double_t SxPrime = 0.;
- Double_t SxxPrime = 0.;
- Double_t SyPrime = 0.;
- Double_t SxyPrime = 0.;
+ Double_t sPrime = 0.;
+ Double_t sxPrime = 0.;
+ Double_t sxxPrime = 0.;
+ Double_t syPrime = 0.;
+ Double_t sxyPrime = 0.;
Double_t chi2 = 0.;
- Int_t num_of_hits = fTrack->GetNumberOfPoints();
+ Int_t numOfHits = fTrack->GetNumberOfPoints();
Int_t co=0;
// ** maybe not necessary, already done in ConfMapPoint
currentHit->SetXYWeight( 1./ (Double_t)(currentHit->GetXerr()*currentHit->GetXerr() + currentHit->GetYerr()*currentHit->GetYerr()) );
// **
- S += currentHit->GetXYWeight();
- Sx += currentHit->GetXYWeight() * currentHit->GetX();
+ s += currentHit->GetXYWeight();
+ sx += currentHit->GetXYWeight() * currentHit->GetX();
}
- if(co!=num_of_hits)
- LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapFit::FitLineXY","TrackFit") << "Mismatch of hits. Counter: "<<co<<" nHits: "<<num_of_hits<<ENDLOG;
+ if(co!=numOfHits)
+ LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapFit::FitLineXY","TrackFit") << "Mismatch of hits. Counter: "<<co<<" nHits: "<<numOfHits<<ENDLOG;
- Double_t xav = (Double_t)Sx / S;
+ Double_t xav = (Double_t)sx / s;
// Calculate weighted means
for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit()) {
AliHLTTPCConfMapPoint *currentHit = (AliHLTTPCConfMapPoint*)fTrack->GetCurrentHit();
Double_t xPrime = currentHit->GetX() - xav;
- SPrime += currentHit->GetXYWeight();
- SxPrime += currentHit->GetXYWeight() * xPrime;
- SxxPrime += currentHit->GetXYWeight() * xPrime * xPrime;
- SyPrime += currentHit->GetXYWeight() * currentHit->GetY();
- SxyPrime += currentHit->GetXYWeight() * xPrime * currentHit->GetY();
+ sPrime += currentHit->GetXYWeight();
+ sxPrime += currentHit->GetXYWeight() * xPrime;
+ sxxPrime += currentHit->GetXYWeight() * xPrime * xPrime;
+ syPrime += currentHit->GetXYWeight() * currentHit->GetY();
+ sxyPrime += currentHit->GetXYWeight() * xPrime * currentHit->GetY();
}
- Double_t det = SPrime*SxxPrime + SxPrime*SxPrime;
+ Double_t det = sPrime*sxxPrime + sxPrime*sxPrime;
if (fabs(det) < 1e-20) {
LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapFit::FitLineXY","TrackFit") << "Determinant == 0" << ENDLOG;
return -1 ;
}
- Double_t b = (Double_t)(SPrime*SxyPrime - SxPrime*SyPrime) / det; // line parameter b
- Double_t aPrime = (Double_t)(SxxPrime*SyPrime - SxPrime*SxyPrime) / det; // line parameter a
+ Double_t b = (Double_t)(sPrime*sxyPrime - sxPrime*syPrime) / det; // line parameter b
+ Double_t aPrime = (Double_t)(sxxPrime*syPrime - sxPrime*sxyPrime) / det; // line parameter a
- Double_t sigma2b = (Double_t)1. / SxxPrime;
- //-- Double_t sigma2aprime = (Double_t)1. /SPrime;
+ Double_t sigma2b = (Double_t)1. / sxxPrime;
+ //-- Double_t sigma2aprime = (Double_t)1. /sPrime;
// Get gradient angle psi of line in xy plane
Double_t psi = (Double_t) atan(b) ;
Double_t S = 0.;
Double_t Ss = 0.;
- Double_t SPrime = 0.;
- Double_t SsPrime = 0.;
- Double_t SssPrime = 0.;
- Double_t SzPrime = 0.;
- Double_t SszPrime = 0.;
+ Double_t sPrime = 0.;
+ Double_t ssPrime = 0.;
+ Double_t sssPrime = 0.;
+ Double_t szPrime = 0.;
+ Double_t sszPrime = 0.;
Double_t chi2 = 0.;
- Double_t s = 0.;
+
+ // Matthias 16.10.2007
+ // what's that!!! local variables 's' and 'S'
+ // change Double_t s = 0.; -> slength
+ Double_t slength = 0.;
AliHLTTPCConfMapPoint *previousHit = NULL;
if(currentHit != fTrack->GetFirstHit()) {
Double_t dx = currentHit->GetX() - previousHit->GetX() ;
Double_t dy = currentHit->GetY() - previousHit->GetY() ;
- s = previousHit->GetS() - (Double_t)sqrt ( dx*dx + dy*dy );
+ slength = previousHit->GetS() - (Double_t)sqrt ( dx*dx + dy*dy );
}
else{
Double_t dx = ((AliHLTTPCConfMapPoint *)fTrack->GetFirstHit())->GetX() - ((AliHLTTPCConfMapPoint *)fTrack->GetLastHit())->GetX();
Double_t dy = ((AliHLTTPCConfMapPoint *)fTrack->GetFirstHit())->GetY() - ((AliHLTTPCConfMapPoint *)fTrack->GetLastHit())->GetY();
- s = (Double_t)sqrt ( dx*dx + dy*dy );
+ slength = (Double_t)sqrt ( dx*dx + dy*dy );
}
- currentHit->SetS(s);
+ currentHit->SetS(slength);
S += currentHit->GetZWeight();
Ss += currentHit->GetZWeight() * currentHit->GetS();
for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit()) {
AliHLTTPCConfMapPoint *currentHit = (AliHLTTPCConfMapPoint*)fTrack->GetCurrentHit();
- Double_t sPrime = currentHit->GetS() - sav;
- SPrime += currentHit->GetZWeight();
- SsPrime += currentHit->GetZWeight() * sPrime;
- SssPrime += currentHit->GetZWeight() * sPrime * sPrime;
- SzPrime += currentHit->GetZWeight() * currentHit->GetZ();
- SszPrime += currentHit->GetZWeight() * sPrime * currentHit->GetZ();
+ // Matthias 20.05.2008
+ // here was a shadowed variable, sPrime is formerly defined
+ // renamed it to lsPrime ('local')
+ Double_t lsPrime = currentHit->GetS() - sav;
+ lsPrime += currentHit->GetZWeight();
+ ssPrime += currentHit->GetZWeight() * lsPrime;
+ sssPrime += currentHit->GetZWeight() * lsPrime * lsPrime;
+ szPrime += currentHit->GetZWeight() * currentHit->GetZ();
+ sszPrime += currentHit->GetZWeight() * lsPrime * currentHit->GetZ();
}
- Double_t det = SPrime*SssPrime + SsPrime*SsPrime;
+ Double_t det = sPrime*sssPrime + ssPrime*ssPrime;
if (fabs(det) < 1e-20) {
LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapFit::FitLineSZ","TrackFit") << "Determinant == 0" << ENDLOG;
return -1 ;
}
- Double_t b = (Double_t)(SPrime*SszPrime - SsPrime*SzPrime) / det; // line parameter b
- Double_t aPrime = (Double_t)(SssPrime*SzPrime - SsPrime*SszPrime) / det; // line parameter a
+ Double_t b = (Double_t)(sPrime*sszPrime - ssPrime*szPrime) / det; // line parameter b
+ Double_t aPrime = (Double_t)(sssPrime*szPrime - ssPrime*sszPrime) / det; // line parameter a
Double_t a = aPrime - b*sav;
- Double_t sigma2b = (Double_t) 1. / SssPrime;
- Double_t sigma2aprime = (Double_t) 1. /SPrime;
+ Double_t sigma2b = (Double_t) 1. / sssPrime;
+ Double_t sigma2aprime = (Double_t) 1. /sPrime;
Double_t sigma2a = sigma2aprime + sav*sav * sigma2b*sigma2b;
return 0;
}
-// #### -B0-CHANGE-END == JMT