From 354cc456d9afe5a33e22846925e29859da14b0e5 Mon Sep 17 00:00:00 2001 From: richterm Date: Mon, 6 Jun 2011 07:11:27 +0000 Subject: [PATCH] calculation of helix parameters and crossing points; added overloaded TObject::Draw; reading tracks from multiple binary files supposted --- HLT/BASE/util/AliHLTGlobalBarrelTrack.cxx | 495 +++++++++++++++++++++- HLT/BASE/util/AliHLTGlobalBarrelTrack.h | 36 +- 2 files changed, 528 insertions(+), 3 deletions(-) diff --git a/HLT/BASE/util/AliHLTGlobalBarrelTrack.cxx b/HLT/BASE/util/AliHLTGlobalBarrelTrack.cxx index fb790aa1dd1..de463e49b88 100644 --- a/HLT/BASE/util/AliHLTGlobalBarrelTrack.cxx +++ b/HLT/BASE/util/AliHLTGlobalBarrelTrack.cxx @@ -30,10 +30,16 @@ #include #include +#include #include "AliHLTGlobalBarrelTrack.h" +#include "AliHLTSpacePointContainer.h" +#include "AliHLTMisc.h" #include "TClonesArray.h" #include "TFile.h" #include "TArrayC.h" +#include "TMath.h" +#include "TMarker.h" +#include "TArc.h" /** ROOT macro for the implementation of ROOT specific class methods */ ClassImp(AliHLTGlobalBarrelTrack) @@ -44,6 +50,10 @@ AliHLTGlobalBarrelTrack::AliHLTGlobalBarrelTrack() , fLastX(0.0) , fLastY(0.0) , fTrackID(-1) + , fHelixRadius(0.0) + , fHelixCenterX(0.0) + , fHelixCenterY(0.0) + , fSpacePoints(NULL) { // see header file for class documentation // or @@ -58,6 +68,10 @@ AliHLTGlobalBarrelTrack::AliHLTGlobalBarrelTrack(const AliHLTGlobalBarrelTrack& , fLastX(t.GetLastPointX()) , fLastY(t.GetLastPointY()) , fTrackID(t.TrackID()) + , fHelixRadius(t.fHelixRadius) + , fHelixCenterX(t.fHelixCenterX) + , fHelixCenterY(t.fHelixCenterY) + , fSpacePoints(NULL) { // see header file for class documentation fPoints.assign(t.fPoints.begin(), t.fPoints.end()); @@ -69,6 +83,10 @@ AliHLTGlobalBarrelTrack::AliHLTGlobalBarrelTrack(const AliHLTExternalTrackParam& , fLastX(p.fLastX) , fLastY(p.fLastY) , fTrackID(p.fTrackID) + , fHelixRadius(0.0) + , fHelixCenterX(0.0) + , fHelixCenterY(0.0) + , fSpacePoints(NULL) { // see header file for class documentation @@ -87,6 +105,10 @@ AliHLTGlobalBarrelTrack::AliHLTGlobalBarrelTrack(const AliExternalTrackParam& p , fLastX(0) , fLastY(0) , fTrackID(0) + , fHelixRadius(0.0) + , fHelixCenterX(0.0) + , fHelixCenterY(0.0) + , fSpacePoints(NULL) { // see header file for class documentation *(dynamic_cast(this))=p; @@ -184,6 +206,74 @@ int AliHLTGlobalBarrelTrack::SetPoints(const UInt_t* pArray, UInt_t arraySize) return fPoints.size(); } +int AliHLTGlobalBarrelTrack::CalculateHelixParams() +{ + // calculate radius and center of the helix + float bfield=AliHLTMisc::Instance().GetBz(); + if (TMath::Abs(bfield) straight lines + fHelixRadius=kVeryBig; + fHelixCenterX=kVeryBig; + fHelixCenterY=kVeryBig; + } else { + fHelixRadius = GetSignedPt()/(-kB2C*bfield); + Double_t trackPhi = Phi()-GetAlpha(); + + //cout << "Helix: phi=" << trackPhi << " x=" << GetX() << " y=" << GetY() << endl; + fHelixCenterX = GetX() + fHelixRadius * sin(trackPhi); + fHelixCenterY = GetY() - fHelixRadius * cos(trackPhi); + //cout << "Helix: center" << " x=" << fHelixCenterX << " y=" << fHelixCenterY << endl; + } + return 0; +} + +int AliHLTGlobalBarrelTrack::CalculateCrossingPoint(float xPlane, float alphaPlane, float& y, float& z) +{ + // calculate crossing point of helix with a plane in yz + // in the local coordinates of the plane + int iResult=0; + if (TMath::Abs(fHelixRadius)=kVeryBig) { + // no magnetic field -> straight lines + } else { + // rotate helix center to local coordinates of the plane reference frame + float cosa=TMath::Cos(alphaPlane-GetAlpha()); + float sina=TMath::Sin(alphaPlane-GetAlpha()); + float cx= fHelixCenterX * cosa + fHelixCenterY * sina; + float cy=-fHelixCenterX * sina + fHelixCenterY * cosa; + + // crossing point of helix with plane + Double_t aa = (cx - xPlane)*(cx - xPlane); + Double_t r2 = fHelixRadius*fHelixRadius; + if(aa > r2) // no crossing + return 0; + + Double_t aa2 = sqrt(r2 - aa); + Double_t y1 = cy + aa2; + Double_t y2 = cy - aa2; + y = y1; + if(TMath::Abs(y2) < TMath::Abs(y1)) y = y2; + + Double_t angle1 = atan2((y - cy),(xPlane - cx)); + if(angle1 < 0) angle1 += 2.*TMath::Pi(); + Double_t angle2 = atan2((GetY() - fHelixCenterY),(GetX() - fHelixCenterX)); + if(angle2 < 0) angle2 += TMath::TwoPi(); + + Double_t diffangle = angle1 - angle2; + diffangle = fmod(diffangle,TMath::TwoPi()); + if((GetSign()*diffangle) > 0) diffangle = diffangle - GetSign()*TMath::TwoPi(); + + Double_t stot = TMath::Abs(diffangle)*fHelixRadius; + + z = GetZ() + stot*GetTgl(); + } + return 1; +} + void AliHLTGlobalBarrelTrack::Print(Option_t* option) const { // see header file for class documentation @@ -198,6 +288,386 @@ void AliHLTGlobalBarrelTrack::Print(Option_t* option) const // cout << " Signed1Pt " << GetSigned1Pt() << endl; } +void AliHLTGlobalBarrelTrack::Draw(Option_t *option) +{ + /// Inherited from TObject, draw the track + float scale=250; + float center[2]={0.5,0.5}; + + if (TMath::Abs(fHelixRadius) tokens(strOption.Tokenize(" ")); + if (!tokens.get()) return; + for (int i=0; iGetEntriesFast(); i++) { + if (!tokens->At(i)) continue; + const char* key=""; + TString arg=tokens->At(i)->GetName(); + + key="scale="; + if (arg.BeginsWith(key)) { + arg.ReplaceAll(key, ""); + scale=arg.Atof(); + continue; + } + key="centerx="; + if (arg.BeginsWith(key)) { + arg.ReplaceAll(key, ""); + center[0]=arg.Atof(); + continue; + } + key="centery="; + if (arg.BeginsWith(key)) { + arg.ReplaceAll(key, ""); + center[1]=arg.Atof(); + continue; + } + key="spacepoints"; + if (arg.CompareTo(key)==0) { + if (fSpacePoints) DrawProjXYSpacePoints(option, fSpacePoints, scale, center); + continue; + } + key="trackarc"; + if (arg.CompareTo(key)==0) { + DrawProjXYTrack(option, scale, center); + continue; + } + } +} + +int AliHLTGlobalBarrelTrack::DrawProjXYSpacePoints(Option_t */*option*/, const AliHLTSpacePointContainer* spacePoints, const float scale, float center[2]) +{ + /// draw space points + int markerColor=3; + + if (!spacePoints) return -EINVAL; + + const UInt_t* pointids=GetPoints(); + for (unsigned i=0; iGetPhi(pointids[i]); + float cosphi=TMath::Cos(clusterphi); + float sinphi=TMath::Sin(clusterphi); + float clusterx=spacePoints->GetX(pointids[i]); + float clustery=spacePoints->GetY(pointids[i]); + // rotate + float pointx= clusterx*sinphi + clustery*cosphi; + float pointy=-clusterx*cosphi + clustery*sinphi; + + // FIXME: cleanup of marker objects + TMarker* m=new TMarker(pointx/(2*scale)+center[0], pointy/(2*scale)+center[1], 3); + m->SetMarkerColor(markerColor); + m->Draw("same"); + } + return 0; +} + +// FIXME: make this a general geometry definition +// through an abstract class interface +const Double_t gkTPCX[159] = { + 85.195, + 85.945, + 86.695, + 87.445, + 88.195, + 88.945, + 89.695, + 90.445, + 91.195, + 91.945, + 92.695, + 93.445, + 94.195, + 94.945, + 95.695, + 96.445, + 97.195, + 97.945, + 98.695, + 99.445, + 100.195, + 100.945, + 101.695, + 102.445, + 103.195, + 103.945, + 104.695, + 105.445, + 106.195, + 106.945, + 107.695, + 108.445, + 109.195, + 109.945, + 110.695, + 111.445, + 112.195, + 112.945, + 113.695, + 114.445, + 115.195, + 115.945, + 116.695, + 117.445, + 118.195, + 118.945, + 119.695, + 120.445, + 121.195, + 121.945, + 122.695, + 123.445, + 124.195, + 124.945, + 125.695, + 126.445, + 127.195, + 127.945, + 128.695, + 129.445, + 130.195, + 130.945, + 131.695, + 135.180, + 136.180, + 137.180, + 138.180, + 139.180, + 140.180, + 141.180, + 142.180, + 143.180, + 144.180, + 145.180, + 146.180, + 147.180, + 148.180, + 149.180, + 150.180, + 151.180, + 152.180, + 153.180, + 154.180, + 155.180, + 156.180, + 157.180, + 158.180, + 159.180, + 160.180, + 161.180, + 162.180, + 163.180, + 164.180, + 165.180, + 166.180, + 167.180, + 168.180, + 169.180, + 170.180, + 171.180, + 172.180, + 173.180, + 174.180, + 175.180, + 176.180, + 177.180, + 178.180, + 179.180, + 180.180, + 181.180, + 182.180, + 183.180, + 184.180, + 185.180, + 186.180, + 187.180, + 188.180, + 189.180, + 190.180, + 191.180, + 192.180, + 193.180, + 194.180, + 195.180, + 196.180, + 197.180, + 198.180, + 199.430, + 200.930, + 202.430, + 203.930, + 205.430, + 206.930, + 208.430, + 209.930, + 211.430, + 212.930, + 214.430, + 215.930, + 217.430, + 218.930, + 220.430, + 221.930, + 223.430, + 224.930, + 226.430, + 227.930, + 229.430, + 230.930, + 232.430, + 233.930, + 235.430, + 236.930, + 238.430, + 239.930, + 241.430, + 242.930, + 244.430, + 245.930 +}; + +int AliHLTGlobalBarrelTrack::DrawProjXYTrack(Option_t *option, const float scale, float center[2]) +{ + /// draw track + bool bDrawArc=false; // draw TArc + bool bNoTrackPoints=false; // don't draw track points + TString strOption(option); + if (strOption.IsNull()) strOption="spacepoints trackarc"; + std::auto_ptr tokens(strOption.Tokenize(" ")); + if (!tokens.get()) return 0; + for (int i=0; iGetEntriesFast(); i++) { + if (!tokens->At(i)) continue; + const char* key=""; + TString arg=tokens->At(i)->GetName(); + + key="drawarc"; + if (arg.BeginsWith(key)) { + bDrawArc=true; + continue; + } + key="notrackpoints"; + if (arg.BeginsWith(key)) { + bNoTrackPoints=true; + continue; + } + } + + float cosa=TMath::Cos(GetAlpha()); + float sina=TMath::Sin(GetAlpha()); + + // first point + float firstpoint[2]; + firstpoint[0]= GetX()*sina + GetY()*cosa; + firstpoint[1]=-GetX()*cosa + GetY()*sina; + { + //cout << " first point alpha=" << GetAlpha() << " x: " << firstpoint[0] << " y: " << firstpoint[1] << endl; + // FIXME: cleanup of marker objects + TMarker* m=new TMarker(firstpoint[0]/(2*scale)+center[0], firstpoint[1]/(2*scale)+center[1], 29); + m->SetMarkerSize(2); + m->SetMarkerColor(2); + m->Draw("same"); + } + + // draw points in step width and remember the last point + float lastpoint[2]={0.0, 0.0}; + int firstpadrow=0; + for (; firstpadrow<159 && gkTPCX[firstpadrow]=kVeryBig) { + // no magnetic field -> straight lines + } else { + // rotate helix center to local coordinates of the plane reference frame + float cx= fHelixCenterX * sina + fHelixCenterY * cosa; + float cy=-fHelixCenterX * cosa + fHelixCenterY * sina; + + float diffx=cx-firstpoint[0]; + float diffy=cy-firstpoint[1]; + float phimin=0.0; + float phimax=2*TMath::Pi(); + if (TMath::Abs(diffx)0) phimin+=TMath::Pi(); + + diffx=cx-lastpoint[0]; + diffy=cy-lastpoint[1]; + if (TMath::Abs(diffx)=0 && padrow<159; padrow+=step) { + float x=gkTPCX[padrow]; + float y=0.0; + float z=0.0; + + int maxshift=9; + int shift=0; + int result=0; + do { + if ((result=CalculateCrossingPoint(x, GetAlpha()-offsetAlpha, y, z))<1) break; + float pointAlpha=TMath::ATan(y/x); + if (TMath::Abs(pointAlpha)>TMath::Pi()/18) { + offsetAlpha+=(pointAlpha>0?-1:1)*TMath::Pi()/9; + result=0; + markerColor++; + cosa=TMath::Cos(GetAlpha()-offsetAlpha); + sina=TMath::Sin(GetAlpha()-offsetAlpha); + } + } while (result==0 && shift++ tracks; iResult=ConvertTrackDataArray(pTracks, buffer->GetSize(), tracks); if (iResult>=0) { - tgt.ExpandCreate(tracks.size()); + int offset=tgt.GetEntriesFast(); + tgt.ExpandCreate(offset+tracks.size()); for (unsigned i=0; i effC++ warning, + AliHLTGlobalBarrelTrack& operator=(const AliHLTGlobalBarrelTrack& t) { + if (this==&t) return *this; + this->~AliHLTGlobalBarrelTrack(); new (this) AliHLTGlobalBarrelTrack(t); + return *this; + } template AliHLTGlobalBarrelTrack& operator=(const c& t) { this->~AliHLTGlobalBarrelTrack(); new (this) AliHLTGlobalBarrelTrack(t); @@ -75,6 +84,9 @@ class AliHLTGlobalBarrelTrack : public AliKalmanTrack /// Get the list of associated points const UInt_t* GetPoints() const; + /// Set the space point data + void SetSpacePointContainer(AliHLTSpacePointContainer* points) {fSpacePoints=points;} + /// Set the list of associated points int SetPoints(const UInt_t* pArray, UInt_t arraySize); @@ -92,11 +104,25 @@ class AliHLTGlobalBarrelTrack : public AliKalmanTrack /// Inherited from TObject, prints the track parameters virtual void Print(Option_t* option = "") const; + /// Inherited from TObject, draw the track + virtual void Draw(Option_t *option=""); + + int DrawProjXYSpacePoints(Option_t *option, const AliHLTSpacePointContainer* fSpacePoints, const float scale, float center[2]); + int DrawProjXYTrack(Option_t *option, const float scale, float center[2]); + int DrawProjXYTrackPoints(Option_t *option, const float scale, const float center[2], int firstpadrow, int step, float lastpoint[2]); + Double_t GetPathLengthTo( Double_t x, Double_t b ) const; static int ReadTracks(const char* filename, TClonesArray& tgt, AliHLTComponentDataType dt=kAliHLTVoidDataType, unsigned specification=kAliHLTVoidDataSpec); + static int ReadTrackList(const char* listfile, TClonesArray& tgt, AliHLTComponentDataType dt=kAliHLTVoidDataType, unsigned specification=kAliHLTVoidDataSpec); protected: + /// calculate and set internal helix parameters + int CalculateHelixParams(); + + /// calculate crossing point with a plane parallel to z axis, at distance x + /// and phi + int CalculateCrossingPoint(float xPlane, float phiPlane, float& u, float& v); private: @@ -111,6 +137,14 @@ class AliHLTGlobalBarrelTrack : public AliKalmanTrack /// track Id for identification during the reconstruction Int_t fTrackID; // + /// helix parameters + Float_t fHelixRadius; // + Float_t fHelixCenterX; // + Float_t fHelixCenterY; // + + /// the space points assigned to the track + AliHLTSpacePointContainer* fSpacePoints; //! + ClassDef(AliHLTGlobalBarrelTrack, 0) }; #endif -- 2.43.0