[u/mrichter/AliRoot.git] / ITS / ITSOccupancy.C

/******************************************************************************

  "ITSOccupancy.C"
  
  this macro calculates the mean occupancy of each ITS layer, counting the
  "fired" digits of each module, and making two overall calculi:
     1) the calculus of the mean overall layer occupancy (as the ratio
	     between the total number of active digits and the total number of 
		  channels in the module;
	  2) a distribution of the occupancies as a funcion of z, to obtain some
		  kind of most significand data for this value, along the z axis
	
  The macro creates also the GIF of the TCanvas where the histograms are plotted
			  
  The arguments that must be inserted are the following:
     1) the answer to the following question (as a Bool_t)
		  "Can I get the total channels distribution from the file totals.root?"
		  if the answer is NO (false) these histograms are re-created and re-stored 
		  in that file.
			 
	  2) the name of the file to analyze; 
	     NOTE: in this argument you MUST omit the ".root" extension, because the 
		        macro uses this argument to recall the file, but also to name the
				  GIF & PS of the output canvas, in order to save the plot of the 
				  distributions to the disk. If you put, as argument, the string
				  "file.root", the macro will look for an unexisting file named
				  "file.root.root"...
                                                                         
  It is also possible (and recommended) to compile this macro, 
  to execute it faster. To do this, it is necessary to:
  1) make sure that de "#define __COMPILED__" instruction below is UNcommented
     and the "#define __NOCOMPILED__" instruction is commented (remember that
	  if you comment or uncomment both these instructions, the macro can't work)
  2) type, at the root prompt, the instruction 
     gSystem->SetIncludePath("-I- -I$ALICE_ROOT/include -I$ALICE_ROOT/ITS -g")
	  to adjoin the ALICE include libraries to the main include directory
  3) load the function via the ".L ITSOccupancy.C++" statement (only the first
     time, because the first time the macro must be compiled). Frome the 
	  second time you use the macro, you must only load it via the 
	  ".L ITSOccupancy.C+" instruction (note the number of '+' signs in each case
  4) call the function with only its name, es: ITSOccupancy("HP", "galice.root")

 Author: Alberto Pulvirenti                                             
******************************************************************************/

#define __COMPILED__
//#define __NOCOMPILED__

#ifdef __COMPILED__
	#include <iostream.h>
	#include <TROOT.h>
	#include <TArrayI.h>
	#include <TCanvas.h>
	#include <TClassTable.h>
	#include <TClonesArray.h>
	#include <TFile.h>
	#include <TObject.h>
	#include <TObjArray.h>
	#include <TTree.h>
	#include <TMath.h>
	#include <TString.h>
	#include <TH1.h>
	#include <AliRun.h>
	#include <AliITS.h>
	#include <AliITSgeom.h>
	#include <AliITSDetType.h>
	#include <AliITSRecPoint.h>
	#include <AliITSdigit.h>
	#include <AliITShit.h>
	#include <AliITSmodule.h> 
	#include <AliITSsegmentation.h>
	#include <AliITSsegmentationSPD.h> 
	#include <AliITSsegmentationSDD.h>
	#include <AliITSsegmentationSSD.h>
#endif

Int_t ITSOccupancy(char *name = "galice", Int_t evNum = 0) {       
		 
	// Evaluate the name of the file to open...
	Int_t filelen = strlen(name);
	Text_t *fileroot = new Text_t[filelen + 6];
	// ...and copy this name to the data ROOT file
	strcpy (fileroot, name);
	strcat (fileroot, ".root");
	
	// Open the Root input file containing the AliRun data and 
	// the Root output data file
	TFile *froot = (TFile*)gROOT->GetListOfFiles()->FindObject(fileroot);
	if (!froot) froot = new TFile(fileroot, "READ");
	// Load the gAlice shared libs if not already in memory and
	// clears an eventually existing AliRun object
#ifdef __NOCOMPILED__
	if (gClassTable->GetID("AliRun") < 0) {
    	gROOT->LoadMacro("loadlibs.C");
    	loadlibs();
  	}
#endif
   if (gAlice) {
      delete gAlice;
	   gAlice = 0;
   }
	// Get the AliRun object from file (if there's none, the macro is stopped)
	Int_t nparticles = 0;
	gAlice = (AliRun*)froot->Get("gAlice");
	if (gAlice) {
		cout << "Found an AliRun object in `" << fileroot << "'" << endl;
		// In this case, select the event number specified (default is 0)
		nparticles = gAlice->GetEvent(evNum);
		if (nparticles)
			cout << "\nNumber of particles   = " << nparticles << endl;
		else {
			cout << "NO PARTICLES!!!!" << endl;
			return 0;
		}
	}
	else {
		cout<<"Not found any AliRun object in `" << fileroot << "'" << endl;
		return 0;
	}
	// Initialize the ITS and its geometry
	AliITS *ITS = (AliITS*)gAlice->GetModule("ITS");
	AliITSgeom *gm = ITS->GetITSgeom();

	// Variables and objects definition and setting
	Float_t zmax[] = {20.0,20.0,25.0,35.0,49.5,54.0}; // z edges for TH1F (cm)
	Int_t nbins[]  = {20,20,20,28,22,24};             // bins number for TH1Fs
	
	// Histos for plotting occupancy distributions
	TH1F *above[6], *below[6];
	
	for (Int_t lay = 0; lay < 6; lay++) {
		Int_t nlads = gm->GetNladders(lay+1);
		Int_t ndets = gm->GetNdetectors(lay+1);
		Int_t dtype = lay / 2;
		Int_t minID = gm->GetModuleIndex(lay+1, 1, 1);
		Int_t maxID = gm->GetModuleIndex(lay+1, nlads, ndets);
		Text_t ffname[20];
		sprintf(ffname, "h_%d", lay+1);
		below[lay] = new TH1F(ffname, "Total z distribution of digits", nbins[lay], -zmax[lay], zmax[lay]);
		cout << "Evaluating total channels number of layer " << lay+1 << endl;
		for (Int_t I = minID; I <= maxID; I++) {		
			AliITSDetType *det = ITS->DetType(dtype);
			AliITSsegmentation *seg = det->GetSegmentationModel();
			Int_t NX = seg->Npx();
			Int_t NZ = seg->Npz();
			for (Int_t ix = 0; ix <= NX; ix++) {
				for (Int_t iz = 0; iz <= NZ; iz++) {
					Float_t lx[] = {0.0,0.0,0.0}, gx[] = {0.0,0.0,0.0};
					seg->DetToLocal(ix, iz, lx[0], lx[2]);
					gm->LtoG(I, lx, gx);
					below[lay]->Fill(gx[2]);
				}
			}
		}
		sprintf(ffname, "H_%d", lay+1);
		above[lay] = new TH1F(ffname, "histo", nbins[lay], -zmax[lay], zmax[lay]);
	}
		
	// Counting the hits, digits and recpoints contained in the ITS
	TTree *TD = gAlice->TreeD();
	ITS->ResetDigits();
	Float_t mean[6];

	for (Int_t L = 0; L < 6; L++) {
		
		cout << "Layer " << L + 1 << ": " << flush;				

		// To avoid two nested loops, are calculated 
		// the ID of the first and last module of the L
		// (remember the L goest from 0 to 5 (not from 1 to 6)
		Int_t first, last;
		first = gm->GetModuleIndex(L + 1, 1, 1);
		last = gm->GetModuleIndex(L + 1, gm->GetNladders(L + 1), gm->GetNdetectors(L + 1));
				
		// Start loop on modules
		for (Int_t ID = first; ID <= last; ID++) {
			Int_t dtype = L / 2;
			AliITSDetType *det = ITS->DetType(dtype);
			AliITSsegmentation *seg = det->GetSegmentationModel();
			if (dtype == 2) seg->SetLayer(L+1);
			
			TD->GetEvent(ID);
			TClonesArray *digits_array = ITS->DigitsAddress(dtype);
			Int_t digits_num = digits_array->GetEntries();
			// Get the coordinates of the module
		  	for (Int_t j = 0; j < digits_num; j++) {
				Float_t lx[] = {0.0,0.0,0.0}, gx[] = {0.0,0.0,0.0};
				AliITSdigit *digit = (AliITSdigit*)digits_array->UncheckedAt(j);
				Int_t iz=digit->fCoord1;  // cell number z
				Int_t ix=digit->fCoord2;  // cell number x
    			// Get local coordinates of the element (microns)
 		   	seg->DetToLocal(ix, iz, lx[0], lx[2]);
				gm->LtoG(ID, lx, gx);
				above[L]->Fill(gx[2]);
			}
		}
		
		above[L]->Divide(above[L], below[L], 100.0, 1.0);
		mean[L] = above[L]->GetSumOfWeights() / above[L]->GetNbinsX(); 
		cout.setf(ios::fixed);
		cout.precision(2);
		cout << " digits occupancy = " << mean[L] << "%" << endl;
	}
	
	TCanvas *view = new TCanvas("view", "Digits occupancy distributions", 600, 900);
	view->Divide(2, 3);
	
	for (Int_t I = 1; I < 7; I++) {
		view->cd(I);
		Text_t title[50];
		sprintf(title, "Layer %d: %4.2f%c", I, mean[I-1], '%');
		title[6] = (Text_t)I + '0';
		above[I-1]->SetTitle(title);
		above[I-1]->SetStats(kFALSE);
		above[I-1]->SetXTitle("z (cm)");
		above[I-1]->SetYTitle("%");
		above[I-1]->Draw();
		view->Update();
	}
	Text_t *filegif  = new Text_t[filelen + 23];
	Text_t *fileps  = new Text_t[filelen + 23];
	strcpy (filegif, name);
	strcpy (fileps, name);
	strcat (filegif, "_digit_occupancy.gif");
	strcat (fileps, "_digit_occupancy.eps");
	view->SaveAs(filegif);
	view->SaveAs(fileps);
	

	return 1;
}
Commit	Line	Data
2bba24f4	1	/******************************************************************************
	2
	3	"ITSOccupancy.C"
	4
faafc0eb	5	this macro calculates the mean occupancy of each ITS layer, counting the
	6	"fired" digits of each module, and making two overall calculi:
	7	1) the calculus of the mean overall layer occupancy (as the ratio
	8	between the total number of active digits and the total number of
	9	channels in the module;
	10	2) a distribution of the occupancies as a funcion of z, to obtain some
	11	kind of most significand data for this value, along the z axis
	12
	13	The macro creates also the GIF of the TCanvas where the histograms are plotted
	14
	15	The arguments that must be inserted are the following:
	16	1) the answer to the following question (as a Bool_t)
	17	"Can I get the total channels distribution from the file totals.root?"
	18	if the answer is NO (false) these histograms are re-created and re-stored
	19	in that file.
	20
	21	2) the name of the file to analyze;
	22	NOTE: in this argument you MUST omit the ".root" extension, because the
	23	macro uses this argument to recall the file, but also to name the
	24	GIF & PS of the output canvas, in order to save the plot of the
	25	distributions to the disk. If you put, as argument, the string
	26	"file.root", the macro will look for an unexisting file named
	27	"file.root.root"...
2bba24f4	28
faafc0eb	29	It is also possible (and recommended) to compile this macro,
faafc0eb	30	to execute it faster. To do this, it is necessary to:
2bba24f4	31	1) make sure that de "#define __COMPILED__" instruction below is UNcommented
	32	and the "#define __NOCOMPILED__" instruction is commented (remember that
	33	if you comment or uncomment both these instructions, the macro can't work)
	34	2) type, at the root prompt, the instruction
	35	gSystem->SetIncludePath("-I- -I$ALICE_ROOT/include -I$ALICE_ROOT/ITS -g")
	36	to adjoin the ALICE include libraries to the main include directory
	37	3) load the function via the ".L ITSOccupancy.C++" statement (only the first
	38	time, because the first time the macro must be compiled). Frome the
	39	second time you use the macro, you must only load it via the
	40	".L ITSOccupancy.C+" instruction (note the number of '+' signs in each case
	41	4) call the function with only its name, es: ITSOccupancy("HP", "galice.root")
	42
2bba24f4	43	Author: Alberto Pulvirenti
	44	******************************************************************************/
	45
	46	#define __COMPILED__
	47	//#define __NOCOMPILED__
	48
	49	#ifdef __COMPILED__
	50	#include <iostream.h>
	51	#include <TROOT.h>
	52	#include <TArrayI.h>
	53	#include <TCanvas.h>
	54	#include <TClassTable.h>
	55	#include <TClonesArray.h>
	56	#include <TFile.h>
	57	#include <TObject.h>
	58	#include <TObjArray.h>
	59	#include <TTree.h>
	60	#include <TMath.h>
	61	#include <TString.h>
	62	#include <TH1.h>
	63	#include <AliRun.h>
	64	#include <AliITS.h>
	65	#include <AliITSgeom.h>
	66	#include <AliITSDetType.h>
	67	#include <AliITSRecPoint.h>
	68	#include <AliITSdigit.h>
	69	#include <AliITShit.h>
	70	#include <AliITSmodule.h>
	71	#include <AliITSsegmentation.h>
	72	#include <AliITSsegmentationSPD.h>
	73	#include <AliITSsegmentationSDD.h>
	74	#include <AliITSsegmentationSSD.h>
	75	#endif
	76
6782062e	77	Int_t ITSOccupancy(char *name = "galice", Int_t evNum = 0) {
2bba24f4	78
faafc0eb	79	// Evaluate the name of the file to open...
	80	Int_t filelen = strlen(name);
	81	Text_t *fileroot = new Text_t[filelen + 6];
	82	// ...and copy this name to the data ROOT file
	83	strcpy (fileroot, name);
	84	strcat (fileroot, ".root");
	85
	86	// Open the Root input file containing the AliRun data and
	87	// the Root output data file
	88	TFile froot = (TFile)gROOT->GetListOfFiles()->FindObject(fileroot);
	89	if (!froot) froot = new TFile(fileroot, "READ");
	90	// Load the gAlice shared libs if not already in memory and
	91	// clears an eventually existing AliRun object
2bba24f4	92	#ifdef __NOCOMPILED__
	93	if (gClassTable->GetID("AliRun") < 0) {
	94	gROOT->LoadMacro("loadlibs.C");
	95	loadlibs();
	96	}
	97	#endif
faafc0eb	98	if (gAlice) {
2bba24f4	99	delete gAlice;
faafc0eb	100	gAlice = 0;
2bba24f4	101	}
2bba24f4	102	// Get the AliRun object from file (if there's none, the macro is stopped)
faafc0eb	103	Int_t nparticles = 0;
	104	gAlice = (AliRun*)froot->Get("gAlice");
	105	if (gAlice) {
	106	cout << "Found an AliRun object in `" << fileroot << "'" << endl;
	107	// In this case, select the event number specified (default is 0)
	108	nparticles = gAlice->GetEvent(evNum);
	109	if (nparticles)
	110	cout << "\nNumber of particles = " << nparticles << endl;
	111	else {
	112	cout << "NO PARTICLES!!!!" << endl;
	113	return 0;
	114	}
2bba24f4	115	}
faafc0eb	116	else {
faafc0eb	117	cout<<"Not found any AliRun object in `" << fileroot << "'" << endl;
2bba24f4	118	return 0;
	119	}
	120	// Initialize the ITS and its geometry
	121	AliITS ITS = (AliITS)gAlice->GetModule("ITS");
	122	AliITSgeom *gm = ITS->GetITSgeom();
2bba24f4	123
faafc0eb	124	// Variables and objects definition and setting
	125	Float_t zmax[] = {20.0,20.0,25.0,35.0,49.5,54.0}; // z edges for TH1F (cm)
	126	Int_t nbins[] = {20,20,20,28,22,24}; // bins number for TH1Fs
	127
	128	// Histos for plotting occupancy distributions
	129	TH1F above[6], below[6];
	130
6782062e	131	for (Int_t lay = 0; lay < 6; lay++) {
	132	Int_t nlads = gm->GetNladders(lay+1);
	133	Int_t ndets = gm->GetNdetectors(lay+1);
	134	Int_t dtype = lay / 2;
	135	Int_t minID = gm->GetModuleIndex(lay+1, 1, 1);
	136	Int_t maxID = gm->GetModuleIndex(lay+1, nlads, ndets);
	137	Text_t ffname[20];
	138	sprintf(ffname, "h_%d", lay+1);
	139	below[lay] = new TH1F(ffname, "Total z distribution of digits", nbins[lay], -zmax[lay], zmax[lay]);
	140	cout << "Evaluating total channels number of layer " << lay+1 << endl;
	141	for (Int_t I = minID; I <= maxID; I++) {
	142	AliITSDetType *det = ITS->DetType(dtype);
	143	AliITSsegmentation *seg = det->GetSegmentationModel();
	144	Int_t NX = seg->Npx();
	145	Int_t NZ = seg->Npz();
	146	for (Int_t ix = 0; ix <= NX; ix++) {
	147	for (Int_t iz = 0; iz <= NZ; iz++) {
	148	Float_t lx[] = {0.0,0.0,0.0}, gx[] = {0.0,0.0,0.0};
	149	seg->DetToLocal(ix, iz, lx[0], lx[2]);
	150	gm->LtoG(I, lx, gx);
	151	below[lay]->Fill(gx[2]);
faafc0eb	152	}
faafc0eb	153	}
2bba24f4	154	}
6782062e	155	sprintf(ffname, "H_%d", lay+1);
6782062e	156	above[lay] = new TH1F(ffname, "histo", nbins[lay], -zmax[lay], zmax[lay]);
2bba24f4	157	}
6782062e	158
faafc0eb	159	// Counting the hits, digits and recpoints contained in the ITS
	160	TTree *TD = gAlice->TreeD();
	161	ITS->ResetDigits();
	162	Float_t mean[6];
2bba24f4	163
	164	for (Int_t L = 0; L < 6; L++) {
	165
faafc0eb	166	cout << "Layer " << L + 1 << ": " << flush;
faafc0eb	167
2bba24f4	168	// To avoid two nested loops, are calculated
	169	// the ID of the first and last module of the L
	170	// (remember the L goest from 0 to 5 (not from 1 to 6)
	171	Int_t first, last;
	172	first = gm->GetModuleIndex(L + 1, 1, 1);
	173	last = gm->GetModuleIndex(L + 1, gm->GetNladders(L + 1), gm->GetNdetectors(L + 1));
	174
	175	// Start loop on modules
2bba24f4	176	for (Int_t ID = first; ID <= last; ID++) {
faafc0eb	177	Int_t dtype = L / 2;
	178	AliITSDetType *det = ITS->DetType(dtype);
	179	AliITSsegmentation *seg = det->GetSegmentationModel();
	180	if (dtype == 2) seg->SetLayer(L+1);
2bba24f4	181
faafc0eb	182	TD->GetEvent(ID);
	183	TClonesArray *digits_array = ITS->DigitsAddress(dtype);
	184	Int_t digits_num = digits_array->GetEntries();
	185	// Get the coordinates of the module
	186	for (Int_t j = 0; j < digits_num; j++) {
	187	Float_t lx[] = {0.0,0.0,0.0}, gx[] = {0.0,0.0,0.0};
	188	AliITSdigit digit = (AliITSdigit)digits_array->UncheckedAt(j);
	189	Int_t iz=digit->fCoord1; // cell number z
	190	Int_t ix=digit->fCoord2; // cell number x
	191	// Get local coordinates of the element (microns)
	192	seg->DetToLocal(ix, iz, lx[0], lx[2]);
	193	gm->LtoG(ID, lx, gx);
	194	above[L]->Fill(gx[2]);
2bba24f4	195	}
2bba24f4	196	}
faafc0eb	197
	198	above[L]->Divide(above[L], below[L], 100.0, 1.0);
	199	mean[L] = above[L]->GetSumOfWeights() / above[L]->GetNbinsX();
	200	cout.setf(ios::fixed);
	201	cout.precision(2);
	202	cout << " digits occupancy = " << mean[L] << "%" << endl;
	203	}
	204
	205	TCanvas *view = new TCanvas("view", "Digits occupancy distributions", 600, 900);
	206	view->Divide(2, 3);
	207
	208	for (Int_t I = 1; I < 7; I++) {
	209	view->cd(I);
	210	Text_t title[50];
	211	sprintf(title, "Layer %d: %4.2f%c", I, mean[I-1], '%');
	212	title[6] = (Text_t)I + '0';
	213	above[I-1]->SetTitle(title);
	214	above[I-1]->SetStats(kFALSE);
	215	above[I-1]->SetXTitle("z (cm)");
	216	above[I-1]->SetYTitle("%");
	217	above[I-1]->Draw();
	218	view->Update();
	219	}
6782062e	220	Text_t *filegif = new Text_t[filelen + 23];
6782062e	221	Text_t *fileps = new Text_t[filelen + 23];
faafc0eb	222	strcpy (filegif, name);
faafc0eb	223	strcpy (fileps, name);
6782062e	224	strcat (filegif, "_digit_occupancy.gif");
6782062e	225	strcat (fileps, "_digit_occupancy.eps");
faafc0eb	226	view->SaveAs(filegif);
faafc0eb	227	view->SaveAs(fileps);
2bba24f4	228
2bba24f4	229
faafc0eb	230	return 1;
2bba24f4	231	}