1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
17 ///////////////////////////////////////////////////////////////////////////////
19 // This class is a container class for the laser beam positions of //
20 // the TPC laser system (ALICE-INT-2002-022 v.1 ). //
21 // The system consits of 6 Laser Rods, each contains four mirror bundels //
22 // on each side of the TPC. Each bundle has seven mirrors. //
24 // The TPC side 0 corresponds to the "Shaft Side (A)", side 1 corresponds //
25 // to the "Muon Side (C)".
26 // The laser rods are counted counter clockwise starting with the one //
27 // with the smalles phi angle. //
28 // The mirror bundles in one rod for each side are counted starting from //
29 // form the greatest z value to the smalles. //
30 // The beams in each bundel are counted counter clockwise. Having all beams //
31 // pointing downward, the first beam will be the leftmost. //
33 ///////////////////////////////////////////////////////////////////////////////
36 Create File with design positions
38 .L AliTPCLaserTrack.cxx+
39 AliTPCLaserTracks ltracks
40 ltracks.WriteTreeDesignData()
41 TFile f("LaserTracksDesing.root");
42 TTree * tree =(TTree*)f.Get("LaserTracks");
46 #include <Riostream.h>
48 #include <TPolyLine3D.h>
52 #include <TEventList.h>
60 #include "AliTPCLaserTracks.h"
62 ////////////////////////////////////////////////////////////////////////
63 // Class AliTPCLaserTracks
64 ////////////////////////////////////////////////////////////////////////
66 ClassImp(AliTPCLaserTracks)
68 AliTPCLaserTracks::AliTPCLaserTracks():TObject(),
88 // AliTPCLaserTracks default constructor
92 //_______________________________________________________________________
93 AliTPCLaserTracks::AliTPCLaserTracks(Int_t npoints):TObject(),
113 // AliTPCLaserTracks constructor to initialise array of points
119 //_______________________________________________________________________
120 Double_t AliTPCLaserTracks::FindBeamLength(TVector3 vF, TVector3 vP)
123 // Find distance between mirror and the intersection between
124 // inner and outer field cage respectively
126 // use the pq formula to find the intersection point between
127 // the beam and the inner and outer fieldcage cylinders
128 // The formulae solved are
131 // xPoint = vF.X() + n*vP.X() = r*cos(phi)
132 // yPoint = vF.Y() + n*vP.Y() = r*sin(phi)
133 // zPoint = vF.Y() + n*vP.Y()
135 // get n from the first two equations, calculate x,y,zPoint
137 // vF is the mirror position and vP the beam direction
138 Double_t r[2]; //smalles r[0] and largest r[1] pad radius
139 r[0] = 80.0; //smalles pad radius
140 r[1] = 260.0; //largest pad radius
142 Double_t fxpxfypy = vF.X()*vP.X()+vF.Y()*vP.Y();
143 Double_t px2py2 = vP.X()*vP.X()+vP.Y()*vP.Y();
145 for (Int_t i=0; i<2; i++){
146 Double_t rad = (fxpxfypy/px2py2)*(fxpxfypy/px2py2)-
147 (vF.X()*vF.X()+vF.Y()*vF.Y()-r[i]*r[i])/px2py2;
150 Double_t n1 = -(fxpxfypy)/px2py2+TMath::Sqrt(rad);
151 TVector3 vI1(vF.X()+n1*vP.X(),
155 Double_t n2 = -(fxpxfypy)/px2py2-TMath::Sqrt(rad);
156 TVector3 vI2(vF.X()+n2*vP.X(),
161 //if we cross two boarders on our way return the closer boarder
162 if ( (n1>0) && (n2>0) ) {
163 if ( (vF-vI1).Mag() <= (vF-vI2).Mag() )
164 return (vF-vI1).Mag();
166 return (vF-vI2).Mag();
170 return (vF-vI1).Mag();
173 return (vF-vI2).Mag();
178 return 3.4e38; //not found!!
181 //_______________________________________________________________________
182 void AliTPCLaserTracks::WriteTreeDesignData()
185 // Write a tree with the design data of the laser track positions
188 TFile *f = TFile::Open("LaserTracksDesing.root","recreate");
190 AliTPCLaserTracks *ltp = new AliTPCLaserTracks(2);
192 TTree *t = new TTree("LaserTracks","LaserTracks");
193 t->Branch("LaserTracks","AliTPCLaserTracks",<p);
196 Double_t phiBeam[7] = {-31.8,-16.,-9.2,2.5,9.2,16.,31.8};
197 Double_t phiRod[6] = {40.,100.,160.,220.,280.,340.}; //-20. for the muon side
198 // Double_t zBundleCorse[2][4] = {{10,79,163,241},{13,85,169,247}};
199 Double_t zBundleCorse[2][4] = {{241.,163.,79.,10.},{247.,169.,85.,13.}};
200 Double_t zBeamFine[7] = {1.45,1.3,1.15,1.3,1.15,1.3,1.45};
203 //loop over TPC side -- 0:A (Shaft) side -- 1:C (Muon) side
204 for (Int_t side=0; side<2; side++){
205 //loop over laser rods -- counterclockwise
206 for (Int_t rod=0; rod<6; rod++){
210 Double_t phiRod2 = phiRod[rod]-side*20;
211 vRod.SetMagPhi(254.25, //center of the rod at 254.25cm
213 (phiRod2) //-20 deg on C-Side
216 //loop over bundle -- counted from large z to small z
217 for (Int_t bundle=0; bundle<4; bundle++){
218 //center of the bundle; not yet rotated
219 Int_t bundleS = bundle;
220 // if ( side == 1 ) bundleS = 4-bundle;
221 if ( side == 1 ) bundleS = 3-bundle;
224 vBundle.SetMagPhi(.65,
226 //? TMath::Power(-1,side)*
227 (phiRod2+180.-90.*bundleS)
232 for (Int_t beam=0; beam<7; beam++){
237 vBeam.SetMagPhi(1.,TMath::DegToRad()*(phiRod2+i*60));
241 TVector2 xyMirror = vRod+vBundle+vBeam;
242 Double_t zMirror = (zBundleCorse[rod%2][bundleS]+zBeamFine[beam])*TMath::Power(-1,side);
243 TVector3 v3Mirror(xyMirror.X(),xyMirror.Y(),zMirror);
245 Double_t phi = TMath::DegToRad()*(phiRod2+180+phiBeam[beam]*TMath::Power(-1,side));
246 Double_t theta = TMath::DegToRad()*90;
249 v3Beam.SetMagThetaPhi(1,theta,phi); //Direction Vector
250 v3Beam.SetMagThetaPhi(FindBeamLength(v3Mirror,v3Beam), theta, phi);
257 ltp->SetBundle(bundleS);
260 ltp->SetX(xyMirror.X());
261 ltp->SetY(xyMirror.Y());
265 ltp->SetTheta(theta);
267 // ltp->InitPoints(2);
268 ltp->SetPoint(0,xyMirror.X(),xyMirror.Y(),zMirror);
269 ltp->SetPoint(1,xyMirror.X()+v3Beam.X(),xyMirror.Y()+v3Beam.Y(),zMirror+v3Beam.Z());
274 << " Bundle: " << bundleS
282 // cout << "Filled!!!" << endl;
295 //_______________________________________________________________________
296 TPolyLine3D* AliTPCLaserTracks::GetLine()
299 // Make a polilyne object oout of the points
301 if ( fNpoints ) return new TPolyLine3D(fNpoints,fXarr,fYarr,fZarr);
306 //_______________________________________________________________________
307 TObjArray* AliTPCLaserTracks::GetLines(const Char_t* file, const Char_t *cuts)
310 // Read the input files with the laser track
311 // Make a array of polylines for visualization
313 TObjArray *array = new TObjArray;
315 TFile *f = TFile::Open(file,"read");
316 TTree *tree = (TTree*)f->Get("LaserTracks");
318 AliTPCLaserTracks *ltp = new AliTPCLaserTracks();
319 // TEventList *evList = new TEventList("evList");
320 TEventList evList("evList");
322 tree->SetBranchAddress("LaserTracks",<p);
324 tree->Draw(">>evList",cuts,"goff");
326 // cout << "N: " << evList.GetN() << endl;
328 for (Int_t ev = 0; ev < evList.GetN(); ev++){
329 // cout << ev << endl;
330 tree->GetEntry( evList.GetEntry(ev) );
331 array->Add(ltp->GetLine());
333 // cout << "delte f"<< endl;
336 // cout << "evlist" << endl;
341 //_______________________________________________________________________
342 void AliTPCLaserTracks::InitPoints()
349 fXarr = new Double_t[fMaxSize];
350 fYarr = new Double_t[fMaxSize];
351 fZarr = new Double_t[fMaxSize];
354 //_______________________________________________________________________
355 Int_t AliTPCLaserTracks::SetPoint(Int_t point, Double_t x, Double_t y, Double_t z)
358 // Set point to point arrays
361 if ( !fXarr || !fYarr || !fZarr ) return 1;
362 if ( point > fNpoints-1 ) return 1;
370 //_______________________________________________________________________
371 Int_t AliTPCLaserTracks::FindMirror(Char_t *file, Double_t x, Double_t y, Double_t z, Double_t phi)
374 // return the id of the mirror in 'file' that maches best the given parameters
377 TFile *f = TFile::Open(file,"read");
378 TTree *tree = (TTree*)f->Get("LaserTracks");
380 AliTPCLaserTracks *ltp = new AliTPCLaserTracks();
381 TEventList evList("evList");
383 tree->SetBranchAddress("LaserTracks",<p);
386 s = "abs(fX-"; s+= x;
387 s+= ")<3&&abs(fY-"; s+=y;
388 s+= ")<3&&abs(fZ-"; s+= z;
390 // s+= "&&((abs(fPhi-";s+= phi;
391 // s+= ")<.06)||(abs(fPhi-";s+= (phi-TMath::DegToRad()*180);
393 cout << s.Data() << endl;
395 tree->Draw(">>evList",s.Data());
397 Double_t dphiMin=TMath::Pi();
400 cout << "nev: " << evList.GetN() << endl;
401 for (Int_t i=0; i<evList.GetN(); i++){
402 tree->GetEntry( evList.GetEntry(i) );
403 Double_t dphi = TMath::Abs(ltp->GetPhi() - phi);
404 if ( dphi > TMath::DegToRad()*180 ) dphi-=TMath::DegToRad()*180;
405 if ( dphi<dphiMin ) {
414 //_______________________________________________________________________
415 AliTPCLaserTracks::~AliTPCLaserTracks()
418 // standard destructor
421 // if ( fP ) delete fP;
422 if ( fXarr ) delete [] fXarr;
423 if ( fYarr ) delete [] fYarr;
424 if ( fZarr ) delete [] fZarr;