// @(#) $Id$
+// Original: AliHLTConfMapFit.cxx,v 1.14 2005/06/14 10:55:21 cvetan
+
+/**************************************************************************
+ * This file is property of and copyright by the ALICE HLT Project *
+ * ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * 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 *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+
+/** @file AliHLTTPCConfMapFit.cxx
+ @author Anders Vestbo, maintained by Matthias Richter
+ @date
+ @brief Fit class for conformal mapping tracking.
+*/
-// Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
-//*-- Copyright © ALICE HLT Group
+// 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 "AliHLTTPCTransform.h"
#include "AliHLTTPCConfMapFit.h"
-/** \class AliHLTTPCConfMapFit
-<pre>
-//_____________________________________________________________
-// AliHLTTPCConfMapFit
-//
-// Fit class for conformal mapping tracking
-</pre>
-*/
-
#if __GNUC__ >= 3
using namespace std;
#endif
-ClassImp(AliHLTTPCConfMapFit)
+ClassImp(AliHLTTPCConfMapFit);
+AliHLTTPCConfMapFit::AliHLTTPCConfMapFit()
+ :
+ fTrack(NULL),
+ fVertex(NULL)
+{
+ //constructor
+}
AliHLTTPCConfMapFit::AliHLTTPCConfMapFit(AliHLTTPCConfMapTrack *track,AliHLTTPCVertex *vertex)
+ :
+ fTrack(track),
+ fVertex(vertex)
+
{
//constructor
- fTrack = track;
- fVertex = vertex;
+}
+
+AliHLTTPCConfMapFit::~AliHLTTPCConfMapFit()
+{
+ // destructor
}
Int_t AliHLTTPCConfMapFit::FitHelix()
return 0;
}
+Int_t AliHLTTPCConfMapFit::FitStraightLine() {
+ //fit the straight line
+ if(FitLineXY()) {
+ LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapFit::FitStraightLine","TrackFit")<<AliHLTTPCLog::kDec<<
+ "Problems during stright line fit in XY plane"<<ENDLOG;
+ return 1;
+ }
+ if(FitLineSZ()){
+ LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapFit::FitStraightLine","TrackFit")<<AliHLTTPCLog::kDec<<
+ "Problems during stright line fit in SZ plane"<<ENDLOG;
+ return 1;
+ }
+ return 0;
+}
+
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 !
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);
return 0 ;
}
+
+
+//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+// Straight Line Fit in x-y plane
+//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+Int_t AliHLTTPCConfMapFit::FitLineXY ( ){
+ // -----------------------------------------------------------------------------
+ // Implementation after Numerical Recipes in C, 2nd Edtion, chapter 15.2, p. 661
+ // with y = b*x + a
+ // and Data Analysis for Physical Science Students, Luis Lyons, chapter 2.4 p.51
+ // with y = a' + bx' , x' = x - <x>
+ // -----------------------------------------------------------------------------
+
+ 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 chi2 = 0.;
+
+ Int_t numOfHits = fTrack->GetNumberOfPoints();
+
+ Int_t co=0;
+
+ // - Loop over hits calculating average : xav
+ for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit()) {
+ co++;
+ AliHLTTPCConfMapPoint *currentHit = (AliHLTTPCConfMapPoint*)fTrack->GetCurrentHit();
+ // ** 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();
+ }
+
+ if(co!=numOfHits)
+ LOG(AliHLTTPCLog::kError,"AliHLTTPCConfMapFit::FitLineXY","TrackFit") << "Mismatch of hits. Counter: "<<co<<" nHits: "<<numOfHits<<ENDLOG;
+
+ 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();
+ }
+
+ Double_t det = sPrime*sxxPrime + sxPrime*sxPrime;
+
+ if (fabs(det) < 1e-20) {
+ LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapFit::FitLineXY","TrackFit") << "Determinant == 0" << ENDLOG;
+ chi2 = 99999.F ;
+ fTrack->SetChiSq1(chi2);
+ 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 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) ;
+
+ // Calculate chi2
+ for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit()) {
+ AliHLTTPCConfMapPoint *currentHit = (AliHLTTPCConfMapPoint*)fTrack->GetCurrentHit();
+ Double_t tempchi = currentHit->GetY() - aPrime - b*(currentHit->GetX() - xav);
+ chi2 += tempchi*tempchi*currentHit->GetXYWeight() ;
+ }
+
+ Double_t a = aPrime - b*xav;
+
+
+ // Set TrackParameter
+ fTrack->SetChiSq1(chi2);
+ fTrack->SetPsi(psi);
+ fTrack->SetPsierr(sigma2b);
+ fTrack->SetCenterX(0.); // Set to point on the track (for UpdateToFirstPoint)
+ fTrack->SetCenterY(a); // Set to point on the track (for UpdateToFirstPoint)
+
+ //Set the first point on the track to the space point coordinates of the innermost track
+ //This will be updated to lie on the fit later on (AliHLTTPCTrack::UpdateToFirstPoint).
+ AliHLTTPCConfMapPoint *lastHit = (AliHLTTPCConfMapPoint*)fTrack->GetLastHit();
+ Double_t x0 = lastHit->GetX();
+ Double_t y0 = lastHit->GetY();
+ fTrack->SetFirstPoint(x0,y0,0); //Z-value is set in FitLineSZ
+
+
+ //Set Defaults
+ fTrack->SetRadius(-1.);
+ fTrack->SetCharge(1);
+ fTrack->SetPt(-1.);
+
+
+ return 0;
+}
+
+
+//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+// Straight Line Fit in s-z plane
+//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+Int_t AliHLTTPCConfMapFit::FitLineSZ ( ){
+ // -----------------------------------------------------------------------------
+ // Implementation after Numerical Recipes in C, 2nd Edtion, chapter 15.2, p. 661
+ // with z = b*s + a
+ // and Data Analysis for Physical Science Students, Luis Lyons, chapter 2.4 p.51
+ // with z = a' + bs' , s' = s - <s>
+ // -----------------------------------------------------------------------------
+
+ 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 chi2 = 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;
+
+ // - Loop over hits calculating length in xy-plane: s
+ // - Loop over hits calculating average : sav
+ for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit()) {
+ AliHLTTPCConfMapPoint *currentHit = (AliHLTTPCConfMapPoint*)fTrack->GetCurrentHit();
+ if(currentHit != fTrack->GetFirstHit()) {
+ Double_t dx = currentHit->GetX() - previousHit->GetX() ;
+ Double_t dy = currentHit->GetY() - previousHit->GetY() ;
+ 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();
+ slength = (Double_t)sqrt ( dx*dx + dy*dy );
+ }
+
+ currentHit->SetS(slength);
+
+ S += currentHit->GetZWeight();
+ Ss += currentHit->GetZWeight() * currentHit->GetS();
+ }
+
+ Double_t sav = (Double_t)Ss / S;
+
+ // Calculate weighted means
+ 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();
+ }
+
+ Double_t det = sPrime*sssPrime + ssPrime*ssPrime;
+
+ if (fabs(det) < 1e-20) {
+ LOG(AliHLTTPCLog::kDebug,"AliHLTTPCConfMapFit::FitLineSZ","TrackFit") << "Determinant == 0" << ENDLOG;
+ chi2 = 99999.F ;
+ fTrack->SetChiSq2(chi2);
+ 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 a = aPrime - b*sav;
+
+ Double_t sigma2b = (Double_t) 1. / sssPrime;
+ Double_t sigma2aprime = (Double_t) 1. /sPrime;
+
+ Double_t sigma2a = sigma2aprime + sav*sav * sigma2b*sigma2b;
+
+ // Calculate chi2
+ for(fTrack->StartLoop(); fTrack->LoopDone(); fTrack->GetNextHit()) {
+ AliHLTTPCConfMapPoint *currentHit = (AliHLTTPCConfMapPoint*)fTrack->GetCurrentHit();
+ Double_t tempchi = currentHit->GetZ() - aPrime - b*(currentHit->GetS() - sav);
+ chi2 += tempchi*tempchi*currentHit->GetZWeight() ;
+ }
+
+ // Set TrackParameter
+ fTrack->SetChiSq2(chi2);
+ fTrack->SetTgl(b);
+ fTrack->SetZ0(a);
+ fTrack->SetTglerr(sigma2b);
+// fTrack->SetZ0err(sigma2aprime); // maybe subject to check
+ fTrack->SetZ0err(sigma2a); // maybe subject to check
+ return 0;
+}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff