[u/mrichter/AliRoot.git] / TRD / AliTRDtrackerDebug.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * 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.                  *
 **************************************************************************/

/* $Id: AliTRDtrackerDebug.cxx 23810 2008-02-08 09:00:27Z hristov $ */

///////////////////////////////////////////////////////////////////////////////
//                                                                           //
//  Tracker debug streamer                                                   //
//                                                                           //
//  Authors:                                                                 //
//    Alex Bercuci <A.Bercuci@gsi.de>                                        //
//    Markus Fasel <M.Fasel@gsi.de>                                          //
//                                                                           //
///////////////////////////////////////////////////////////////////////////////

#include "AliTRDtrackerDebug.h"

#include "TFile.h"
#include "TTree.h"
#include "TTreeStream.h"
#include "TLinearFitter.h"

#include "AliLog.h"
#include "AliTRDgeometry.h"
#include "AliTRDtrackV1.h"
#include "AliTRDseedV1.h"
#include "AliTRDseed.h"
#include "AliTRDcluster.h"

ClassImp(AliTRDtrackerDebug)

//____________________________________________________
AliTRDtrackerDebug::AliTRDtrackerDebug() : AliTRDtrackerV1()
	,fOutputStreamer(0x0)
	,fTree(0x0)
	,fTracklet(0x0)
	,fTrack(0x0)
	,fNClusters(0)
	,fAlpha(0.)
{
        //
	// Default constructor
	//
	fOutputStreamer = new TTreeSRedirector("TRD.Debug.root");
}

//____________________________________________________
AliTRDtrackerDebug::~AliTRDtrackerDebug()
{
	// destructor
	
	delete fOutputStreamer;
}


//____________________________________________________
void AliTRDtrackerDebug::Draw(const Option_t *)
{
// steer draw function
}


//____________________________________________________
Bool_t AliTRDtrackerDebug::Init()
{
// steer linking data for various debug streams	
	
	fTrack = new AliTRDtrackV1();
	fTree->SetBranchAddress("ncl", &fNClusters);
	fTree->SetBranchAddress("track.", &fTrack);
	return kTRUE;
}

//____________________________________________________
Bool_t AliTRDtrackerDebug::Open(const char *method)
{
	// Connect to the tracker debug file
	
	TDirectory *savedir = gDirectory; 
	TFile::Open("TRD.TrackerDebugger.root");
	fTree = (TTree*)gFile->Get(method);
	if(!fTree){
		AliInfo(Form("Can not find debug stream for the %s method.\n", method));
		savedir->cd();
		return kFALSE;
	}
	savedir->cd();
	return kTRUE;
}

//____________________________________________________
Int_t AliTRDtrackerDebug::Process()
{
// steer debug process threads
	
	for(int it = 0; it<fTree->GetEntries(); it++){
		if(!fTree->GetEntry(it)) continue;
		if(!fNClusters) continue;
		fAlpha = fTrack->GetAlpha();
		//printf("Processing track %d [%d] ...\n", it, fNClusters);
		ResidualsTrackletsTrack();

		const AliTRDseedV1 *tracklet = 0x0;
		for(int ip = 5; ip>=0; ip--){
			if(!(tracklet = fTrack->GetTracklet(ip))) continue;
			if(!tracklet->GetN()) continue;
			
			ResidualsClustersTrack(tracklet);
			ResidualsClustersTracklet(tracklet);
			ResidualsClustersParametrisation(tracklet);
		}
	}
	return kTRUE;
}


//____________________________________________________
void AliTRDtrackerDebug::ResidualsClustersTrack(const AliTRDseedV1 *tracklet)
{
// Calculate averange distances from clusters to the TRD track	
	
	Double_t x[3]; 
	AliTRDcluster *c = 0x0;
	for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
		if(!(c = tracklet->GetClusters(ic))) continue;
		Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();

		// propagate track to cluster 
		PropagateToX(*fTrack, xc, 2.); 
		fTrack->GetXYZ(x);
		
		// transform to local tracking coordinates
		//Double_t xg =  x[0] * TMath::Cos(fAlpha) + x[1] * TMath::Sin(fAlpha); 
		Double_t yg = -x[0] * TMath::Sin(fAlpha) + x[1] * TMath::Cos(fAlpha);

		// apply tilt pad correction
		yc+= (zc - x[2]) * tracklet->GetTilt();
		
		Double_t dy = yc-yg;

		TTreeSRedirector &cstreamer = *fOutputStreamer;
		cstreamer << "ResidualsClustersTrack"
			<< "c.="   << c
			<< "dy="   << dy
			<< "\n";
	}
}

//____________________________________________________
void AliTRDtrackerDebug::ResidualsClustersTracklet(const AliTRDseedV1 *tracklet) const
{
// Calculates distances from clusters to tracklets
	
	Double_t x0 = tracklet->GetX0(), 
	         y0 = tracklet->GetYfit(0), 
	         ys = tracklet->GetYfit(1);
	         //z0 = tracklet->GetZfit(0), 
	         //zs = tracklet->GetZfit(1);
	
	AliTRDcluster *c = 0x0;
	for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
		if(!(c = tracklet->GetClusters(ic))) continue;
		Double_t xc = c->GetX(), yc = c->GetY()/*, zc = c->GetZ()*/;
		Double_t dy = yc- (y0-(x0-xc)*ys);

		//To draw  use : 
    //ResidualsClustersTracklet->Draw("TMath::Abs(10.*dy):TMath::ATan(ys)*TMath::RadToDeg()>>h(20, -40, 40)", "", "prof");
		TTreeSRedirector &cstreamer = *fOutputStreamer;
		cstreamer << "ResidualsClustersTracklet"
			<< "c.="   << c
			<< "ys="   << ys
			<< "dy="   << dy
			<< "\n";
	}
}

//____________________________________________________
void AliTRDtrackerDebug::ResidualsClustersParametrisation(const AliTRDseedV1 *tracklet) const
{
// Calculates distances from clusters to Rieman fit.
	
	// store cluster positions
	Double_t x0 = tracklet->GetX0();
	AliTRDcluster *c = 0x0;
	
	Double_t x[2]; Int_t ncl, mcl, jc;
	TLinearFitter fitter(3, "hyp2");
	for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
		if(!(c = tracklet->GetClusters(ic))) continue;
		Double_t xc = c->GetX(), yc = c->GetY()/*, zc = c->GetZ()*/;
		
		jc = ic; ncl = 0; mcl=0; fitter.ClearPoints();
		while(ncl<6){
			// update index
			mcl++;
			jc = ic + ((mcl&1)?-1:1)*(mcl>>1);

			if(jc<0 || jc>=35) continue;
			if(!(c = tracklet->GetClusters(jc))) continue;

			x[0] = c->GetX()-x0;
			x[1] = x[0]*x[0];
			fitter.AddPoint(x, c->GetY(), c->GetSigmaY2());
			ncl++;
		}
		fitter.Eval();
		Double_t dy = yc - fitter.GetParameter(0) -fitter.GetParameter(1) * (xc-x0) - fitter.GetParameter(2)* (xc-x0)*(xc-x0); 
	
		TTreeSRedirector &cstreamer = *fOutputStreamer;
		cstreamer << "ResidualsClustersParametrisation"
			<< "dy="   << dy
			<< "\n";
	}
}


//____________________________________________________
void AliTRDtrackerDebug::ResidualsTrackletsTrack() const
{
// Calculates distances from tracklets to the TRD track.
	
	if(fTrack->GetNumberOfTracklets() < 6) return;

	// build a working copy of the tracklets attached to the track 
	// and initialize working variables fX, fY and fZ
	AliTRDcluster *c = 0x0;
	AliTRDseedV1 tracklet[6] = {0x0, 0x0, 0x0, 0x0, 0x0, 0x0};
	const AliTRDseedV1 *ctracklet = 0x0;
	for(int ip = 0; ip<6; ip++){
		if(!(ctracklet = fTrack->GetTracklet(ip))) continue;
		tracklet[ip] = (*ctracklet); 
		Double_t x0 = tracklet[ip].GetX0();
		for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
			if(!(c = tracklet[ip].GetClusters(ic))) continue;
			Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
			tracklet[ip].SetX(ic, xc-x0);
			tracklet[ip].SetY(ic, yc);
			tracklet[ip].SetZ(ic, zc);
		}
	}
	
	// Do a Rieman fit (with tilt correction) for all tracklets 
	// except the one which is tested. 
	// (Based on AliTRDseed::IsOK() return false)
	for(int ip=0; ip<6; ip++){
		// reset tracklet to be tested
		Double_t x0 = tracklet[ip].GetX0();
		new(&tracklet[ip]) AliTRDseedV1();
		tracklet[ip].SetX0(x0);

		// fit Rieman with tilt correction
		AliTRDseedV1::FitRiemanTilt(&tracklet[0], kTRUE);

		// make a copy of the fit result
		Double_t 
			y0   = tracklet[ip].GetYref(0),
			dydx = tracklet[ip].GetYref(1),
			z0   = tracklet[ip].GetZref(0),
			dzdx = tracklet[ip].GetZref(1);

		// restore tracklet
		tracklet[ip] = (*fTrack->GetTracklet(ip)); 
		for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
			if(!(c = tracklet[ip].GetClusters(ic))) continue;
			Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
			tracklet[ip].SetX(ic, xc-x0);
			tracklet[ip].SetY(ic, yc);
			tracklet[ip].SetZ(ic, zc);
		}		
		
		// fit clusters
		AliTRDseedV1 ts(tracklet[ip]);
		ts.SetYref(0, y0); ts.SetYref(1, dydx);
		ts.SetZref(0, z0); ts.SetZref(1, dzdx);
		ts.Update();

		// save results for plotting
		Int_t plane   = tracklet[ip].GetPlane();
		Double_t dy   = tracklet[ip].GetYfit(0) - ts.GetYfit(0);
		Double_t tgy  = tracklet[ip].GetYfit(1);
		Double_t dtgy = (tracklet[ip].GetYfit(1) - ts.GetYfit(1))/(1. + tracklet[ip].GetYfit(1) * ts.GetYfit(1));
		Double_t dz   = tracklet[ip].GetZfit(0) - ts.GetZfit(0);
		Double_t tgz  = tracklet[ip].GetZfit(1);
		Double_t dtgz = (tracklet[ip].GetZfit(1) - ts.GetZfit(1))/(1. + tracklet[ip].GetZfit(1) * ts.GetZfit(1));
		TTreeSRedirector &cstreamer = *fOutputStreamer;
		cstreamer << "ResidualsTrackletsTrack"
			<< "ref.="   << &tracklet[ip]
			<< "fit.="   << &ts
			<< "plane="  << plane
			<< "dy="     << dy
			<< "tgy="    << tgy
			<< "dtgy="   << dtgy
			<< "dz="     << dz
			<< "tgz="    << tgz
			<< "dtgz="   << dtgz
			<< "\n";
	}
}
Commit	Line	Data
eb38ed55	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	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	**************************************************************************/
	15
	16	/* $Id: AliTRDtrackerDebug.cxx 23810 2008-02-08 09:00:27Z hristov $ */
	17
	18	///////////////////////////////////////////////////////////////////////////////
	19	// //
	20	// Tracker debug streamer //
	21	// //
	22	// Authors: //
	23	// Alex Bercuci <A.Bercuci@gsi.de> //
	24	// Markus Fasel <M.Fasel@gsi.de> //
	25	// //
	26	///////////////////////////////////////////////////////////////////////////////
	27
	28	#include "AliTRDtrackerDebug.h"
	29
	30	#include "TFile.h"
	31	#include "TTree.h"
	32	#include "TTreeStream.h"
	33	#include "TLinearFitter.h"
	34
	35	#include "AliLog.h"
	36	#include "AliTRDgeometry.h"
	37	#include "AliTRDtrackV1.h"
	38	#include "AliTRDseedV1.h"
	39	#include "AliTRDseed.h"
	40	#include "AliTRDcluster.h"
	41
	42	ClassImp(AliTRDtrackerDebug)
	43
	44	//____________________________________________________
	45	AliTRDtrackerDebug::AliTRDtrackerDebug() : AliTRDtrackerV1()
	46	,fOutputStreamer(0x0)
	47	,fTree(0x0)
	48	,fTracklet(0x0)
	49	,fTrack(0x0)
	50	,fNClusters(0)
	51	,fAlpha(0.)
	52	{
	53	//
	54	// Default constructor
	55	//
	56	fOutputStreamer = new TTreeSRedirector("TRD.Debug.root");
	57	}
	58
	59	//____________________________________________________
	60	AliTRDtrackerDebug::~AliTRDtrackerDebug()
	61	{
	62	// destructor
	63
	64	delete fOutputStreamer;
65	}
66
67
68	//____________________________________________________
69	void AliTRDtrackerDebug::Draw(const Option_t *)
70	{
71	// steer draw function
72	}
73
74
75	//____________________________________________________
76	Bool_t AliTRDtrackerDebug::Init()
77	{
78	// steer linking data for various debug streams
79
80	fTrack = new AliTRDtrackV1();
81	fTree->SetBranchAddress("ncl", &fNClusters);
82	fTree->SetBranchAddress("track.", &fTrack);
83	return kTRUE;
84	}
85
86	//____________________________________________________
87	Bool_t AliTRDtrackerDebug::Open(const char *method)
88	{
89	// Connect to the tracker debug file
90
91	TDirectory *savedir = gDirectory;
92	TFile::Open("TRD.TrackerDebugger.root");
93	fTree = (TTree*)gFile->Get(method);
94	if(!fTree){
95	AliInfo(Form("Can not find debug stream for the %s method.\n", method));
96	savedir->cd();
97	return kFALSE;
98	}
99	savedir->cd();
100	return kTRUE;
101	}
102
103	//____________________________________________________
104	Int_t AliTRDtrackerDebug::Process()
105	{
106	// steer debug process threads
107
108	for(int it = 0; it<fTree->GetEntries(); it++){
109	if(!fTree->GetEntry(it)) continue;
110	if(!fNClusters) continue;
111	fAlpha = fTrack->GetAlpha();
112	//printf("Processing track %d [%d] ...\n", it, fNClusters);
113	ResidualsTrackletsTrack();
114
115	const AliTRDseedV1 *tracklet = 0x0;
116	for(int ip = 5; ip>=0; ip--){
117	if(!(tracklet = fTrack->GetTracklet(ip))) continue;
118	if(!tracklet->GetN()) continue;
119
120	ResidualsClustersTrack(tracklet);
121	ResidualsClustersTracklet(tracklet);
122	ResidualsClustersParametrisation(tracklet);
123	}
124	}
125	return kTRUE;
126	}
127
128
129	//____________________________________________________
130	void AliTRDtrackerDebug::ResidualsClustersTrack(const AliTRDseedV1 *tracklet)
131	{
132	// Calculate averange distances from clusters to the TRD track
133
134	Double_t x[3];
135	AliTRDcluster *c = 0x0;
136	for(int ic=0; ic<35/AliTRDseed:knTimebins/; ic++){
137	if(!(c = tracklet->GetClusters(ic))) continue;
138	Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
139
140	// propagate track to cluster
141	PropagateToX(*fTrack, xc, 2.);
142	fTrack->GetXYZ(x);
143
144	// transform to local tracking coordinates
145	//Double_t xg = x[0] * TMath::Cos(fAlpha) + x[1] * TMath::Sin(fAlpha);
146	Double_t yg = -x[0] * TMath::Sin(fAlpha) + x[1] * TMath::Cos(fAlpha);
147
148	// apply tilt pad correction
149	yc+= (zc - x[2]) * tracklet->GetTilt();
150
151	Double_t dy = yc-yg;
152
153	TTreeSRedirector &cstreamer = *fOutputStreamer;
154	cstreamer << "ResidualsClustersTrack"
155	<< "c.=" << c
156	<< "dy=" << dy
157	<< "\n";
158	}
159	}
160
161	//____________________________________________________
162	void AliTRDtrackerDebug::ResidualsClustersTracklet(const AliTRDseedV1 *tracklet) const
163	{
164	// Calculates distances from clusters to tracklets
165
166	Double_t x0 = tracklet->GetX0(),
167	y0 = tracklet->GetYfit(0),
168	ys = tracklet->GetYfit(1);
169	//z0 = tracklet->GetZfit(0),
170	//zs = tracklet->GetZfit(1);
171
172	AliTRDcluster *c = 0x0;
173	for(int ic=0; ic<35/AliTRDseed:knTimebins/; ic++){
174	if(!(c = tracklet->GetClusters(ic))) continue;
175	Double_t xc = c->GetX(), yc = c->GetY()/, zc = c->GetZ()/;
176	Double_t dy = yc- (y0-(x0-xc)*ys);
177
178	//To draw use :
179	//ResidualsClustersTracklet->Draw("TMath::Abs(10.dy):TMath::ATan(ys)TMath::RadToDeg()>>h(20, -40, 40)", "", "prof");
180	TTreeSRedirector &cstreamer = *fOutputStreamer;
181	cstreamer << "ResidualsClustersTracklet"
182	<< "c.=" << c
183	<< "ys=" << ys
184	<< "dy=" << dy
185	<< "\n";
186	}
187	}
188
189	//____________________________________________________
190	void AliTRDtrackerDebug::ResidualsClustersParametrisation(const AliTRDseedV1 *tracklet) const
191	{
192	// Calculates distances from clusters to Rieman fit.
193
194	// store cluster positions
195	Double_t x0 = tracklet->GetX0();
196	AliTRDcluster *c = 0x0;
197
198	Double_t x[2]; Int_t ncl, mcl, jc;
199	TLinearFitter fitter(3, "hyp2");
200	for(int ic=0; ic<35/AliTRDseed:knTimebins/; ic++){
201	if(!(c = tracklet->GetClusters(ic))) continue;
202	Double_t xc = c->GetX(), yc = c->GetY()/, zc = c->GetZ()/;
203
204	jc = ic; ncl = 0; mcl=0; fitter.ClearPoints();
205	while(ncl<6){
206	// update index
207	mcl++;
208	jc = ic + ((mcl&1)?-1:1)*(mcl>>1);
209
210	if(jc<0 \|\| jc>=35) continue;
211	if(!(c = tracklet->GetClusters(jc))) continue;
212
213	x[0] = c->GetX()-x0;
214	x[1] = x[0]*x[0];
215	fitter.AddPoint(x, c->GetY(), c->GetSigmaY2());
216	ncl++;
217	}
218	fitter.Eval();
219	Double_t dy = yc - fitter.GetParameter(0) -fitter.GetParameter(1) * (xc-x0) - fitter.GetParameter(2)* (xc-x0)*(xc-x0);
220
221	TTreeSRedirector &cstreamer = *fOutputStreamer;
222	cstreamer << "ResidualsClustersParametrisation"
223	<< "dy=" << dy
224	<< "\n";
225	}
226	}
227
228
229	//____________________________________________________
230	void AliTRDtrackerDebug::ResidualsTrackletsTrack() const
231	{
232	// Calculates distances from tracklets to the TRD track.
233
234	if(fTrack->GetNumberOfTracklets() < 6) return;
235
236	// build a working copy of the tracklets attached to the track
237	// and initialize working variables fX, fY and fZ
238	AliTRDcluster *c = 0x0;
239	AliTRDseedV1 tracklet[6] = {0x0, 0x0, 0x0, 0x0, 0x0, 0x0};
240	const AliTRDseedV1 *ctracklet = 0x0;
241	for(int ip = 0; ip<6; ip++){
242	if(!(ctracklet = fTrack->GetTracklet(ip))) continue;
243	tracklet[ip] = (*ctracklet);
244	Double_t x0 = tracklet[ip].GetX0();
245	for(int ic=0; ic<35/AliTRDseed:knTimebins/; ic++){
246	if(!(c = tracklet[ip].GetClusters(ic))) continue;
247	Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
248	tracklet[ip].SetX(ic, xc-x0);
249	tracklet[ip].SetY(ic, yc);
250	tracklet[ip].SetZ(ic, zc);
251	}
252	}
253
254	// Do a Rieman fit (with tilt correction) for all tracklets
255	// except the one which is tested.
256	// (Based on AliTRDseed::IsOK() return false)
257	for(int ip=0; ip<6; ip++){
258	// reset tracklet to be tested
259	Double_t x0 = tracklet[ip].GetX0();
260	new(&tracklet[ip]) AliTRDseedV1();
261	tracklet[ip].SetX0(x0);
262
263	// fit Rieman with tilt correction
264	AliTRDseedV1::FitRiemanTilt(&tracklet[0], kTRUE);
265
266	// make a copy of the fit result
267	Double_t
268	y0 = tracklet[ip].GetYref(0),
269	dydx = tracklet[ip].GetYref(1),
270	z0 = tracklet[ip].GetZref(0),
271	dzdx = tracklet[ip].GetZref(1);
272
273	// restore tracklet
274	tracklet[ip] = (*fTrack->GetTracklet(ip));
275	for(int ic=0; ic<35/AliTRDseed:knTimebins/; ic++){
276	if(!(c = tracklet[ip].GetClusters(ic))) continue;
277	Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
278	tracklet[ip].SetX(ic, xc-x0);
279	tracklet[ip].SetY(ic, yc);
280	tracklet[ip].SetZ(ic, zc);
281	}
282
283	// fit clusters
284	AliTRDseedV1 ts(tracklet[ip]);
285	ts.SetYref(0, y0); ts.SetYref(1, dydx);
286	ts.SetZref(0, z0); ts.SetZref(1, dzdx);
287	ts.Update();
288
289	// save results for plotting
290	Int_t plane = tracklet[ip].GetPlane();
291	Double_t dy = tracklet[ip].GetYfit(0) - ts.GetYfit(0);
292	Double_t tgy = tracklet[ip].GetYfit(1);
293	Double_t dtgy = (tracklet[ip].GetYfit(1) - ts.GetYfit(1))/(1. + tracklet[ip].GetYfit(1) * ts.GetYfit(1));
294	Double_t dz = tracklet[ip].GetZfit(0) - ts.GetZfit(0);
295	Double_t tgz = tracklet[ip].GetZfit(1);
296	Double_t dtgz = (tracklet[ip].GetZfit(1) - ts.GetZfit(1))/(1. + tracklet[ip].GetZfit(1) * ts.GetZfit(1));
297	TTreeSRedirector &cstreamer = *fOutputStreamer;
298	cstreamer << "ResidualsTrackletsTrack"
299	<< "ref.=" << &tracklet[ip]
300	<< "fit.=" << &ts
301	<< "plane=" << plane
302	<< "dy=" << dy
303	<< "tgy=" << tgy
304	<< "dtgy=" << dtgy
305	<< "dz=" << dz
306	<< "tgz=" << tgz
307	<< "dtgz=" << dtgz
308	<< "\n";
309	}
310	}
311