[u/mrichter/AliRoot.git] / PWGPP / TPC / AliTPCtaskQA.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.                  *
 **************************************************************************/
//
//
// 
//Task for analysis if TPC base QA 

// QA histogram
// Dimensions
// 0 - chi2
// 1 - number of clusters
// 2 - number of findable clusters
// 3 - number of clusters/ findable clusters  
// 4 - pt          - at the entrance of the TPC
// 5 - eta         - at the entrance of the TPC
// 6 - phi         - at the entrance of the TPC

//-----------------------------------------------------------------------
// Author : M.Ivanov  marian.ivanov@cern.ch - 
//-----------------------------------------------------------------------


//
// ROOT includes
#include <TChain.h>
#include <TMath.h>
#include <TVectorD.h>
#include <TSystem.h>
#include <TFile.h>
#include <TParticle.h>

// ALIROOT includes
#include <AliAnalysisManager.h>
#include <AliESDInputHandler.h>
#include "AliStack.h"
#include "AliMCEvent.h"
#include "AliMCEventHandler.h"
#include "AliMathBase.h"

#include <AliESD.h>
#include "AliExternalTrackParam.h"
#include "AliTracker.h"
#include "AliTPCseed.h"
//
#include "AliTPCtaskQA.h"
//
#include <THnSparse.h>

//

// STL includes
#include <iostream>

using namespace std;


ClassImp(AliTPCtaskQA)

//________________________________________________________________________
AliTPCtaskQA::AliTPCtaskQA() : 
  AliAnalysisTask(), 
  fMCinfo(0),     //! MC event handler
  fESD(0),
  fList(0),
  fTPCqa(0)
{
  //
  // Default constructor (should not be used)
  //
}

AliTPCtaskQA::AliTPCtaskQA(const AliTPCtaskQA& info) : 
  AliAnalysisTask(info), 
  fMCinfo(info.fMCinfo),     //! MC event handler
  fESD(info.fESD),        //!
  fList(0),
  fTPCqa(0)
{
  //
  // Dummy Copy  constructor - no copy constructor for THnSparse 
  //
  fList = (TObjArray*)(info.fList->Clone());
}


//________________________________________________________________________
AliTPCtaskQA::AliTPCtaskQA(const char *name) : 
  AliAnalysisTask(name, "AliTPCtaskQA"), 
  fMCinfo(0),     //! MC event handler
  fESD(0),
  fList(0),
  fTPCqa(0)
{
  //
  // Normal constructor
  //
  // Input slot #0 works with a TChain
  DefineInput(0, TChain::Class());
  // Output slot #0 writes into a TList
  DefineOutput(0, TObjArray::Class());
  //
  //make histos
  Init(); 
}

void AliTPCtaskQA::Init(){
  //
  // Init qa histogram
  // Dimensions
  //
  // 0 - chi2
  // 1 - number of clusters
  // 2 - number of findable clusters
  // 3 - number of clusters/ findable clusters  
  // 4 - pt          - at the entrance of the TPC
  // 5 - eta         - at the entrance of the TPC
  // 6 - phi         - at the entrance of the TPC
  

  Double_t xmin[7],  xmax[7];
  Int_t    nbins[7];
  // 
  nbins[0]=100;                // chi2
  xmin[0]=0; xmax[0]=10;
  //
  nbins[1]=80;                 // ncls
  xmin[1]=0; xmax[1]=160;
  //
  nbins[2]=80;                 // nclsf
  xmin[2]=0; xmax[2]=160;
  //
  nbins[3]=40;                 // ncls/nclsf
  xmin[3] =-0.1; xmax[3]=1.1;
  //
  nbins[4]=50;                 // pt
  xmin[4] =0.1; xmax[4]=100;

  nbins[5]=40;                 // eta
  xmin[5] =-2; xmax[5]=2;

  nbins[6]= 360;                 // phi - 10 bins per sector
  xmin[6] = -TMath::Pi(); xmax[6]=TMath::Pi();


  fTPCqa = new THnSparseS("TPC qa","TPC qa",7,nbins,xmin,xmax);
  //
  //
  BinLogX(fTPCqa->GetAxis(4));
  
  const char *hisAxisName[7] ={"chi2/N_{cl}","N_{cl}","N_{clF}","N_{clR}","p_{t}","#eta","#phi"};
  //  
  for (Int_t i=0;i<7;i++) {
    fTPCqa->GetAxis(i)->SetTitle(hisAxisName[i]);
    fTPCqa->GetAxis(i)->SetName(hisAxisName[i]);
  }
  fList = new TObjArray(3);
  fList->AddAt(fTPCqa,0);
}


AliTPCtaskQA::~AliTPCtaskQA(){
  //
  //
  //
  delete fTPCqa;
}


//________________________________________________________________________
void AliTPCtaskQA::ConnectInputData(Option_t *) 
{
  //
  // Connect the input data
  //

  TTree* tree=dynamic_cast<TTree*>(GetInputData(0));
  if (!tree) {
    //Printf("ERROR: Could not read chain from input slot 0");
  }
  else {
    AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
    if (!esdH) {
      //Printf("ERROR: Could not get ESDInputHandler");
    }
    else {
      fESD = esdH->GetEvent();
      //Printf("*** CONNECTED NEW EVENT ****");
    }  
  }
  AliMCEventHandler* mcinfo = (AliMCEventHandler*) (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());  
  mcinfo->SetReadTR(kTRUE);
  
  fMCinfo = mcinfo->MCEvent();


}


//________________________________________________________________________
void AliTPCtaskQA::CreateOutputObjects() 
{
  //
  // Connect the output objects
  //

}


//________________________________________________________________________
void AliTPCtaskQA::Exec(Option_t *) {
  //
  // Execute analysis for current event 
  //


  // If MC has been connected   

  if (!fMCinfo){
    cout << "Not MC info\n" << endl;
  }else{
    ProcessMCInfo();
  }
  //
  PostData(0, fList);
}      


void  AliTPCtaskQA::ProcessMCInfo(){
  //
  //
  //
  //

  if (!fTPCqa) Init();
  Int_t npart   = fMCinfo->GetNumberOfTracks();
  Int_t ntracks = fESD->GetNumberOfTracks(); 
  if (npart<=0) return;
  if (ntracks<=0) return;
  //
  //
  TParticle * particle= new TParticle;
  TClonesArray * trefs = new TClonesArray("AliTrackReference");
  
  for (Int_t itrack=0;itrack<ntracks;itrack++){
    AliESDtrack *track = fESD->GetTrack(itrack);
    //
    if ((track->GetStatus()&AliESDtrack::kTPCrefit)==0) continue;   // only refited tracks
    Int_t ipart = TMath::Abs(track->GetLabel());
    //
    Int_t status = fMCinfo->GetParticleAndTR(ipart, particle, trefs);
    if (status<0) continue;
    if (!particle) continue;
    if (!trefs) continue;
    //
    //
    AliTrackReference *tpcRef=0;
    for (Int_t iref=0; iref<trefs->GetEntries(); iref++){
      AliTrackReference *ref = (AliTrackReference*)trefs->At(iref);
      if (ref->DetectorId()== AliTrackReference::kTPC){
	tpcRef=ref;
	break;
      }
    }
    if (!tpcRef) continue;
    
    //
    // Fill histos
    //
    Double_t x[7];
    x[0]= track->GetTPCchi2()/track->GetTPCNcls();
    x[1]= track->GetTPCNcls();
    x[2]= track->GetTPCNclsF();
    x[3]= Float_t(track->GetTPCNcls())/Float_t(track->GetTPCNclsF());
    x[4]= tpcRef->Pt();
    x[5]= -0.5*TMath::Log((tpcRef->P()+tpcRef->Pz())/(tpcRef->P()-tpcRef->Pz()));
    x[6]= TMath::ATan2(tpcRef->Y(),tpcRef->X());
    //
    fTPCqa->Fill(x);
  }
}


void AliTPCtaskQA::BinLogX(TAxis *axis) {
  //
  //
  //
  Int_t bins = axis->GetNbins();
  
  Double_t from = axis->GetXmin();
  Double_t to   = axis->GetXmax();
  Double_t *new_bins = new Double_t[bins + 1];
  
  new_bins[0] = from;
  Double_t factor = pow(to/from, 1./bins);
  
  for (Int_t i = 1; i <= bins; i++) {
    new_bins[i] = factor * new_bins[i-1];
  }
  axis->Set(bins, new_bins);
  delete [] new_bins;
}


AliTPCtaskQA* AliTPCtaskQA::ReadFromFile(const char *fname){
  //
  //
  //
  AliTPCtaskQA *qa = new AliTPCtaskQA;
  TFile *f = new TFile(fname);
  TObjArray *array= (TObjArray*)f->Get("tpcTaskQA");
  qa->fTPCqa =  (THnSparse*)array->At(0);  
  delete f;
  return qa;
}
Commit	Line	Data
7cc34f08	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	//
	17	//
	18	//Task for analysis if TPC base QA
	19
	20	// QA histogram
	21	// Dimensions
	22	// 0 - chi2
	23	// 1 - number of clusters
	24	// 2 - number of findable clusters
	25	// 3 - number of clusters/ findable clusters
	26	// 4 - pt - at the entrance of the TPC
	27	// 5 - eta - at the entrance of the TPC
	28	// 6 - phi - at the entrance of the TPC
	29
	30	//-----------------------------------------------------------------------
	31	// Author : M.Ivanov marian.ivanov@cern.ch -
	32	//-----------------------------------------------------------------------
	33
	34
	35
	36	//
	37	// ROOT includes
	38	#include <TChain.h>
	39	#include <TMath.h>
	40	#include <TVectorD.h>
	41	#include <TSystem.h>
	42	#include <TFile.h>
	43	#include <TParticle.h>
	44
	45	// ALIROOT includes
	46	#include <AliAnalysisManager.h>
	47	#include <AliESDInputHandler.h>
	48	#include "AliStack.h"
	49	#include "AliMCEvent.h"
	50	#include "AliMCEventHandler.h"
	51	#include "AliMathBase.h"
	52
	53	#include <AliESD.h>
	54	#include "AliExternalTrackParam.h"
	55	#include "AliTracker.h"
	56	#include "AliTPCseed.h"
	57	//
	58	#include "AliTPCtaskQA.h"
	59	//
	60	#include <THnSparse.h>
	61
	62	//
	63
	64	// STL includes
65	#include <iostream>
66
67	using namespace std;
68
69
70	ClassImp(AliTPCtaskQA)
71
72	//________________________________________________________________________
73	AliTPCtaskQA::AliTPCtaskQA() :
74	AliAnalysisTask(),
75	fMCinfo(0), //! MC event handler
76	fESD(0),
77	fList(0),
78	fTPCqa(0)
79	{
80	//
81	// Default constructor (should not be used)
82	//
83	}
84
85	AliTPCtaskQA::AliTPCtaskQA(const AliTPCtaskQA& info) :
86	AliAnalysisTask(info),
87	fMCinfo(info.fMCinfo), //! MC event handler
88	fESD(info.fESD), //!
89	fList(0),
90	fTPCqa(0)
91	{
92	//
93	// Dummy Copy constructor - no copy constructor for THnSparse
94	//
95	fList = (TObjArray*)(info.fList->Clone());
96	}
97
98
99
100	//________________________________________________________________________
101	AliTPCtaskQA::AliTPCtaskQA(const char *name) :
102	AliAnalysisTask(name, "AliTPCtaskQA"),
103	fMCinfo(0), //! MC event handler
104	fESD(0),
105	fList(0),
106	fTPCqa(0)
107	{
108	//
109	// Normal constructor
110	//
111	// Input slot #0 works with a TChain
112	DefineInput(0, TChain::Class());
113	// Output slot #0 writes into a TList
114	DefineOutput(0, TObjArray::Class());
115	//
116	//make histos
117	Init();
118	}
119
120	void AliTPCtaskQA::Init(){
121	//
122	// Init qa histogram
123	// Dimensions
124	//
125	// 0 - chi2
126	// 1 - number of clusters
127	// 2 - number of findable clusters
128	// 3 - number of clusters/ findable clusters
129	// 4 - pt - at the entrance of the TPC
130	// 5 - eta - at the entrance of the TPC
131	// 6 - phi - at the entrance of the TPC
132
133
134
135	Double_t xmin[7], xmax[7];
136	Int_t nbins[7];
137	//
138	nbins[0]=100; // chi2
139	xmin[0]=0; xmax[0]=10;
140	//
141	nbins[1]=80; // ncls
142	xmin[1]=0; xmax[1]=160;
143	//
144	nbins[2]=80; // nclsf
145	xmin[2]=0; xmax[2]=160;
146	//
147	nbins[3]=40; // ncls/nclsf
148	xmin[3] =-0.1; xmax[3]=1.1;
149	//
150	nbins[4]=50; // pt
151	xmin[4] =0.1; xmax[4]=100;
152
153	nbins[5]=40; // eta
154	xmin[5] =-2; xmax[5]=2;
155
156	nbins[6]= 360; // phi - 10 bins per sector
157	xmin[6] = -TMath::Pi(); xmax[6]=TMath::Pi();
158
159
160
161	fTPCqa = new THnSparseS("TPC qa","TPC qa",7,nbins,xmin,xmax);
162	//
163	//
164	BinLogX(fTPCqa->GetAxis(4));
165
c6d45264	166	const char *hisAxisName[7] ={"chi2/N_{cl}","N_{cl}","N_{clF}","N_{clR}","p_{t}","#eta","#phi"};
7cc34f08	167	//
	168	for (Int_t i=0;i<7;i++) {
	169	fTPCqa->GetAxis(i)->SetTitle(hisAxisName[i]);
	170	fTPCqa->GetAxis(i)->SetName(hisAxisName[i]);
	171	}
	172	fList = new TObjArray(3);
	173	fList->AddAt(fTPCqa,0);
	174	}
	175
	176
	177
	178
	179	AliTPCtaskQA::~AliTPCtaskQA(){
	180	//
	181	//
	182	//
	183	delete fTPCqa;
	184	}
	185
	186
	187	//________________________________________________________________________
	188	void AliTPCtaskQA::ConnectInputData(Option_t *)
	189	{
	190	//
	191	// Connect the input data
	192	//
	193
	194	TTree* tree=dynamic_cast<TTree*>(GetInputData(0));
	195	if (!tree) {
	196	//Printf("ERROR: Could not read chain from input slot 0");
	197	}
	198	else {
	199	AliESDInputHandler esdH = dynamic_cast<AliESDInputHandler> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
	200	if (!esdH) {
	201	//Printf("ERROR: Could not get ESDInputHandler");
	202	}
	203	else {
	204	fESD = esdH->GetEvent();
	205	//Printf("* CONNECTED NEW EVENT **");
	206	}
	207	}
	208	AliMCEventHandler* mcinfo = (AliMCEventHandler*) (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
	209	mcinfo->SetReadTR(kTRUE);
	210
	211	fMCinfo = mcinfo->MCEvent();
	212
	213
	214	}
	215
	216
	217
	218
	219
	220
	221	//________________________________________________________________________
	222	void AliTPCtaskQA::CreateOutputObjects()
	223	{
	224	//
	225	// Connect the output objects
	226	//
	227
	228	}
	229
	230
231	//________________________________________________________________________
232	void AliTPCtaskQA::Exec(Option_t *) {
233	//
234	// Execute analysis for current event
235	//
236
237
238	// If MC has been connected
239
240	if (!fMCinfo){
241	cout << "Not MC info\n" << endl;
242	}else{
243	ProcessMCInfo();
244	}
245	//
246	PostData(0, fList);
247	}
248
249
250
251
252
253
254
255
256
257
258	void AliTPCtaskQA::ProcessMCInfo(){
259	//
260	//
261	//
262	//
263
264	if (!fTPCqa) Init();
265	Int_t npart = fMCinfo->GetNumberOfTracks();
266	Int_t ntracks = fESD->GetNumberOfTracks();
267	if (npart<=0) return;
268	if (ntracks<=0) return;
269	//
270	//
271	TParticle * particle= new TParticle;
272	TClonesArray * trefs = new TClonesArray("AliTrackReference");
273
274	for (Int_t itrack=0;itrack<ntracks;itrack++){
275	AliESDtrack *track = fESD->GetTrack(itrack);
276	//
277	if ((track->GetStatus()&AliESDtrack::kTPCrefit)==0) continue; // only refited tracks
278	Int_t ipart = TMath::Abs(track->GetLabel());
279	//
280	Int_t status = fMCinfo->GetParticleAndTR(ipart, particle, trefs);
281	if (status<0) continue;
282	if (!particle) continue;
283	if (!trefs) continue;
284	//
285	//
286	AliTrackReference *tpcRef=0;
287	for (Int_t iref=0; iref<trefs->GetEntries(); iref++){
288	AliTrackReference ref = (AliTrackReference)trefs->At(iref);
289	if (ref->DetectorId()== AliTrackReference::kTPC){
290	tpcRef=ref;
291	break;
292	}
293	}
294	if (!tpcRef) continue;
295
296	//
297	// Fill histos
298	//
299	Double_t x[7];
300	x[0]= track->GetTPCchi2()/track->GetTPCNcls();
301	x[1]= track->GetTPCNcls();
302	x[2]= track->GetTPCNclsF();
303	x[3]= Float_t(track->GetTPCNcls())/Float_t(track->GetTPCNclsF());
304	x[4]= tpcRef->Pt();
305	x[5]= -0.5*TMath::Log((tpcRef->P()+tpcRef->Pz())/(tpcRef->P()-tpcRef->Pz()));
306	x[6]= TMath::ATan2(tpcRef->Y(),tpcRef->X());
307	//
308	fTPCqa->Fill(x);
309	}
310	}
311
312
313
314
315
316
317	void AliTPCtaskQA::BinLogX(TAxis *axis) {
318	//
319	//
320	//
321	Int_t bins = axis->GetNbins();
322
323	Double_t from = axis->GetXmin();
324	Double_t to = axis->GetXmax();
325	Double_t *new_bins = new Double_t[bins + 1];
326
327	new_bins[0] = from;
328	Double_t factor = pow(to/from, 1./bins);
329
330	for (Int_t i = 1; i <= bins; i++) {
331	new_bins[i] = factor * new_bins[i-1];
332	}
333	axis->Set(bins, new_bins);
4ce766eb	334	delete [] new_bins;
7cc34f08	335	}
	336
	337
	338	AliTPCtaskQA* AliTPCtaskQA::ReadFromFile(const char *fname){
	339	//
	340	//
	341	//
	342	AliTPCtaskQA *qa = new AliTPCtaskQA;
	343	TFile *f = new TFile(fname);
	344	TObjArray array= (TObjArray)f->Get("tpcTaskQA");
	345	qa->fTPCqa = (THnSparse*)array->At(0);
	346	delete f;
	347	return qa;
	348	}
	349