[u/mrichter/AliRoot.git] / TRD / qaRec / AliTRDtrackingResolution.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: AliTRDtrackingResolution.cxx 27496 2008-07-22 08:35:45Z cblume $ */

////////////////////////////////////////////////////////////////////////////
//                                                                        //
//  Reconstruction QA                                                     //
//                                                                        //
//  Authors:                                                              //
//    Markus Fasel <M.Fasel@gsi.de>                                       //
//                                                                        //
////////////////////////////////////////////////////////////////////////////

#include <cstring>


#include <TObjArray.h>
#include <TList.h>
#include <TH2.h>
#include <TH1.h>
#include <TF1.h>
#include <TProfile.h>
#include <TGraphErrors.h>
#include <TMath.h>
#include "TTreeStream.h"
#include "TGeoManager.h"

#include "AliAnalysisManager.h"
#include "AliTrackReference.h"
#include "AliTrackPointArray.h"
#include "AliCDBManager.h"

#include "AliTRDcluster.h"
#include "AliTRDseedV1.h"
#include "AliTRDtrackV1.h"
#include "AliTRDtrackerV1.h"
#include "AliTRDReconstructor.h"
#include "AliTRDrecoParam.h"

#include "AliTRDtrackInfo/AliTRDtrackInfo.h"
#include "AliTRDtrackingResolution.h"

ClassImp(AliTRDtrackingResolution)

//________________________________________________________
AliTRDtrackingResolution::AliTRDtrackingResolution(const char * name):
  AliAnalysisTask(name, "")
  ,fTrackObjects(0x0)
  ,fOutputHistograms(0x0)
  ,fYRes(0x0)
  ,fPhiRes(0x0)
  ,fReconstructor(0x0)
  ,fDebugLevel(0)
  ,fDebugStream(0x0)
{
  AliCDBManager *cdb = AliCDBManager::Instance();
  cdb->SetDefaultStorage("local://$ALICE_ROOT");
  cdb->SetRun(0);
  TGeoManager::Import("geometry.root");
  fReconstructor = new AliTRDReconstructor();
  fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
  AliTRDtrackerV1::SetNTimeBins(24);

  DefineInput(0, TObjArray::Class());
  DefineOutput(0, TList::Class());
}

//________________________________________________________
void AliTRDtrackingResolution::ConnectInputData(Option_t *){
  fTrackObjects = dynamic_cast<TObjArray *>(GetInputData(0));
}

//________________________________________________________
void AliTRDtrackingResolution::CreateOutputObjects()
{
  // spatial resolution
  OpenFile(0, "RECREATE");
  fOutputHistograms = new TList();
  fYRes = new TH2I("fY", "", 21, -21., 21., 100, -.5, .5);
  fOutputHistograms->Add(fYRes);
  fPhiRes = new TH2I("fPhi", "", 21, -21., 21., 100, -10., 10.);
  fOutputHistograms->Add(fPhiRes);
}

//________________________________________________________
void AliTRDtrackingResolution::Exec(Option_t *)
{
  // spatial Resolution: res = pos_{Tracklet}(x = x_{Anode wire}) - pos_{TrackRef}(x = x_{Anode wire})
  // angular Resolution: res = Tracklet angle - TrackRef Angle

  Int_t nTrackInfos = fTrackObjects->GetEntriesFast();
  if(fDebugLevel>=2) printf("Number of Histograms: %d\n", fOutputHistograms->GetEntries());

  Double_t dy, dz;
  AliTRDtrackV1 *fTrack = 0x0;
  AliTRDtrackInfo *fInfo = 0x0;
  if(fDebugLevel>=2) printf("Number of TrackInfos: %d\n", nTrackInfos);
  for(Int_t iTI = 0; iTI < nTrackInfos; iTI++){
    // check if ESD and MC-Information are available
    if(!(fInfo = dynamic_cast<AliTRDtrackInfo *>(fTrackObjects->UncheckedAt(iTI)))) continue;
    if(!(fTrack = fInfo->GetTRDtrack())) continue;
    if(fDebugLevel>=3) printf("\tDoing track[%d] NTrackRefs[%d]\n", iTI, fInfo->GetNTrackRefs());

    if(fInfo->GetNTrackRefs() < 2) continue; 
  
    AliTRDseedV1 *fTracklet = 0x0;
    AliTrackReference *fTrackRefs[2];
    for(Int_t iplane = 0; iplane < kNLayers; iplane++){
      if(!(fTracklet = fTrack->GetTracklet(iplane))) continue;
      if(!fTracklet->IsOK()) continue;

      if(fDebugLevel>=4) printf("\t\tLy[%d] X0[%6.3f] Ncl[%d]\n", iplane, fTracklet->GetX0(), fTracklet->GetN());

      // check for 2 track ref where the radial position has a distance less than 3.7mm
      Int_t nFound = 0;
      memset(fTrackRefs, 0, sizeof(AliTrackReference*) * 2);
      AliTrackReference *tempTrackRef = 0;
      for(Int_t itr = 0; itr < fInfo->GetNTrackRefs(); itr++){
        if(fDebugLevel>=5) printf("nFound = %d\n", nFound);
        if(nFound >= 2) break;
        tempTrackRef = fInfo->GetTrackRef(itr);
        if(!tempTrackRef) continue;
        if(fDebugLevel>=5) printf("TrackRef %d: x = %f\n", itr, tempTrackRef->LocalX());
        if(fTracklet->GetX0() - tempTrackRef->LocalX() > 3.7) continue;
        if(tempTrackRef->LocalX() - fTracklet->GetX0() > 3.7) break;
        if(fDebugLevel>=5) printf("accepted\n");
        if(nFound == 1)
          if(fTrackRefs[0]->LocalX() >= tempTrackRef->LocalX()) continue;
        fTrackRefs[nFound++] = tempTrackRef;
      }
      if(fDebugLevel>=4) printf("\t\tFound track crossing [%d]\n", nFound);
      if(nFound < 2) continue;
      // We found 2 track refs for the tracklet, get y and z at the anode wire by a linear approximation


      // RESOLUTION
      Double_t dx = fTrackRefs[1]->LocalX() - fTrackRefs[0]->LocalX();
      Double_t dydx = (fTrackRefs[1]->LocalY() - fTrackRefs[0]->LocalY()) / dx;
      Double_t dzdx = (fTrackRefs[1]->Z() - fTrackRefs[0]->Z()) / dx;
      Double_t dx0 = fTrackRefs[1]->LocalX() - fTracklet->GetX0();
      Double_t ymc =  fTrackRefs[1]->LocalY() - dydx*dx0;
      Double_t zmc =  fTrackRefs[1]->Z() - dzdx*dx0;
      
      // recalculate tracklet based on the MC info
      fTracklet->SetYref(0, zmc);
      fTracklet->SetYref(1, dzdx);
      fTracklet->Fit();
      dy = fTracklet->GetYfit(0) - ymc;
      dz = fTracklet->GetZfit(0) - zmc;
      //res_y *= 100; // in mm
      Double_t momentum = fTrackRefs[0]->P();

//       Double_t phi     = fTrackRefs[0]->Phi();
//       Double_t theta   = fTrackRefs[0]->Theta();
      Double_t phi   = TMath::ATan(dydx);
      Double_t theta = TMath::ATan(dzdx);
      Double_t dphi   = TMath::ATan(fTracklet->GetYfit(1)) - phi;
      
      // Fill Histograms
      if(TMath::Abs(dx-3.7)<1.E-3){
        fYRes->Fill(phi*TMath::RadToDeg(), dy);
        fPhiRes->Fill(phi*TMath::RadToDeg(), dphi*TMath::RadToDeg());
      }      

      if(fDebugLevel>=4) printf("\t\tdx[%6.4f] dy[%6.4f] dz[%6.4f] dphi[%6.4f] \n", dx, dy, dz, dphi);
      
      // Fill Debug Tree
      if(fDebugLevel>=1){
        (*fDebugStream) << "Resolution"
          << "plane="	 	<< iplane
          << "p="       << momentum
          << "dx="      << dx
          << "dy="		  << dy
          << "dz="	 	  << dz
          << "phi="			<< phi
          << "theta="		<< theta
          << "dphi="		<< dphi
          << "\n";
      }
    }
    if(fTrack->GetNumberOfTracklets() < 6) continue;

    // RESIDUALS
    Int_t nc = 0;
    AliTrackPoint tr[AliTRDtrackV1::kMAXCLUSTERSPERTRACK], cl[AliTRDtrackV1::kMAXCLUSTERSPERTRACK];
    // prepare
    for(Int_t iplane = 0; iplane < kNLayers; iplane++){
      if(!(fTracklet = fTrack->GetTracklet(iplane))) continue;
      if(!fTracklet->IsOK()) continue;
      AliTRDcluster *c = 0x0;
      for(Int_t ic=AliTRDseed::knTimebins-1; ic>=0; ic--){
        if(!(c = fTracklet->GetClusters(ic))) continue;
        tr[nc].SetXYZ(c->GetX(), 0., 0.);
        cl[nc].SetXYZ(c->GetX(), c->GetY(), c->GetZ());
        nc++;
      }
    }
    AliTRDtrackerV1::FitRiemanTilt(fTrack, 0x0, kTRUE, nc, tr);
    for(Int_t ip=0; ip<nc; ip++){
      dy = cl[ip].GetY() - tr[ip].GetY();
      dz = cl[ip].GetZ() - tr[ip].GetZ();
      if(fDebugLevel>=1){
        (*fDebugStream) << "ResidualsRT"
          << "dy="		  << dy
          << "dz="	 	  << dz
/*          << "phi="			<< phi
          << "theta="		<< theta
          << "dphi="		<< dphi*/
          << "\n";
      }
    }

//     AliTRDtrackerV1::FitKalman(fTrack, 0x0, kFALSE, nc, tr);
//     for(Int_t ip=0; ip<nc; ip++){
//       dy = cl[ip].GetY() - tr[ip].GetY();
//       dz = cl[ip].GetZ() - tr[ip].GetZ();
//       if(fDebugLevel>=1){
//         (*fDebugStream) << "ResidualsKF"
//           << "dy="		  << dy
//           << "dz="	 	  << dz
// /*          << "phi="			<< phi
//           << "theta="		<< theta
//           << "dphi="		<< dphi*/
//           << "\n";
//       }
//     }    
  }
  PostData(0, fOutputHistograms);
}

//________________________________________________________
void AliTRDtrackingResolution::Terminate(Option_t *)
{
  if(fDebugStream) delete fDebugStream;

  //process distributions
  fOutputHistograms = dynamic_cast<TList*>(GetOutputData(0));
  if (!fOutputHistograms) {
    Printf("ERROR: list not available");
    return;
  }
  TH2I *h2 = 0x0;
  TH1D *h = 0x0;

  // y resolution
  TF1 *f = new TF1("f1", "gaus", -.5, .5);  
  h2 = (TH2I*)fOutputHistograms->At(0);
  TGraphErrors *gm = new TGraphErrors(h2->GetNbinsX());
  gm->SetNameTitle("meany", "Mean dy");
  TGraphErrors *gs = new TGraphErrors(h2->GetNbinsX());
  gs->SetNameTitle("sigmy", "Sigma y");
  for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
    Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
    f->SetParameter(1, 0.);f->SetParameter(2, 2.e-2);
    h = h2->ProjectionY("py", iphi, iphi);
    h->Fit(f, "q", "goff", -.5, .5);
    Int_t jphi = iphi -1;
    gm->SetPoint(jphi, phi, f->GetParameter(1));
    gm->SetPointError(jphi, 0., f->GetParError(1));
    gs->SetPoint(jphi, phi, f->GetParameter(2));
    gs->SetPointError(jphi, 0., f->GetParError(2));
  }
  fOutputHistograms->Add(gm);
  fOutputHistograms->Add(gs);

  // phi resolution
  h2 = (TH2I*)fOutputHistograms->At(1);
  gm = new TGraphErrors(h2->GetNbinsX());
  gm->SetNameTitle("meanphi", "Mean Phi");
  gs = new TGraphErrors(h2->GetNbinsX());
  gs->SetNameTitle("sigmphi", "Sigma Phi");
  for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
    Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
    f->SetParameter(1, 0.);f->SetParameter(2, 2.e-2);
    h = h2->ProjectionY("py", iphi, iphi);
    h->Fit(f, "q", "goff", -.5, .5);
    Int_t jphi = iphi -1;
    gm->SetPoint(jphi, phi, f->GetParameter(1));
    gm->SetPointError(jphi, 0., f->GetParError(1));
    gs->SetPoint(jphi, phi, f->GetParameter(2));
    gs->SetPointError(jphi, 0., f->GetParError(2));
  }
  fOutputHistograms->Add(gm);
  fOutputHistograms->Add(gs);


  delete f;
}

//________________________________________________________
void AliTRDtrackingResolution::SetDebugLevel(Int_t level){
  fDebugLevel = level;
  if(!fDebugLevel) return;
  if(fDebugStream) return;
  fDebugStream = new TTreeSRedirector("TRD.Resolution.root");
}

//________________________________________________________
void AliTRDtrackingResolution::SetRecoParam(AliTRDrecoParam *r)
{
  fReconstructor->SetRecoParam(r);
}
Commit	Line	Data
77203477	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: AliTRDtrackingResolution.cxx 27496 2008-07-22 08:35:45Z cblume $ */
	17
	18	////////////////////////////////////////////////////////////////////////////
	19	// //
	20	// Reconstruction QA //
	21	// //
	22	// Authors: //
	23	// Markus Fasel <M.Fasel@gsi.de> //
	24	// //
	25	////////////////////////////////////////////////////////////////////////////
	26
aaf47b30	27	#include <cstring>
	28
	29
77203477	30	#include <TObjArray.h>
77203477	31	#include <TList.h>
7102d1b1	32	#include <TH2.h>
	33	#include <TH1.h>
	34	#include <TF1.h>
77203477	35	#include <TProfile.h>
7102d1b1	36	#include <TGraphErrors.h>
77203477	37	#include <TMath.h>
aaf47b30	38	#include "TTreeStream.h"
aaf47b30	39	#include "TGeoManager.h"
77203477	40
77203477	41	#include "AliAnalysisManager.h"
77203477	42	#include "AliTrackReference.h"
aaf47b30	43	#include "AliTrackPointArray.h"
	44	#include "AliCDBManager.h"
	45
	46	#include "AliTRDcluster.h"
	47	#include "AliTRDseedV1.h"
	48	#include "AliTRDtrackV1.h"
	49	#include "AliTRDtrackerV1.h"
	50	#include "AliTRDReconstructor.h"
	51	#include "AliTRDrecoParam.h"
77203477	52
	53	#include "AliTRDtrackInfo/AliTRDtrackInfo.h"
	54	#include "AliTRDtrackingResolution.h"
	55
77203477	56	ClassImp(AliTRDtrackingResolution)
	57
	58	//________________________________________________________
	59	AliTRDtrackingResolution::AliTRDtrackingResolution(const char * name):
aaf47b30	60	AliAnalysisTask(name, "")
	61	,fTrackObjects(0x0)
	62	,fOutputHistograms(0x0)
	63	,fYRes(0x0)
	64	,fPhiRes(0x0)
	65	,fReconstructor(0x0)
	66	,fDebugLevel(0)
	67	,fDebugStream(0x0)
77203477	68	{
aaf47b30	69	AliCDBManager *cdb = AliCDBManager::Instance();
	70	cdb->SetDefaultStorage("local://$ALICE_ROOT");
	71	cdb->SetRun(0);
	72	TGeoManager::Import("geometry.root");
	73	fReconstructor = new AliTRDReconstructor();
	74	fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
	75	AliTRDtrackerV1::SetNTimeBins(24);
	76
77203477	77	DefineInput(0, TObjArray::Class());
	78	DefineOutput(0, TList::Class());
	79	}
	80
	81	//________________________________________________________
	82	void AliTRDtrackingResolution::ConnectInputData(Option_t *){
	83	fTrackObjects = dynamic_cast<TObjArray *>(GetInputData(0));
	84	}
	85
	86	//________________________________________________________
	87	void AliTRDtrackingResolution::CreateOutputObjects()
	88	{
	89	// spatial resolution
77203477	90	OpenFile(0, "RECREATE");
77203477	91	fOutputHistograms = new TList();
7102d1b1	92	fYRes = new TH2I("fY", "", 21, -21., 21., 100, -.5, .5);
	93	fOutputHistograms->Add(fYRes);
	94	fPhiRes = new TH2I("fPhi", "", 21, -21., 21., 100, -10., 10.);
	95	fOutputHistograms->Add(fPhiRes);
77203477	96	}
	97
	98	//________________________________________________________
7102d1b1	99	void AliTRDtrackingResolution::Exec(Option_t *)
7102d1b1	100	{
77203477	101	// spatial Resolution: res = pos_{Tracklet}(x = x_{Anode wire}) - pos_{TrackRef}(x = x_{Anode wire})
77203477	102	// angular Resolution: res = Tracklet angle - TrackRef Angle
7102d1b1	103
77203477	104	Int_t nTrackInfos = fTrackObjects->GetEntriesFast();
77203477	105	if(fDebugLevel>=2) printf("Number of Histograms: %d\n", fOutputHistograms->GetEntries());
aaf47b30	106
	107	Double_t dy, dz;
	108	AliTRDtrackV1 *fTrack = 0x0;
77203477	109	AliTRDtrackInfo *fInfo = 0x0;
	110	if(fDebugLevel>=2) printf("Number of TrackInfos: %d\n", nTrackInfos);
	111	for(Int_t iTI = 0; iTI < nTrackInfos; iTI++){
77203477	112	// check if ESD and MC-Information are available
7102d1b1	113	if(!(fInfo = dynamic_cast<AliTRDtrackInfo *>(fTrackObjects->UncheckedAt(iTI)))) continue;
aaf47b30	114	if(!(fTrack = fInfo->GetTRDtrack())) continue;
7102d1b1	115	if(fDebugLevel>=3) printf("\tDoing track[%d] NTrackRefs[%d]\n", iTI, fInfo->GetNTrackRefs());
	116
	117	if(fInfo->GetNTrackRefs() < 2) continue;
aaf47b30	118
aaf47b30	119	AliTRDseedV1 *fTracklet = 0x0;
77203477	120	AliTrackReference *fTrackRefs[2];
77203477	121	for(Int_t iplane = 0; iplane < kNLayers; iplane++){
aaf47b30	122	if(!(fTracklet = fTrack->GetTracklet(iplane))) continue;
	123	if(!fTracklet->IsOK()) continue;
	124
7102d1b1	125	if(fDebugLevel>=4) printf("\t\tLy[%d] X0[%6.3f] Ncl[%d]\n", iplane, fTracklet->GetX0(), fTracklet->GetN());
7102d1b1	126
77203477	127	// check for 2 track ref where the radial position has a distance less than 3.7mm
77203477	128	Int_t nFound = 0;
	129	memset(fTrackRefs, 0, sizeof(AliTrackReference) 2);
	130	AliTrackReference *tempTrackRef = 0;
	131	for(Int_t itr = 0; itr < fInfo->GetNTrackRefs(); itr++){
7102d1b1	132	if(fDebugLevel>=5) printf("nFound = %d\n", nFound);
77203477	133	if(nFound >= 2) break;
	134	tempTrackRef = fInfo->GetTrackRef(itr);
	135	if(!tempTrackRef) continue;
7102d1b1	136	if(fDebugLevel>=5) printf("TrackRef %d: x = %f\n", itr, tempTrackRef->LocalX());
77203477	137	if(fTracklet->GetX0() - tempTrackRef->LocalX() > 3.7) continue;
77203477	138	if(tempTrackRef->LocalX() - fTracklet->GetX0() > 3.7) break;
7102d1b1	139	if(fDebugLevel>=5) printf("accepted\n");
77203477	140	if(nFound == 1)
	141	if(fTrackRefs[0]->LocalX() >= tempTrackRef->LocalX()) continue;
	142	fTrackRefs[nFound++] = tempTrackRef;
	143	}
7102d1b1	144	if(fDebugLevel>=4) printf("\t\tFound track crossing [%d]\n", nFound);
77203477	145	if(nFound < 2) continue;
77203477	146	// We found 2 track refs for the tracklet, get y and z at the anode wire by a linear approximation
aaf47b30	147
	148
	149	// RESOLUTION
77203477	150	Double_t dx = fTrackRefs[1]->LocalX() - fTrackRefs[0]->LocalX();
	151	Double_t dydx = (fTrackRefs[1]->LocalY() - fTrackRefs[0]->LocalY()) / dx;
	152	Double_t dzdx = (fTrackRefs[1]->Z() - fTrackRefs[0]->Z()) / dx;
	153	Double_t dx0 = fTrackRefs[1]->LocalX() - fTracklet->GetX0();
	154	Double_t ymc = fTrackRefs[1]->LocalY() - dydx*dx0;
	155	Double_t zmc = fTrackRefs[1]->Z() - dzdx*dx0;
	156
aaf47b30	157	// recalculate tracklet based on the MC info
	158	fTracklet->SetYref(0, zmc);
	159	fTracklet->SetYref(1, dzdx);
	160	fTracklet->Fit();
	161	dy = fTracklet->GetYfit(0) - ymc;
	162	dz = fTracklet->GetZfit(0) - zmc;
77203477	163	//res_y *= 100; // in mm
77203477	164	Double_t momentum = fTrackRefs[0]->P();
77203477	165
	166	// Double_t phi = fTrackRefs[0]->Phi();
	167	// Double_t theta = fTrackRefs[0]->Theta();
	168	Double_t phi = TMath::ATan(dydx);
	169	Double_t theta = TMath::ATan(dzdx);
7102d1b1	170	Double_t dphi = TMath::ATan(fTracklet->GetYfit(1)) - phi;
77203477	171
7102d1b1	172	// Fill Histograms
	173	if(TMath::Abs(dx-3.7)<1.E-3){
	174	fYRes->Fill(phi*TMath::RadToDeg(), dy);
	175	fPhiRes->Fill(phiTMath::RadToDeg(), dphiTMath::RadToDeg());
aaf47b30	176	}
7102d1b1	177
7102d1b1	178	if(fDebugLevel>=4) printf("\t\tdx[%6.4f] dy[%6.4f] dz[%6.4f] dphi[%6.4f] \n", dx, dy, dz, dphi);
77203477	179
	180	// Fill Debug Tree
	181	if(fDebugLevel>=1){
	182	(*fDebugStream) << "Resolution"
	183	<< "plane=" << iplane
	184	<< "p=" << momentum
	185	<< "dx=" << dx
	186	<< "dy=" << dy
	187	<< "dz=" << dz
	188	<< "phi=" << phi
	189	<< "theta=" << theta
	190	<< "dphi=" << dphi
	191	<< "\n";
	192	}
	193	}
aaf47b30	194	if(fTrack->GetNumberOfTracklets() < 6) continue;
	195
	196	// RESIDUALS
	197	Int_t nc = 0;
	198	AliTrackPoint tr[AliTRDtrackV1::kMAXCLUSTERSPERTRACK], cl[AliTRDtrackV1::kMAXCLUSTERSPERTRACK];
	199	// prepare
	200	for(Int_t iplane = 0; iplane < kNLayers; iplane++){
	201	if(!(fTracklet = fTrack->GetTracklet(iplane))) continue;
	202	if(!fTracklet->IsOK()) continue;
	203	AliTRDcluster *c = 0x0;
	204	for(Int_t ic=AliTRDseed::knTimebins-1; ic>=0; ic--){
	205	if(!(c = fTracklet->GetClusters(ic))) continue;
	206	tr[nc].SetXYZ(c->GetX(), 0., 0.);
	207	cl[nc].SetXYZ(c->GetX(), c->GetY(), c->GetZ());
	208	nc++;
	209	}
	210	}
	211	AliTRDtrackerV1::FitRiemanTilt(fTrack, 0x0, kTRUE, nc, tr);
	212	for(Int_t ip=0; ip<nc; ip++){
	213	dy = cl[ip].GetY() - tr[ip].GetY();
	214	dz = cl[ip].GetZ() - tr[ip].GetZ();
	215	if(fDebugLevel>=1){
	216	(*fDebugStream) << "ResidualsRT"
	217	<< "dy=" << dy
	218	<< "dz=" << dz
	219	/* << "phi=" << phi
	220	<< "theta=" << theta
	221	<< "dphi=" << dphi*/
	222	<< "\n";
	223	}
	224	}
	225
	226	// AliTRDtrackerV1::FitKalman(fTrack, 0x0, kFALSE, nc, tr);
	227	// for(Int_t ip=0; ip<nc; ip++){
	228	// dy = cl[ip].GetY() - tr[ip].GetY();
	229	// dz = cl[ip].GetZ() - tr[ip].GetZ();
	230	// if(fDebugLevel>=1){
	231	// (*fDebugStream) << "ResidualsKF"
	232	// << "dy=" << dy
	233	// << "dz=" << dz
	234	// /* << "phi=" << phi
	235	// << "theta=" << theta
	236	// << "dphi=" << dphi*/
	237	// << "\n";
	238	// }
	239	// }
77203477	240	}
	241	PostData(0, fOutputHistograms);
	242	}
	243
	244	//________________________________________________________
7102d1b1	245	void AliTRDtrackingResolution::Terminate(Option_t *)
7102d1b1	246	{
77203477	247	if(fDebugStream) delete fDebugStream;
7102d1b1	248
	249	//process distributions
	250	fOutputHistograms = dynamic_cast<TList*>(GetOutputData(0));
	251	if (!fOutputHistograms) {
	252	Printf("ERROR: list not available");
	253	return;
	254	}
	255	TH2I *h2 = 0x0;
	256	TH1D *h = 0x0;
	257
	258	// y resolution
	259	TF1 *f = new TF1("f1", "gaus", -.5, .5);
	260	h2 = (TH2I*)fOutputHistograms->At(0);
	261	TGraphErrors *gm = new TGraphErrors(h2->GetNbinsX());
	262	gm->SetNameTitle("meany", "Mean dy");
	263	TGraphErrors *gs = new TGraphErrors(h2->GetNbinsX());
	264	gs->SetNameTitle("sigmy", "Sigma y");
	265	for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
	266	Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
	267	f->SetParameter(1, 0.);f->SetParameter(2, 2.e-2);
	268	h = h2->ProjectionY("py", iphi, iphi);
	269	h->Fit(f, "q", "goff", -.5, .5);
	270	Int_t jphi = iphi -1;
	271	gm->SetPoint(jphi, phi, f->GetParameter(1));
	272	gm->SetPointError(jphi, 0., f->GetParError(1));
	273	gs->SetPoint(jphi, phi, f->GetParameter(2));
	274	gs->SetPointError(jphi, 0., f->GetParError(2));
	275	}
	276	fOutputHistograms->Add(gm);
	277	fOutputHistograms->Add(gs);
	278
	279	// phi resolution
	280	h2 = (TH2I*)fOutputHistograms->At(1);
	281	gm = new TGraphErrors(h2->GetNbinsX());
	282	gm->SetNameTitle("meanphi", "Mean Phi");
	283	gs = new TGraphErrors(h2->GetNbinsX());
	284	gs->SetNameTitle("sigmphi", "Sigma Phi");
	285	for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
	286	Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
	287	f->SetParameter(1, 0.);f->SetParameter(2, 2.e-2);
	288	h = h2->ProjectionY("py", iphi, iphi);
	289	h->Fit(f, "q", "goff", -.5, .5);
	290	Int_t jphi = iphi -1;
	291	gm->SetPoint(jphi, phi, f->GetParameter(1));
	292	gm->SetPointError(jphi, 0., f->GetParError(1));
	293	gs->SetPoint(jphi, phi, f->GetParameter(2));
	294	gs->SetPointError(jphi, 0., f->GetParError(2));
	295	}
	296	fOutputHistograms->Add(gm);
	297	fOutputHistograms->Add(gs);
	298
	299
	300	delete f;
77203477	301	}
	302
	303	//________________________________________________________
	304	void AliTRDtrackingResolution::SetDebugLevel(Int_t level){
	305	fDebugLevel = level;
	306	if(!fDebugLevel) return;
	307	if(fDebugStream) return;
	308	fDebugStream = new TTreeSRedirector("TRD.Resolution.root");
	309	}
aaf47b30	310
	311	//________________________________________________________
	312	void AliTRDtrackingResolution::SetRecoParam(AliTRDrecoParam *r)
	313	{
	314	fReconstructor->SetRecoParam(r);
	315	}