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

#include "iostream.h"


  void AliITSComparisonV1(Int_t evNumber1=0,Int_t evNumber2=0) {

  const char *filename="itstracks.root";
  
  ///////////////// Dynamically link some shared libs ////////////////////////////////
  
  if (gClassTable->GetID("AliRun") < 0) {
    gROOT->LoadMacro("loadlibs.C");
    loadlibs();
  } else {
    delete gAlice;
    gAlice=0;
  }

// Connect the Root Galice file containing Geometry, Kine and Hits
   TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
   if (!file) file = new TFile(filename);
   
   ifstream in("itsgood_tracks"); 
//
//   Loop over events 
//

   char tname[30];

     
  struct GoodTrack {
   Int_t fEventN;  
    Int_t lab,code;
    Float_t px,py,pz,x,y,z,pxg,pyg,pzg,ptg;
    Bool_t flag;
  };   
    Int_t i;
  GoodTrack gt[15000]; 
  Int_t ngood=0;  
  cerr<<"Reading itsgood tracks...\n";
  while (in>>gt[ngood].fEventN>>gt[ngood].lab>>gt[ngood].code
	  >>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz
	  >>gt[ngood].x  >>gt[ngood].y >>gt[ngood].z
	  >>gt[ngood].pxg  >>gt[ngood].pyg >>gt[ngood].pzg
	  >>gt[ngood].ptg >>gt[ngood].flag) {
    ngood++;
    cerr<<ngood<<'\r';
    if (ngood==15000) {
      cerr<<"Too many good tracks !\n";
      break;
    }
  }
  if (!in.eof()) cerr<<"Read error (itsgood_tracks) !\n";
   ofstream out1 ("AliITSTrag.out");
  Int_t countpos=0,countneg=0;

   
  //for (i=0; i<ngood; i++) out1 << gt[i].ptg << "\n";
  for (i=0; i<ngood; i++) {
    out1 << gt[i].ptg << "\n";
    Int_t codpar=gt[i].code;
    if(codpar==2212 || codpar==-11 || codpar==-13 || codpar==211 || codpar==321 || codpar==3222
	 || codpar==213 || codpar==323 || codpar==10323 || codpar==3224 || codpar==2224 || codpar==2214
	 || codpar==-1114 || codpar==-3112 || codpar==-3312 || codpar==3224 || codpar==-3114 || codpar==-3314
	 || codpar==411 || codpar==431 || codpar==413 || codpar==433 || codpar==-15 || codpar==4232
	 || codpar==4222 || codpar==4322 || codpar==4422 || codpar==4412 || codpar==4432 || codpar==4224 
	 ||codpar==4214 || codpar==4324 || codpar==4424 || codpar==4414 || codpar==4434 || codpar==4444)
    countpos++;
	 if(codpar==-2212 || codpar==11 || codpar==13 || codpar==-211 || codpar==-321 || codpar==3112
	 || codpar==-213 || codpar==-323 || codpar==-10323 || codpar==3114 || codpar==1114 || codpar==-2224
	 || codpar==-2214 || codpar==33112 || codpar==-3222 || codpar==3114 || codpar==3314 || codpar==3334 
	 || codpar==-3224 || codpar==-411 || codpar==-431 || codpar==-413 || codpar==-433 || codpar==15 
	 || codpar==-4422 || codpar==-4432 || codpar==-4214 || codpar==-4324 || codpar==-4424 || codpar==-4434 
	 || codpar==-444)
	 countneg++; 		   
  }
  out1.close();
   ofstream out ("AliITSTra.out"); 
///   definition of nm of good track for each event
   TVector ntrackevtpc(evNumber2-evNumber1 +1);
   Int_t nchange=-1;
   Int_t nmev=-100;
   
   
    Int_t i;
    for( i =0; i<ngood ; i++){ 
    if(gt[i].fEventN != nmev ){
    nmev=gt[i].fEventN;
    nchange++; 
    }
    if(nmev == gt[i].fEventN) ntrackevtpc(nchange)++;
    }  
//////////////////////////////////////////////////////////////     
    
   for (Int_t nev=evNumber1; nev<= evNumber2; nev++) {


   sprintf(tname,"TreeT%d",nev);
   TTree *tracktree=(TTree*)file->Get(tname);
   TBranch *tbranch=tracktree->GetBranch("ITStracks");
   cout<<" nev = "<<nev<<"\n";   
	//cout<<" open the file \n"; 
	
   Int_t nentr=tracktree->GetEntries();

   TObjArray tarray(nentr);
  // AliITSIOTrack *iotrack=0;
   printf("nentr %d\n",nentr);
	
   for (Int_t i=0; i<nentr; i++) {
      AliITSIOTrack *iotrack=new AliITSIOTrack;
      // tarray.AddAt(new AliITSIOTrack,i);
      // iotrack=(AliITSiotrack*)tarray.UncheckedAt(i);
       tbranch->SetAddress(&iotrack);
       tracktree->GetEvent(i);
	tarray.AddLast(iotrack);
   }
	//file->Close();		 
	
	  AliITSIOTrack *iotrack;
   for (Int_t i=0; i<nentr; i++) {
     AliITSIOTrack *iotrack=new AliITSIOTrack;      	
	 iotrack=(AliITSIOTrack*)tarray.UncheckedAt(i);
	 if(!iotrack) continue;
     Int_t labITS=iotrack->GetLabel();
     Int_t labTPC=iotrack->GetTPCLabel();   
	  Double_t phistate=iotrack->GetStatePhi();
	  Double_t tgl=iotrack->GetStateTgl();	
	  Double_t Zstate=iotrack->GetStateZ();
	  Double_t Dr=iotrack->GetStateD();		  
	  Double_t C=iotrack->GetStateC();
	  Double_t px=iotrack->GetPx();
	  Double_t py=iotrack->GetPy();	  
	  Double_t pz=iotrack->GetPz();	
	  Double_t pt=TMath::Sqrt(px*px+py*py);
	  Int_t c = iotrack->GetCharge(); 
	  Double_t x=iotrack->GetX();
	  Double_t y=iotrack->GetY();
	  Double_t z= iotrack->GetZ();
	  Double_t Dz=iotrack->GetDz(); 

	  
     // cout<<" track label = "<<label<<"\n";
     // cout<<" phi z D tanl C = "<<phistate<<" "<<Zstate<<" "<<Dr<<" "<<tgl<<" "<<C<<"\n"; 	  
	  
  //  Int_t ilab=TMath::Abs(iotrack->GetLabel());
    Int_t flaglab=0;   
    Int_t iii=0;
   Double_t ptg=0.,pxg=0.,pyg=0.,pzg=0.;
	Double_t xo=0., yo=0., zo=0.;  	 

    Int_t mingood=0, maxgood=0, jj=0;
    if(nev==evNumber1) mingood=0;
    else
    {for(jj=evNumber1; jj<nev; jj++) mingood += ntrackevtpc(jj);}
     for(jj=evNumber1; jj<=nev; jj++) maxgood+= ntrackevtpc(jj);    

   Int_t ilab=TMath::Abs(labTPC);
   // for (iii=0;iii<ngood;iii++) {
   for(iii=mingood; iii<maxgood; iii++){
	 //cout<<" ilab, gt[iii].lab = "<<ilab<<" "<<gt[iii].lab<<"\n"; getchar();
      if (ilab==gt[iii].lab) { 
	flaglab=1;
	ptg=gt[iii].ptg; 
	pxg=gt[iii].pxg;
	pyg=gt[iii].pyg;
	pzg=gt[iii].pzg;
	xo=gt[iii].x;
	yo=gt[iii].y;
	zo=gt[iii].z;			
	break;
      }
    }   

    if (!flaglab) continue;  
    //cout<<" j =  " <<j<<"\n";  getchar();  	  		  		    
   TVector dataOut(10);
    Int_t kkk=0;
    
    dataOut(kkk) = ptg; kkk++; dataOut(kkk)=labITS; kkk++; dataOut(kkk)=labTPC; kkk++;	
  
      ///////////////////////////////
      Double_t difpt= (pt-ptg)/ptg*100.;
      //cout<<" difpt = "<<difpt<<"\n"; getchar();  	            	                
      dataOut(kkk)=difpt; kkk++;                                             
      Double_t lambdag=TMath::ATan(pzg/ptg);
      Double_t   lam=TMath::ATan(tgl);      
      Double_t diflam = (lam - lambdag)*1000.;
      dataOut(kkk) = diflam; kkk++;	    	  			    
      Double_t phig=TMath::ATan2(pyg,pxg);  if(phig<0) phig=2.*TMath::Pi()+phig;       
     // Double_t phi=phivertex;
       Double_t phi=TMath::ACos(px/pt);
     Double_t duepi=2.*TMath::Pi();	 
     if(phi>duepi) phi-=duepi;
     if(phi<0.) phi+=duepi;      
      Double_t signC=0.; 
      if(c>0) signC=1.; else signC=-1.;
	  Double_t Dtot= TMath::Sqrt(Dr*Dr+Dz*Dz);       
      Double_t difphi = (phi - phig)*1000.;
      dataOut(kkk)=difphi; kkk++;
      dataOut(kkk)=Dtot*1.e4; kkk++;
      dataOut(kkk)=Dr*1.e4; kkk++;
      dataOut(kkk)=Dz*1.e4; kkk++;
      dataOut(kkk)=signC; kkk++; 
      Int_t r;
      for (r=0; r<10; r++) { out<<dataOut(r)<<" ";}
      out<<"\n";
      kkk=0;        	    

    delete iotrack;		 
   }  
  //out.close(); 
   }   // event loop 
  out.close();
  in.close();   
  file->Close();      
}
Commit	Line	Data
ae0cd643	1	#include "iostream.h"
	2
	3
	4	void AliITSComparisonV1(Int_t evNumber1=0,Int_t evNumber2=0) {
	5
	6	const char *filename="itstracks.root";
	7
	8	///////////////// Dynamically link some shared libs ////////////////////////////////
	9
	10	if (gClassTable->GetID("AliRun") < 0) {
	11	gROOT->LoadMacro("loadlibs.C");
	12	loadlibs();
	13	} else {
	14	delete gAlice;
	15	gAlice=0;
	16	}
	17
	18	// Connect the Root Galice file containing Geometry, Kine and Hits
	19	TFile file = (TFile)gROOT->GetListOfFiles()->FindObject(filename);
	20	if (!file) file = new TFile(filename);
	21
	22	ifstream in("itsgood_tracks");
	23	//
	24	// Loop over events
	25	//
	26
	27	char tname[30];
	28
	29
	30	struct GoodTrack {
	31	Int_t fEventN;
	32	Int_t lab,code;
	33	Float_t px,py,pz,x,y,z,pxg,pyg,pzg,ptg;
	34	Bool_t flag;
	35	};
	36	Int_t i;
	37	GoodTrack gt[15000];
	38	Int_t ngood=0;
	39	cerr<<"Reading itsgood tracks...\n";
	40	while (in>>gt[ngood].fEventN>>gt[ngood].lab>>gt[ngood].code
	41	>>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz
	42	>>gt[ngood].x >>gt[ngood].y >>gt[ngood].z
	43	>>gt[ngood].pxg >>gt[ngood].pyg >>gt[ngood].pzg
	44	>>gt[ngood].ptg >>gt[ngood].flag) {
	45	ngood++;
	46	cerr<<ngood<<'\r';
	47	if (ngood==15000) {
	48	cerr<<"Too many good tracks !\n";
	49	break;
	50	}
	51	}
	52	if (!in.eof()) cerr<<"Read error (itsgood_tracks) !\n";
	53	ofstream out1 ("AliITSTrag.out");
	54	Int_t countpos=0,countneg=0;
	55
	56
	57	//for (i=0; i<ngood; i++) out1 << gt[i].ptg << "\n";
	58	for (i=0; i<ngood; i++) {
	59	out1 << gt[i].ptg << "\n";
	60	Int_t codpar=gt[i].code;
	61	if(codpar==2212 \|\| codpar==-11 \|\| codpar==-13 \|\| codpar==211 \|\| codpar==321 \|\| codpar==3222
	62	\|\| codpar==213 \|\| codpar==323 \|\| codpar==10323 \|\| codpar==3224 \|\| codpar==2224 \|\| codpar==2214
	63	\|\| codpar==-1114 \|\| codpar==-3112 \|\| codpar==-3312 \|\| codpar==3224 \|\| codpar==-3114 \|\| codpar==-3314
	64	\|\| codpar==411 \|\| codpar==431 \|\| codpar==413 \|\| codpar==433 \|\| codpar==-15 \|\| codpar==4232
65	\|\| codpar==4222 \|\| codpar==4322 \|\| codpar==4422 \|\| codpar==4412 \|\| codpar==4432 \|\| codpar==4224
66	\|\|codpar==4214 \|\| codpar==4324 \|\| codpar==4424 \|\| codpar==4414 \|\| codpar==4434 \|\| codpar==4444)
67	countpos++;
68	if(codpar==-2212 \|\| codpar==11 \|\| codpar==13 \|\| codpar==-211 \|\| codpar==-321 \|\| codpar==3112
69	\|\| codpar==-213 \|\| codpar==-323 \|\| codpar==-10323 \|\| codpar==3114 \|\| codpar==1114 \|\| codpar==-2224
70	\|\| codpar==-2214 \|\| codpar==33112 \|\| codpar==-3222 \|\| codpar==3114 \|\| codpar==3314 \|\| codpar==3334
71	\|\| codpar==-3224 \|\| codpar==-411 \|\| codpar==-431 \|\| codpar==-413 \|\| codpar==-433 \|\| codpar==15
72	\|\| codpar==-4422 \|\| codpar==-4432 \|\| codpar==-4214 \|\| codpar==-4324 \|\| codpar==-4424 \|\| codpar==-4434
73	\|\| codpar==-444)
74	countneg++;
75	}
76	out1.close();
77	ofstream out ("AliITSTra.out");
78	/// definition of nm of good track for each event
79	TVector ntrackevtpc(evNumber2-evNumber1 +1);
80	Int_t nchange=-1;
81	Int_t nmev=-100;
82
83
84	Int_t i;
85	for( i =0; i<ngood ; i++){
86	if(gt[i].fEventN != nmev ){
87	nmev=gt[i].fEventN;
88	nchange++;
89	}
90	if(nmev == gt[i].fEventN) ntrackevtpc(nchange)++;
91	}
92	//////////////////////////////////////////////////////////////
93
94	for (Int_t nev=evNumber1; nev<= evNumber2; nev++) {
95
96
97	sprintf(tname,"TreeT%d",nev);
98	TTree tracktree=(TTree)file->Get(tname);
99	TBranch *tbranch=tracktree->GetBranch("ITStracks");
100	cout<<" nev = "<<nev<<"\n";
101	//cout<<" open the file \n";
102
103	Int_t nentr=tracktree->GetEntries();
104
105	TObjArray tarray(nentr);
106	// AliITSIOTrack *iotrack=0;
107	printf("nentr %d\n",nentr);
108
109	for (Int_t i=0; i<nentr; i++) {
110	AliITSIOTrack *iotrack=new AliITSIOTrack;
111	// tarray.AddAt(new AliITSIOTrack,i);
112	// iotrack=(AliITSiotrack*)tarray.UncheckedAt(i);
113	tbranch->SetAddress(&iotrack);
114	tracktree->GetEvent(i);
115	tarray.AddLast(iotrack);
116	}
117	//file->Close();
118
119	AliITSIOTrack *iotrack;
120	for (Int_t i=0; i<nentr; i++) {
121	AliITSIOTrack *iotrack=new AliITSIOTrack;
122	iotrack=(AliITSIOTrack*)tarray.UncheckedAt(i);
123	if(!iotrack) continue;
124	Int_t labITS=iotrack->GetLabel();
61231d74	125	Int_t labTPC=iotrack->GetTPCLabel();
ae0cd643	126	Double_t phistate=iotrack->GetStatePhi();
	127	Double_t tgl=iotrack->GetStateTgl();
	128	Double_t Zstate=iotrack->GetStateZ();
	129	Double_t Dr=iotrack->GetStateD();
	130	Double_t C=iotrack->GetStateC();
	131	Double_t px=iotrack->GetPx();
	132	Double_t py=iotrack->GetPy();
	133	Double_t pz=iotrack->GetPz();
	134	Double_t pt=TMath::Sqrt(pxpx+pypy);
	135	Int_t c = iotrack->GetCharge();
	136	Double_t x=iotrack->GetX();
	137	Double_t y=iotrack->GetY();
61231d74	138	Double_t z= iotrack->GetZ();
	139	Double_t Dz=iotrack->GetDz();
	140
ae0cd643	141
	142	// cout<<" track label = "<<label<<"\n";
	143	// cout<<" phi z D tanl C = "<<phistate<<" "<<Zstate<<" "<<Dr<<" "<<tgl<<" "<<C<<"\n";
	144
	145	// Int_t ilab=TMath::Abs(iotrack->GetLabel());
	146	Int_t flaglab=0;
	147	Int_t iii=0;
335d71a8	148	Double_t ptg=0.,pxg=0.,pyg=0.,pzg=0.;
335d71a8	149	Double_t xo=0., yo=0., zo=0.;
ae0cd643	150
	151	Int_t mingood=0, maxgood=0, jj=0;
	152	if(nev==evNumber1) mingood=0;
	153	else
	154	{for(jj=evNumber1; jj<nev; jj++) mingood += ntrackevtpc(jj);}
	155	for(jj=evNumber1; jj<=nev; jj++) maxgood+= ntrackevtpc(jj);
	156
	157	Int_t ilab=TMath::Abs(labTPC);
	158	// for (iii=0;iii<ngood;iii++) {
	159	for(iii=mingood; iii<maxgood; iii++){
	160	//cout<<" ilab, gt[iii].lab = "<<ilab<<" "<<gt[iii].lab<<"\n"; getchar();
	161	if (ilab==gt[iii].lab) {
	162	flaglab=1;
	163	ptg=gt[iii].ptg;
	164	pxg=gt[iii].pxg;
	165	pyg=gt[iii].pyg;
335d71a8	166	pzg=gt[iii].pzg;
	167	xo=gt[iii].x;
	168	yo=gt[iii].y;
	169	zo=gt[iii].z;
ae0cd643	170	break;
	171	}
	172	}
	173
	174	if (!flaglab) continue;
	175	//cout<<" j = " <<j<<"\n"; getchar();
	176	TVector dataOut(10);
	177	Int_t kkk=0;
	178
	179	dataOut(kkk) = ptg; kkk++; dataOut(kkk)=labITS; kkk++; dataOut(kkk)=labTPC; kkk++;
	180
	181	///////////////////////////////
	182	Double_t difpt= (pt-ptg)/ptg*100.;
	183	//cout<<" difpt = "<<difpt<<"\n"; getchar();
	184	dataOut(kkk)=difpt; kkk++;
	185	Double_t lambdag=TMath::ATan(pzg/ptg);
	186	Double_t lam=TMath::ATan(tgl);
	187	Double_t diflam = (lam - lambdag)*1000.;
	188	dataOut(kkk) = diflam; kkk++;
	189	Double_t phig=TMath::ATan2(pyg,pxg); if(phig<0) phig=2.*TMath::Pi()+phig;
	190	// Double_t phi=phivertex;
	191	Double_t phi=TMath::ACos(px/pt);
	192	Double_t duepi=2.*TMath::Pi();
	193	if(phi>duepi) phi-=duepi;
	194	if(phi<0.) phi+=duepi;
	195	Double_t signC=0.;
	196	if(c>0) signC=1.; else signC=-1.;
335d71a8	197	Double_t Dtot= TMath::Sqrt(DrDr+DzDz);
ae0cd643	198	Double_t difphi = (phi - phig)*1000.;
	199	dataOut(kkk)=difphi; kkk++;
	200	dataOut(kkk)=Dtot*1.e4; kkk++;
	201	dataOut(kkk)=Dr*1.e4; kkk++;
	202	dataOut(kkk)=Dz*1.e4; kkk++;
	203	dataOut(kkk)=signC; kkk++;
	204	Int_t r;
	205	for (r=0; r<10; r++) { out<<dataOut(r)<<" ";}
	206	out<<"\n";
	207	kkk=0;
	208
	209	delete iotrack;
	210	}
	211	//out.close();
	212	} // event loop
	213	out.close();
	214	in.close();
	215	file->Close();
	216	}
	217