[u/mrichter/AliRoot.git] / TRD / AliTRDmcmTracklet.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.                  *
 **************************************************************************/

///////////////////////////////////////////////////////
//                                                   //
//                                                   //
//  Tracklet object created in the local tracking    //
//                                                   //
//                                                   //
///////////////////////////////////////////////////////

#include <TGraph.h>
#include <TMath.h>
#include <TF1.h>

#include "AliTRDcalibDB.h"
#include "AliTRDCommonParam.h"
#include "AliTRDpadPlane.h"
#include "AliTRDgeometry.h"

#include "AliTRDmcmTracklet.h"

ClassImp(AliTRDmcmTracklet)

//_____________________________________________________________________________
AliTRDmcmTracklet::AliTRDmcmTracklet() 
{
  //
  // AliTRDmcmTracklet default constructor
  //

  fDetector = -1;
  fRow      = -1;

  for (Int_t time = 0; time < kNtimeBins; time++) {
    for (Int_t icl = 0; icl < kNclsPads; icl++) {
      fADC[time][icl] = 0;
    }
    for (Int_t it = 0; it < kNdict; it++) {
      fTrack[time][it] = -1;
    }
    fTime[time]   = 0;
    fCol[time]    = 0;
  }

  fNclusters  =  0;
  fN          =  0;
  fTrackLabel = -1;

  fGPos = 0;
  fGAmp = 0;

  fSlope  = 0.0;
  fOffset = 0.0;
  fTime0  = 0.0;
  fRowz   = 0.0;
  fPt     = 0.0;
  fdQdl   = 0.0;

}

//_____________________________________________________________________________
AliTRDmcmTracklet::AliTRDmcmTracklet(Int_t det, Int_t row, Int_t n) 
{
  //
  // AliTRDmcmTracklet default constructor
  //

  fDetector = det;
  fRow      = row;

  for (Int_t time = 0; time < kNtimeBins; time++) {
    for (Int_t icl = 0; icl < kNclsPads; icl++) {
      fADC[time][icl] = 0;
    }
    for (Int_t it = 0; it < kNdict; it++) {
      fTrack[time][it] = -1;
    }
    fTime[time]   = 0;
    fCol[time]    = 0;
  }

  fNclusters = 0;

  fN = n;

  fTrackLabel = -1;

  fGPos = new TGraph(0);
  fGAmp = new TGraph(0);

  fSlope  = 0.0;
  fOffset = 0.0;
  fTime0  = 0.0;
  fRowz   = 0.0;
  fPt     = 0.0;
  fdQdl   = 0.0;

}

//_____________________________________________________________________________
AliTRDmcmTracklet::~AliTRDmcmTracklet() 
{
  //
  // AliTRDmcmTracklet destructor
  //

  if (fGPos != 0) delete fGPos;
  if (fGAmp != 0) delete fGAmp;

}

//_____________________________________________________________________________
AliTRDmcmTracklet &AliTRDmcmTracklet::operator=(const AliTRDmcmTracklet &t)
{
  //
  // assignment operator
  //

  if (this != &t) ((AliTRDmcmTracklet &) t).Copy(*this); 
  return *this;

}

//_____________________________________________________________________________
void AliTRDmcmTracklet::Copy(TObject &t) const
{
  //
  // copy function
  //

  ((AliTRDmcmTracklet &) t).fDetector    = fDetector;
  ((AliTRDmcmTracklet &) t).fRow         = fRow;
  ((AliTRDmcmTracklet &) t).fTrackLabel  = fTrackLabel;
  ((AliTRDmcmTracklet &) t).fNclusters   = fNclusters;
  ((AliTRDmcmTracklet &) t).fN           = fN;
  ((AliTRDmcmTracklet &) t).fGPos        = NULL;
  ((AliTRDmcmTracklet &) t).fGAmp        = NULL;
  ((AliTRDmcmTracklet &) t).fTime0       = fTime0;
  ((AliTRDmcmTracklet &) t).fRowz        = fRowz;
  ((AliTRDmcmTracklet &) t).fSlope       = fSlope;
  ((AliTRDmcmTracklet &) t).fOffset      = fOffset;
  ((AliTRDmcmTracklet &) t).fPt          = fPt;
  ((AliTRDmcmTracklet &) t).fdQdl        = fdQdl;

  for (Int_t time = 0; time < kNtimeBins; time++) {
    for (Int_t icl = 0; icl < kNclsPads; icl++) {
      ((AliTRDmcmTracklet &) t).fADC[time][icl] = 0;
    }
    for (Int_t it = 0; it < kNdict; it++) {
      ((AliTRDmcmTracklet &) t).fTrack[time][it] = -1;
    }
    ((AliTRDmcmTracklet &) t).fTime[time]   = 0;
    ((AliTRDmcmTracklet &) t).fCol[time]    = 0;
  }

}

//_____________________________________________________________________________
void AliTRDmcmTracklet::Reset() 
{
  //
  // Reset the tracklet information
  //

  fDetector = -1;
  fRow      = -1;

  for (Int_t time = 0; time < kNtimeBins; time++) {
    for (Int_t icl = 0; icl < kNclsPads; icl++) {
      fADC[time][icl] = 0;
    }
    for (Int_t it = 0; it < kNdict; it++) {
      fTrack[time][it] = -1;
    }
    fTime[time]   = 0;
    fCol[time]    = 0;
  }

  fNclusters  =  0;
  fN          =  0;
  fTrackLabel = -1;

  fGPos->Set(0);
  fGAmp->Set(0);

  fSlope  = 0.0;
  fOffset = 0.0;
  fTime0  = 0.0;
  fRowz   = 0.0;
  fPt     = 0.0;
  fdQdl   = 0.0;

}

//_____________________________________________________________________________
void AliTRDmcmTracklet::AddCluster(Int_t icol, Int_t itb, Float_t *adc, Int_t *track) 
{
  //
  // Add a cluster to the tracklet
  //
 
  if (fNclusters >= kNtimeBins) return;

  for (Int_t icl = 0; icl < kNclsPads; icl++) {
    //fADC[fNclusters][icl] = (Int_t)adc[icl]; 
    fADC[fNclusters][icl] = adc[icl]; 
  }

  fTrack[fNclusters][0] = track[0];
  fTrack[fNclusters][1] = track[1];
  fTrack[fNclusters][2] = track[2];
  fTime[fNclusters] = itb;
  fCol[fNclusters] = icol;

  fNclusters++;

}

//_____________________________________________________________________________
void AliTRDmcmTracklet::MakeTrackletGraph(AliTRDgeometry *geo, Float_t field) 
{
  //
  // Tracklet graph of positions (global coordinates, rotated [cm])
  //
  
  if (!geo) {
    Error("MakeTrackletGraph","No geometry.");
    return;
  }
  
  AliTRDCommonParam* commonParam = AliTRDCommonParam::Instance();
  if (!commonParam)
  {
    Error("MakeTrackletGraph","No common params.");
    return;
  }

  AliTRDcalibDB* calibration = AliTRDcalibDB::Instance();
  if (!calibration)
  {
    Error("MakeTrackletGraph","No instance of AliTRDcalibDB.");
    return;
  }

  Int_t iplan, icham;

  iplan = geo->GetPlane(fDetector);
  icham = geo->GetChamber(fDetector);

  AliTRDpadPlane *padPlane = commonParam->GetPadPlane(iplan,icham);

  Float_t samplFreq = calibration->GetSamplingFrequency();

  Int_t time, col;
  Float_t amp[3];
  Float_t xpos, ypos, xzero, colSize, timeBinSize;
  Float_t vDrift, omegaTau, lorentzAngle, thetaSlope;
  Float_t tiltingAngle;
  Int_t npg = 0;

  for (Int_t icl = 0; icl < fNclusters; icl++) {

    time   = GetClusterTime(icl);

    amp[0] = GetClusterADC(icl)[0];
    amp[1] = GetClusterADC(icl)[1];
    amp[2] = GetClusterADC(icl)[2];

    col    = GetClusterCol(icl);

    if (amp[0] < 0.0 || amp[1] < 0.0 || amp[2] < 0.0) continue;

    ypos = GetClusY(amp,iplan);

    colSize = padPlane->GetColSize(col);
    vDrift = calibration->GetVdrift(fDetector,col,fRow);
    timeBinSize = vDrift/samplFreq;
    
    // From v4-03-Release to HEAD28Mar06 the sign has changed from "-" to "+" 
    // due to a change in the digitizer
    omegaTau = +TMath::Sign(1.0,(Double_t)field)*GetOmegaTau(vDrift,TMath::Abs(field));
    lorentzAngle = TMath::ATan(omegaTau)*180.0/TMath::Pi();
    
    xpos = (time+0.5) * timeBinSize;
    xpos = geo->GetTime0(iplan) - xpos;

    ypos = padPlane->GetColPos(col) - (ypos + 0.5) * colSize;

    // ExB correction
    xzero = geo->GetTime0(iplan);
    ypos = ypos + (xpos-xzero) * omegaTau;

    // tilted pads correction
    thetaSlope = - padPlane->GetRowPos(fRow)/geo->GetTime0(iplan);
    tiltingAngle = padPlane->GetTiltingAngle()/180.0*TMath::Pi();
    ypos = ypos - (xpos-xzero) * thetaSlope * TMath::Sin(tiltingAngle);

    fGPos->SetPoint(npg,(Double_t)xpos,(Double_t)ypos);
    npg++;

  }

  fGPos->Set(npg);
  
  fTime0 = geo->GetTime0(iplan) - AliTRDgeometry::CdrHght() - 0.5*AliTRDgeometry::CamHght();
  fRowz = 0.5*(padPlane->GetRowPos(fRow) + padPlane->GetRowPos(fRow+1));

  Double_t xMin = 0, xMax = 0, x, y;
  fGPos->GetPoint(0    ,x,y); xMax = x + 0.1;
  fGPos->GetPoint(npg-1,x,y); xMin = x - 0.1;
  
  TF1 *line = new TF1("line","[0]+x*[1]",xMin,xMax);
  fGPos->Fit(line,"WRQ0");

  fOffset = line->Eval(fTime0);
  fSlope  = TMath::ATan(line->GetParameter(1))*180.0/TMath::Pi();

  line->Delete();
  
  Float_t fx = fTime0;
  Float_t fy = fOffset;
  
  Float_t infSlope = TMath::ATan(fy/fx)/TMath::Pi()*180.0;    
  Float_t alpha = fSlope - infSlope;
  Float_t r = TMath::Sqrt(fx*fx + fy*fy)/(2.0*TMath::Sin(alpha/180.0*TMath::Pi()));

  fPt = 0.3 * field * 0.01 * r;
  
  return;

}

//_____________________________________________________________________________
void AliTRDmcmTracklet::MakeClusAmpGraph() 
{
  //
  // Tracklet graph of cluster charges
  //

  Int_t time;
  Float_t amp[3];
  Int_t npg = 0;
  fdQdl = 0.0;
  for (Int_t icl = 0; icl < fNclusters; icl++) {

    time   = GetClusterTime(icl);

    amp[0] = GetClusterADC(icl)[0];
    amp[1] = GetClusterADC(icl)[1];
    amp[2] = GetClusterADC(icl)[2];

    fGAmp->SetPoint(npg,(Double_t)(time+0.5),(Double_t)(amp[0]+amp[1]+amp[2]));
    npg++;

    fdQdl += amp[0]+amp[1]+amp[2];

  }

  fGAmp->Set(npg);

  fdQdl /= (Float_t)npg;

  return;

}

//_____________________________________________________________________________
Float_t AliTRDmcmTracklet::GetClusY(Float_t *adc, Int_t pla) const
{
  //
  // Cluster position in the phi direction in pad units (relative to the pad border)
  //

  Float_t ypos = 0.0;

  Float_t a0 = adc[0];
  Float_t a1 = adc[1];
  Float_t a2 = adc[2];

  Float_t t1, t2, w1 ,w2;

  Float_t w = 1.0;  // pad units

  Float_t sigma = 0.0;

  switch(pla) {
  case 0:
    sigma = 0.515; break;
  case 1:
    sigma = 0.501; break;
  case 2:
    sigma = 0.491; break;
  case 3:
    sigma = 0.481; break;
  case 4:
    sigma = 0.471; break;
  case 5:
    sigma = 0.463; break;
  default:
    Error("GetClusY","Wrong plane number.");
    return 0.0;
  }

  sigma *= w;

  t1 = 0.0;
  w1 = 0.0;
  if( a0 > 0 ) {
    w1 = a0*a0;
    //w1 = a0;
    t1 = w1*((sigma*sigma)/w*TMath::Log(a1/a0)-0.5*w);
  }
  t2 = 0.0;
  w2 = 0.0;
  if( a2 > 0 ) {
    w2 = a2*a2;
    //w2 = a2;
    t2 = w2*((sigma*sigma)/w*TMath::Log(a2/a1)+0.5*w);
  }

  ypos = w*(t1+t2)/(w1+w2);  // range: -0.5*w ... +0.5*w

  return ypos;

}

//_____________________________________________________________________________
void AliTRDmcmTracklet::CookLabel(Float_t frac) 
{
  //
  // Cook the track label from cluster labels
  //

  const Int_t kMaxTracks = 10;
  Int_t trackLabel[kMaxTracks];
  Int_t trackCount[kMaxTracks];
  for (Int_t it = 0; it < kMaxTracks; it++) {
    trackLabel[it] = -1;
    trackCount[it] = 0;
  }

  Bool_t counted;
  Int_t label, nTracks = 0;
  for (Int_t icl = 0; icl < fNclusters; icl++) {

    for (Int_t id = 0; id < kNdict; id++) {

      if (fTrack[icl][id] == -1) continue;

      label = fTrack[icl][id];

      counted = kFALSE;
      for (Int_t it = 0; it < nTracks; it++) {
	if (label == trackLabel[it]) {
	  trackCount[it]++;
	  counted = kTRUE;
	  break;
	}
      }
      if (!counted) {
	trackLabel[nTracks] = label;
	trackCount[nTracks]++;
	nTracks++;
	if (nTracks == kMaxTracks) {
	  Warning("CookLabel","Too many tracks for this tracklet.");
	  nTracks--;
	  break;
	}
      }

    }

  }

  for (Int_t it = 0; it < kMaxTracks; it++) {
    if (trackCount[it] >= (Int_t)(frac*fNclusters)) {
      fTrackLabel = trackLabel[it];
      break;
    }
  }

}

//_____________________________________________________________________________
Float_t AliTRDmcmTracklet::GetOmegaTau(Float_t vdrift, Float_t field) const
{
  //
  // Returns omega*tau (tan(Lorentz-angle)) for a given drift velocity <vd> 
  // and a B-field <b> for Xe/CO2 (15%).
  // The values are according to a GARFIELD simulation.
  //

  //
  // Copy of the "AliTRDcalibDB" function, taking as argument the magnetic field too
  //
  
  const Int_t kNb = 5;
  Float_t p0[kNb] = {  0.004810,  0.007412,  0.010252,  0.013409,  0.016888 };
  Float_t p1[kNb] = {  0.054875,  0.081534,  0.107333,  0.131983,  0.155455 };
  Float_t p2[kNb] = { -0.008682, -0.012896, -0.016987, -0.020880, -0.024623 };
  Float_t p3[kNb] = {  0.000155,  0.000238,  0.000330,  0.000428,  0.000541 };

  Int_t ib = ((Int_t) (10 * (field - 0.15)));
  ib       = TMath::Max(  0,ib);
  ib       = TMath::Min(kNb,ib);

  Float_t alphaL = p0[ib] 
      + p1[ib] * vdrift
      + p2[ib] * vdrift*vdrift
      + p3[ib] * vdrift*vdrift*vdrift;

  return TMath::Tan(alphaL);

}
Commit	Line	Data
0ee00e25	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
e3b2b5e5	16	///////////////////////////////////////////////////////
	17	// //
	18	// //
	19	// Tracklet object created in the local tracking //
	20	// //
	21	// //
	22	///////////////////////////////////////////////////////
	23
0ee00e25	24	#include <TGraph.h>
	25	#include <TMath.h>
	26	#include <TF1.h>
	27
	28	#include "AliTRDcalibDB.h"
	29	#include "AliTRDCommonParam.h"
	30	#include "AliTRDpadPlane.h"
	31	#include "AliTRDgeometry.h"
	32
	33	#include "AliTRDmcmTracklet.h"
	34
	35	ClassImp(AliTRDmcmTracklet)
	36
	37	//_____________________________________________________________________________
	38	AliTRDmcmTracklet::AliTRDmcmTracklet()
	39	{
0ee00e25	40	//
	41	// AliTRDmcmTracklet default constructor
	42	//
	43
	44	fDetector = -1;
	45	fRow = -1;
	46
	47	for (Int_t time = 0; time < kNtimeBins; time++) {
	48	for (Int_t icl = 0; icl < kNclsPads; icl++) {
	49	fADC[time][icl] = 0;
	50	}
	51	for (Int_t it = 0; it < kNdict; it++) {
	52	fTrack[time][it] = -1;
	53	}
	54	fTime[time] = 0;
	55	fCol[time] = 0;
	56	}
	57
	58	fNclusters = 0;
	59	fN = 0;
	60	fTrackLabel = -1;
	61
	62	fGPos = 0;
	63	fGAmp = 0;
	64
	65	fSlope = 0.0;
	66	fOffset = 0.0;
	67	fTime0 = 0.0;
	68	fRowz = 0.0;
	69	fPt = 0.0;
	70	fdQdl = 0.0;
	71
	72	}
	73
	74	//_____________________________________________________________________________
	75	AliTRDmcmTracklet::AliTRDmcmTracklet(Int_t det, Int_t row, Int_t n)
	76	{
0ee00e25	77	//
	78	// AliTRDmcmTracklet default constructor
	79	//
	80
	81	fDetector = det;
	82	fRow = row;
	83
	84	for (Int_t time = 0; time < kNtimeBins; time++) {
	85	for (Int_t icl = 0; icl < kNclsPads; icl++) {
	86	fADC[time][icl] = 0;
	87	}
	88	for (Int_t it = 0; it < kNdict; it++) {
	89	fTrack[time][it] = -1;
	90	}
	91	fTime[time] = 0;
	92	fCol[time] = 0;
	93	}
	94
	95	fNclusters = 0;
	96
	97	fN = n;
	98
	99	fTrackLabel = -1;
	100
	101	fGPos = new TGraph(0);
	102	fGAmp = new TGraph(0);
	103
	104	fSlope = 0.0;
	105	fOffset = 0.0;
	106	fTime0 = 0.0;
	107	fRowz = 0.0;
	108	fPt = 0.0;
	109	fdQdl = 0.0;
	110
	111	}
	112
	113	//_____________________________________________________________________________
	114	AliTRDmcmTracklet::~AliTRDmcmTracklet()
	115	{
0ee00e25	116	//
	117	// AliTRDmcmTracklet destructor
	118	//
3fcb9908	119
	120	if (fGPos != 0) delete fGPos;
	121	if (fGAmp != 0) delete fGAmp;
	122
e3b2b5e5	123	}
	124
	125	//_____________________________________________________________________________
	126	AliTRDmcmTracklet &AliTRDmcmTracklet::operator=(const AliTRDmcmTracklet &t)
	127	{
	128	//
	129	// assignment operator
	130	//
	131
	132	if (this != &t) ((AliTRDmcmTracklet &) t).Copy(*this);
	133	return *this;
	134
	135	}
	136
	137	//_____________________________________________________________________________
	138	void AliTRDmcmTracklet::Copy(TObject &t) const
	139	{
	140	//
	141	// copy function
	142	//
	143
	144	((AliTRDmcmTracklet &) t).fDetector = fDetector;
	145	((AliTRDmcmTracklet &) t).fRow = fRow;
	146	((AliTRDmcmTracklet &) t).fTrackLabel = fTrackLabel;
	147	((AliTRDmcmTracklet &) t).fNclusters = fNclusters;
	148	((AliTRDmcmTracklet &) t).fN = fN;
	149	((AliTRDmcmTracklet &) t).fGPos = NULL;
	150	((AliTRDmcmTracklet &) t).fGAmp = NULL;
	151	((AliTRDmcmTracklet &) t).fTime0 = fTime0;
	152	((AliTRDmcmTracklet &) t).fRowz = fRowz;
	153	((AliTRDmcmTracklet &) t).fSlope = fSlope;
	154	((AliTRDmcmTracklet &) t).fOffset = fOffset;
	155	((AliTRDmcmTracklet &) t).fPt = fPt;
	156	((AliTRDmcmTracklet &) t).fdQdl = fdQdl;
	157
	158	for (Int_t time = 0; time < kNtimeBins; time++) {
	159	for (Int_t icl = 0; icl < kNclsPads; icl++) {
	160	((AliTRDmcmTracklet &) t).fADC[time][icl] = 0;
	161	}
	162	for (Int_t it = 0; it < kNdict; it++) {
	163	((AliTRDmcmTracklet &) t).fTrack[time][it] = -1;
	164	}
	165	((AliTRDmcmTracklet &) t).fTime[time] = 0;
	166	((AliTRDmcmTracklet &) t).fCol[time] = 0;
	167	}
	168
	169	}
	170
	171	//_____________________________________________________________________________
3fcb9908	172	void AliTRDmcmTracklet::Reset()
e3b2b5e5	173	{
e3b2b5e5	174	//
3fcb9908	175	// Reset the tracklet information
	176	//
	177
	178	fDetector = -1;
	179	fRow = -1;
	180
	181	for (Int_t time = 0; time < kNtimeBins; time++) {
	182	for (Int_t icl = 0; icl < kNclsPads; icl++) {
	183	fADC[time][icl] = 0;
	184	}
	185	for (Int_t it = 0; it < kNdict; it++) {
	186	fTrack[time][it] = -1;
	187	}
	188	fTime[time] = 0;
	189	fCol[time] = 0;
	190	}
	191
	192	fNclusters = 0;
	193	fN = 0;
	194	fTrackLabel = -1;
	195
	196	fGPos->Set(0);
	197	fGAmp->Set(0);
	198
	199	fSlope = 0.0;
	200	fOffset = 0.0;
	201	fTime0 = 0.0;
	202	fRowz = 0.0;
	203	fPt = 0.0;
	204	fdQdl = 0.0;
	205
	206	}
	207
	208	//_____________________________________________________________________________
	209	void AliTRDmcmTracklet::AddCluster(Int_t icol, Int_t itb, Float_t adc, Int_t track)
	210	{
e3b2b5e5	211	//
3fcb9908	212	// Add a cluster to the tracklet
e3b2b5e5	213	//
	214
	215	if (fNclusters >= kNtimeBins) return;
	216
	217	for (Int_t icl = 0; icl < kNclsPads; icl++) {
	218	//fADC[fNclusters][icl] = (Int_t)adc[icl];
	219	fADC[fNclusters][icl] = adc[icl];
	220	}
	221
	222	fTrack[fNclusters][0] = track[0];
	223	fTrack[fNclusters][1] = track[1];
	224	fTrack[fNclusters][2] = track[2];
	225	fTime[fNclusters] = itb;
	226	fCol[fNclusters] = icol;
	227
	228	fNclusters++;
0ee00e25	229
	230	}
	231
	232	//_____________________________________________________________________________
	233	void AliTRDmcmTracklet::MakeTrackletGraph(AliTRDgeometry *geo, Float_t field)
	234	{
0ee00e25	235	//
	236	// Tracklet graph of positions (global coordinates, rotated [cm])
	237	//
	238
	239	if (!geo) {
	240	Error("MakeTrackletGraph","No geometry.");
	241	return;
	242	}
	243
	244	AliTRDCommonParam* commonParam = AliTRDCommonParam::Instance();
	245	if (!commonParam)
	246	{
	247	Error("MakeTrackletGraph","No common params.");
	248	return;
	249	}
	250
	251	AliTRDcalibDB* calibration = AliTRDcalibDB::Instance();
	252	if (!calibration)
	253	{
	254	Error("MakeTrackletGraph","No instance of AliTRDcalibDB.");
	255	return;
	256	}
	257
	258	Int_t iplan, icham;
	259
	260	iplan = geo->GetPlane(fDetector);
	261	icham = geo->GetChamber(fDetector);
	262
	263	AliTRDpadPlane *padPlane = commonParam->GetPadPlane(iplan,icham);
	264
e3b2b5e5	265	Float_t samplFreq = calibration->GetSamplingFrequency();
0ee00e25	266
	267	Int_t time, col;
	268	Float_t amp[3];
e3b2b5e5	269	Float_t xpos, ypos, xzero, colSize, timeBinSize;
	270	Float_t vDrift, omegaTau, lorentzAngle, thetaSlope;
	271	Float_t tiltingAngle;
0ee00e25	272	Int_t npg = 0;
	273
	274	for (Int_t icl = 0; icl < fNclusters; icl++) {
	275
	276	time = GetClusterTime(icl);
	277
	278	amp[0] = GetClusterADC(icl)[0];
	279	amp[1] = GetClusterADC(icl)[1];
	280	amp[2] = GetClusterADC(icl)[2];
	281
	282	col = GetClusterCol(icl);
	283
	284	if (amp[0] < 0.0 \|\| amp[1] < 0.0 \|\| amp[2] < 0.0) continue;
	285
e3b2b5e5	286	ypos = GetClusY(amp,iplan);
0ee00e25	287
e3b2b5e5	288	colSize = padPlane->GetColSize(col);
0ee00e25	289	vDrift = calibration->GetVdrift(fDetector,col,fRow);
e3b2b5e5	290	timeBinSize = vDrift/samplFreq;
0ee00e25	291
	292	// From v4-03-Release to HEAD28Mar06 the sign has changed from "-" to "+"
	293	// due to a change in the digitizer
e3b2b5e5	294	omegaTau = +TMath::Sign(1.0,(Double_t)field)*GetOmegaTau(vDrift,TMath::Abs(field));
e3b2b5e5	295	lorentzAngle = TMath::ATan(omegaTau)*180.0/TMath::Pi();
0ee00e25	296
e3b2b5e5	297	xpos = (time+0.5) * timeBinSize;
e3b2b5e5	298	xpos = geo->GetTime0(iplan) - xpos;
0ee00e25	299
e3b2b5e5	300	ypos = padPlane->GetColPos(col) - (ypos + 0.5) * colSize;
0ee00e25	301
0ee00e25	302	// ExB correction
e3b2b5e5	303	xzero = geo->GetTime0(iplan);
e3b2b5e5	304	ypos = ypos + (xpos-xzero) * omegaTau;
0ee00e25	305
0ee00e25	306	// tilted pads correction
e3b2b5e5	307	thetaSlope = - padPlane->GetRowPos(fRow)/geo->GetTime0(iplan);
	308	tiltingAngle = padPlane->GetTiltingAngle()/180.0*TMath::Pi();
	309	ypos = ypos - (xpos-xzero) * thetaSlope * TMath::Sin(tiltingAngle);
0ee00e25	310
e3b2b5e5	311	fGPos->SetPoint(npg,(Double_t)xpos,(Double_t)ypos);
0ee00e25	312	npg++;
	313
	314	}
	315
	316	fGPos->Set(npg);
	317
	318	fTime0 = geo->GetTime0(iplan) - AliTRDgeometry::CdrHght() - 0.5*AliTRDgeometry::CamHght();
3fcb9908	319	fRowz = 0.5*(padPlane->GetRowPos(fRow) + padPlane->GetRowPos(fRow+1));
0ee00e25	320
	321	Double_t xMin = 0, xMax = 0, x, y;
	322	fGPos->GetPoint(0 ,x,y); xMax = x + 0.1;
	323	fGPos->GetPoint(npg-1,x,y); xMin = x - 0.1;
	324
	325	TF1 line = new TF1("line","[0]+x[1]",xMin,xMax);
	326	fGPos->Fit(line,"WRQ0");
	327
	328	fOffset = line->Eval(fTime0);
	329	fSlope = TMath::ATan(line->GetParameter(1))*180.0/TMath::Pi();
	330
	331	line->Delete();
	332
e3b2b5e5	333	Float_t fx = fTime0;
e3b2b5e5	334	Float_t fy = fOffset;
0ee00e25	335
e3b2b5e5	336	Float_t infSlope = TMath::ATan(fy/fx)/TMath::Pi()*180.0;
0ee00e25	337	Float_t alpha = fSlope - infSlope;
e3b2b5e5	338	Float_t r = TMath::Sqrt(fxfx + fyfy)/(2.0TMath::Sin(alpha/180.0TMath::Pi()));
0ee00e25	339
e3b2b5e5	340	fPt = 0.3 * field * 0.01 * r;
0ee00e25	341
	342	return;
	343
	344	}
	345
	346	//_____________________________________________________________________________
	347	void AliTRDmcmTracklet::MakeClusAmpGraph()
	348	{
	349	//
	350	// Tracklet graph of cluster charges
	351	//
	352
0ee00e25	353	Int_t time;
	354	Float_t amp[3];
	355	Int_t npg = 0;
	356	fdQdl = 0.0;
	357	for (Int_t icl = 0; icl < fNclusters; icl++) {
	358
	359	time = GetClusterTime(icl);
	360
	361	amp[0] = GetClusterADC(icl)[0];
	362	amp[1] = GetClusterADC(icl)[1];
	363	amp[2] = GetClusterADC(icl)[2];
	364
	365	fGAmp->SetPoint(npg,(Double_t)(time+0.5),(Double_t)(amp[0]+amp[1]+amp[2]));
	366	npg++;
	367
	368	fdQdl += amp[0]+amp[1]+amp[2];
	369
	370	}
	371
	372	fGAmp->Set(npg);
	373
	374	fdQdl /= (Float_t)npg;
	375
	376	return;
	377
	378	}
	379
	380	//_____________________________________________________________________________
e3b2b5e5	381	Float_t AliTRDmcmTracklet::GetClusY(Float_t *adc, Int_t pla) const
0ee00e25	382	{
	383	//
	384	// Cluster position in the phi direction in pad units (relative to the pad border)
	385	//
	386
e3b2b5e5	387	Float_t ypos = 0.0;
0ee00e25	388
e3b2b5e5	389	Float_t a0 = adc[0];
	390	Float_t a1 = adc[1];
	391	Float_t a2 = adc[2];
0ee00e25	392
e3b2b5e5	393	Float_t t1, t2, w1 ,w2;
0ee00e25	394
e3b2b5e5	395	Float_t w = 1.0; // pad units
0ee00e25	396
e3b2b5e5	397	Float_t sigma = 0.0;
0ee00e25	398
	399	switch(pla) {
	400	case 0:
e3b2b5e5	401	sigma = 0.515; break;
0ee00e25	402	case 1:
e3b2b5e5	403	sigma = 0.501; break;
0ee00e25	404	case 2:
e3b2b5e5	405	sigma = 0.491; break;
0ee00e25	406	case 3:
e3b2b5e5	407	sigma = 0.481; break;
0ee00e25	408	case 4:
e3b2b5e5	409	sigma = 0.471; break;
0ee00e25	410	case 5:
e3b2b5e5	411	sigma = 0.463; break;
0ee00e25	412	default:
	413	Error("GetClusY","Wrong plane number.");
	414	return 0.0;
	415	}
	416
e3b2b5e5	417	sigma *= w;
0ee00e25	418
e3b2b5e5	419	t1 = 0.0;
	420	w1 = 0.0;
	421	if( a0 > 0 ) {
	422	w1 = a0*a0;
	423	//w1 = a0;
	424	t1 = w1((sigmasigma)/wTMath::Log(a1/a0)-0.5w);
0ee00e25	425	}
e3b2b5e5	426	t2 = 0.0;
	427	w2 = 0.0;
	428	if( a2 > 0 ) {
	429	w2 = a2*a2;
	430	//w2 = a2;
	431	t2 = w2((sigmasigma)/wTMath::Log(a2/a1)+0.5w);
0ee00e25	432	}
0ee00e25	433
e3b2b5e5	434	ypos = w(t1+t2)/(w1+w2); // range: -0.5w ... +0.5*w
0ee00e25	435
e3b2b5e5	436	return ypos;
0ee00e25	437
	438	}
	439
	440	//_____________________________________________________________________________
	441	void AliTRDmcmTracklet::CookLabel(Float_t frac)
	442	{
	443	//
	444	// Cook the track label from cluster labels
	445	//
	446
	447	const Int_t kMaxTracks = 10;
	448	Int_t trackLabel[kMaxTracks];
	449	Int_t trackCount[kMaxTracks];
	450	for (Int_t it = 0; it < kMaxTracks; it++) {
	451	trackLabel[it] = -1;
	452	trackCount[it] = 0;
	453	}
	454
	455	Bool_t counted;
	456	Int_t label, nTracks = 0;
	457	for (Int_t icl = 0; icl < fNclusters; icl++) {
	458
	459	for (Int_t id = 0; id < kNdict; id++) {
	460
	461	if (fTrack[icl][id] == -1) continue;
	462
	463	label = fTrack[icl][id];
	464
	465	counted = kFALSE;
	466	for (Int_t it = 0; it < nTracks; it++) {
	467	if (label == trackLabel[it]) {
	468	trackCount[it]++;
	469	counted = kTRUE;
	470	break;
	471	}
	472	}
	473	if (!counted) {
	474	trackLabel[nTracks] = label;
	475	trackCount[nTracks]++;
	476	nTracks++;
	477	if (nTracks == kMaxTracks) {
	478	Warning("CookLabel","Too many tracks for this tracklet.");
	479	nTracks--;
	480	break;
	481	}
	482	}
	483
	484	}
	485
	486	}
	487
	488	for (Int_t it = 0; it < kMaxTracks; it++) {
	489	if (trackCount[it] >= (Int_t)(frac*fNclusters)) {
	490	fTrackLabel = trackLabel[it];
	491	break;
	492	}
	493	}
	494
	495	}
	496
	497	//_____________________________________________________________________________
e3b2b5e5	498	Float_t AliTRDmcmTracklet::GetOmegaTau(Float_t vdrift, Float_t field) const
0ee00e25	499	{
	500	//
	501	// Returns omega*tau (tan(Lorentz-angle)) for a given drift velocity <vd>
	502	// and a B-field <b> for Xe/CO2 (15%).
	503	// The values are according to a GARFIELD simulation.
	504	//
	505
	506	//
	507	// Copy of the "AliTRDcalibDB" function, taking as argument the magnetic field too
	508	//
	509
	510	const Int_t kNb = 5;
	511	Float_t p0[kNb] = { 0.004810, 0.007412, 0.010252, 0.013409, 0.016888 };
	512	Float_t p1[kNb] = { 0.054875, 0.081534, 0.107333, 0.131983, 0.155455 };
	513	Float_t p2[kNb] = { -0.008682, -0.012896, -0.016987, -0.020880, -0.024623 };
	514	Float_t p3[kNb] = { 0.000155, 0.000238, 0.000330, 0.000428, 0.000541 };
	515
	516	Int_t ib = ((Int_t) (10 * (field - 0.15)));
	517	ib = TMath::Max( 0,ib);
	518	ib = TMath::Min(kNb,ib);
	519
	520	Float_t alphaL = p0[ib]
	521	+ p1[ib] * vdrift
	522	+ p2[ib] * vdrift*vdrift
	523	+ p3[ib] * vdriftvdriftvdrift;
	524
	525	return TMath::Tan(alphaL);
	526
	527	}