Home » PANDA » PandaRoot » Bugs, Fixes, Releases » Fix in lhetrack
Fix in lhetrack [message #6339] |
Thu, 10 April 2008 10:57 |
StefanoSpataro
Messages: 2736 Registered: June 2005 Location: Torino
|
first-grade participant |
From: *physik.uni-giessen.de
|
|
Dear all,
I have found a problem in lhetrack code that gave me a crash while processing many events. In particular, the function PndTpcLheTrackCuts::Circle3pnts(...) which estimates from the first and the last point of the track candidate the centre and the radius of the circle crossing those two points and the origin (0,0).
You can see below the old version in PndTpcLheTrackCuts.cxx
Toggle Spoiler
void PndTpcLheTrackCuts::Circle3pnts(Double_t x[],Double_t y[], Double_t r[]) {
// calc center and R of circle from 3 points
Double_t m_a = (y[1] - y[0]) / (x[1] - x[0]);
Double_t m_b = (y[2] - y[1]) / (x[2] - x[1]);
if (m_a != m_b) {
Double_t x_c = .5*(m_a* m_b*(y[0] - y[2]) +
m_b*(x[0] + x[1]) -
m_a*(x[1] + x[2])) / (m_b - m_a);
Double_t y_c = - (x_c - .5*(x[0] + x[1])) / m_a +
.5*(y[0] + y[1]);
Double_t dely = y_c - y[0];
Double_t delx = x_c - x[0];
Double_t rad =TMath::Sqrt(delx*delx + dely*dely);
r[0] = x_c;
r[1] = y_c;
r[2] = rad;
}
}
As you can see, the function is not protected agains division by zero (so when x[1]-x[0]==0 or x[2]-x[1]==0). Moreover, the case m_a==m_b (the three points lay on a straight line) does not give results.
I have corrected it in the repository in the following way:
Toggle Spoiler
void PndTpcLheTrackCuts::Circle3pnts(Double_t x[],Double_t y[], Double_t r[]) {
// calc center and R of circle from 3 points
Double_t x_c = 0., y_c = 0.;
if ( (y[1] - y[0])*(x[2] - x[1]) == (x[1] - x[0])*(y[2] - y[1]))
{
cout << " -W- PndTpcLheTrackCuts::Circle3pnts: The points lay on a straight line" << endl;
x_c = 10000.;
y_c = 10000.;
}
else
{
if ((x[1] - x[0])!=0 && (x[2] - x[1])!=0)
{
Double_t m_a= (y[1] - y[0]) / (x[1] - x[0]);
Double_t m_b = (y[2] - y[1]) / (x[2] - x[1]);
x_c = .5*(m_a* m_b*(y[0] - y[2]) +
m_b*(x[0] + x[1]) -
m_a*(x[1] + x[2])) / (m_b - m_a);
y_c = - (x_c - .5*(x[0] + x[1])) / m_a +
.5*(y[0] + y[1]);
}
else
{
Double_t m_a = (x[1] - x[0]) / (y[1] - y[0]);
Double_t m_b = (x[2] - x[1]) / (y[2] - y[1]);
y_c = .5*(m_a* m_b*(x[0] - x[2]) +
m_b*(y[0] + y[1]) -
m_a*(x[1] + x[2])) / (m_b - m_a);
x_c = - (y_c - .5*(y[0] + y[1])) / m_a +
.5*(x[0] + x[1]);
}
} // end of not collinear points condition
Double_t dely = y_c - y[0];
Double_t delx = x_c - x[0];
r[0] = x_c;
r[1] = y_c;
r[2] = TMath::Sqrt(delx*delx + dely*dely);
}
If the three points lay on a straight line, the circle centre is set at (10000,10000) just to put a high number, probably the code will kill those tracks or be able to recover them.
If there is the division by zero, I have simply exchanged x->y y->x (in this system there is no division by zero anymore), calculated the centre coordinates and re-transformated them into the original system.
I hope I have not messed up things, however now I do not have the crash anymore (or at least there).
Comments and checks are welcome.
|
|
|
Update of lhetrack [message #6357 is a reply to message #6339] |
Thu, 10 April 2008 19:25 |
StefanoSpataro
Messages: 2736 Registered: June 2005 Location: Torino
|
first-grade participant |
From: *physik.uni-giessen.de
|
|
Dear all,
since version 2466 lhetrack package has several "important" updates.
1) First of all, a small bug that produced analysis crashes was solved (see entry before).
2) Now it is possible to include MVD in the fitting by a task option, without recompiling (as done before)
3) It is possible to run the code with pure MC points, with MCPoints smeared by fixing the spatial resolution, with reconstructed hits (so TpcCluster and MVDHits).
The construction is very easy. You have just to add in your analysis macro:
// ----- LHETRACK ---------------------------------
PndTpcLheHitsMaker* trackMS = new PndTpcLheHitsMaker("Tracking routine");
trackMS->SetTpcMode(2, -1); // 0 OFF, 1 TpcPoint, 2 TpcCluster // TpcPoint smearing [cm], if negative no smearing
trackMS->SetMvdMode(2, -1); // 0 OFF, 1 MVDPoint, 2 MVDHit // MVDPoint smearing [cm], if negative no smearing
fRun->AddTask(trackMS);
PndTpcLheTrackFinder* trackFinder = new PndTpcLheTrackFinder();
fRun->AddTask(trackFinder);
PndTpcLheTrackFitter* trackFitter = new PndTpcLheTrackFitter("fitting");
fRun->AddTask(trackFitter);
With the SetXXXMode(Int_t mode, Float_t res) function for TPC and MVD you can switch OFF your detector (mode==0), or decide to use just Points (mode==1), Points smeared according to a space resolution value res (for TPC sigma_Z = 2*res), or reconstructed hits (mode==2).
Under macro/lhetrack you can run the reconstruction chain:
run_sim_tpcmvd.C will produce MC points
run_track_tpcmvd.C will run digitization+reconstruction for TPC and MVD, and will use the reco points for the fitting
plot_pT.C will plot the results of the fitting.
Known (and unsolved) problems (at the moment)
1) If lhetrack runs in a 2nd step macro, it is crashing. The problem seems to stay in the sigma for MVD hits. The functions MVDHits::GetDx() y and z called by lhetrack code require the geometry which is not in the file, so the macro crashes. Adding the sim file as fRun->AddFriend() does not help. For these reasons everything should run in the same macro, at the moment.
2) It seems that if one turns on the RiemannFit task, the macro goes well. But when launching the plot_pT.C macro to show results, this crashes. It is like it is not able to access properly to the tree. This should be investigated.
3) While for the MCPoints is easy to get the MCTrack index, for reco hits this is not possible. In this case the TrackID is set to -1.
4) lhetrack has some tools to fill track candidated from the MC and from the fit, by using TpcPoints TCA. Now lhetrack works even with other kind of data objects such as MVDPoints and reco hits. So this part should be adjusted. For MC it will be quite simple, but for reco... that's another story.
Some feedback is welcome, of course.
|
|
|
|
|
|
Re: Pid with lhetrack [message #6519 is a reply to message #6516] |
Fri, 18 April 2008 11:27 |
StefanoSpataro
Messages: 2736 Registered: June 2005 Location: Torino
|
first-grade participant |
From: *physik.uni-giessen.de
|
|
Hello,
I attach the plot for MVD dE/dx versus momentum, where:
MVD dE/dx is the total energy deposited in MVD strip/pixel (MvdHits, no MC) divided by the number of MVD strips/pixels fired -> no correction for thickness or angle of the track inside the sensor
momentum is the reconstructed momentum using reco hits, so TpcCluster and MvdHit, calculated by lhetrack.
Here the plot for GEANT3
and for GEANT4
The plots were generated with Box generator, uniform in p [0,1 GeV/c], theta [20°,120°] and phi [0°,360°]. I have generated separately protons, kaons pions and electrons.
You can see that he momentum reconstruction seems to stop below 100-150 MeV/c, but it depends also on the particle type.
Second, you can see that in "kaon" events there is a structure in the pion region, most probably connected to kaon decays.
Third, electrons have different energy loss in G3 and G4!!!
Here you are the direct comparison:
In Geant3 electrons are losing more than twice the energy lost in Geant4. In Geant4 we have a MIPS point similar to the one of the pions, while in G3 it seems that we could use MVD to discriminate between pions and electrons
However, I remember you the dE/dx plot for STT done with geant3 (that I have also showed at CHEP), and even there the electrons stayed much higher than pions.
But from dE/dx plots on tpc of the data particle booklet (page 269), it seems G3 is closer to reality...
So... another problem of G4?
[Updated on: Fri, 18 April 2008 11:27] Report message to a moderator
|
|
|
|
Re: Pid with lhetrack [message #6522 is a reply to message #6520] |
Fri, 18 April 2008 11:52 |
StefanoSpataro
Messages: 2736 Registered: June 2005 Location: Torino
|
first-grade participant |
From: *physik.uni-giessen.de
|
|
Hi Sebastian,
I remember your considerations, but I have supposed they were connected to the bad behaviour of the g4 "stepLimiter" option (as found by Dima) that now is switched off. But it does not seems the case, there is still something strange.
About trusting on g3 or g4... I am not too much sure on what should be the correct response. Are you able to produce the same plots for tpc, just to cross check? In theory the plot in the data particle booklet should be our reference, and in that case electrons stay higher than pions (so in this case g3 should be more realistic).
About the second structure, as I wrote before,
kaons can decay into pions, and these are placed in the region you evidence. This is the reason "I think" of the second structure. "Kaon" means "event where I have created one primary kaon", not "the particle is a kaon".
I am not applying PID because I am not using MC info at all, so at the moment I cannot tell you exactly what they are, but I think my suspect should be true.
Opinions are welcome..
|
|
|
|
|
|
|
|
Re: Update of lhetrack [message #6853 is a reply to message #6802] |
Fri, 06 June 2008 18:03 |
StefanoSpataro
Messages: 2736 Registered: June 2005 Location: Torino
|
first-grade participant |
From: *physik.uni-giessen.de
|
|
Dear all,
in svn you can find the latest version of lhetrack code.
By using the geometrical extrapolation, a correlation with TOF and EMC detectors is now provided: the object PndLhePidTrack contains together the informations from TPC, MVD, EMC and TOF.
The following plot shows the tof vs momentum distribution for protons (the top one), kaons (intermediate), pions, muons, and electrons (the lower bands almost horizontal):
and here you are the E/P (energy loss in EMC divided by momentum) versus momentum plot for electrons (the horizontal line at 1), protons (the band which is going up), pions and muons (the band which is going down) and kaons (mixed with the pions, but a higher momenta)
The bump at E/B~0.8 should be due to secondary pions generated by the decay of a primary kaon.
No MonteCarlo information was used for the previous plots, everything comes from full simulation and full reconstruction.
Still the extrapolation should be improved, and the path lenght to the tof detector has to be calculated.
-
Attachment: TofTof.gif
(Size: 17.75KB, Downloaded 793 times)
-
Attachment: EmcEloss.gif
(Size: 19.04KB, Downloaded 762 times)
|
|
|
|
|
|
|
|
Re: Dirc+TOF [message #7182 is a reply to message #7180] |
Wed, 06 August 2008 14:42 |
asanchez
Messages: 350 Registered: March 2006
|
first-grade participant |
From: *gsi.de
|
|
hi again,
i thought the z0, was the z coordinante at (x0,y0)
where x0,y0 are the coordinates of the center of the circle described by the particle in its movement.
So my point was, that primaries are statrting at 0,0,0, where the extrapolation to tof radius to determine the phi is correct.
But if the particle is not starting at 0,0,0 the extrapolation to tof radius is not valid anymore at least you take also into a count the Phi value for the first hit in the track.
What do you think?
Alicia.
|
|
|
|
|
|
|
|
|
Goto Forum:
Current Time: Tue Oct 15 14:21:13 CEST 2024
Total time taken to generate the page: 0.00722 seconds
|