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

//-----------------------------------------------------//
//                                                     //
//           Date   : August 05 2003                   //
//  This reads the file PMD.digits.root(TreeD),        //
//  calls the Clustering algorithm and stores the      //
//  clustering output in PMD.RecPoints.root(TreeR)     // 
//                                                     //
//-----------------------------------------------------//

#include <Riostream.h>
#include <TMath.h>
#include <TBRIK.h>
#include <TNode.h>
#include <TTree.h>
#include <TGeometry.h>
#include <TObjArray.h>
#include <TClonesArray.h>
#include <TFile.h>
#include <TNtuple.h>
#include <TParticle.h>

#include "AliPMDcluster.h"
#include "AliPMDclupid.h"
#include "AliPMDDiscriminator.h"

ClassImp(AliPMDDiscriminator)

AliPMDDiscriminator::AliPMDDiscriminator() :
  fDiscrim(0),
  fDebug(0)
{
//
// Default Constructor
//

}
// -----------------------------------------------------------------------
AliPMDDiscriminator::~AliPMDDiscriminator()
{
  // Destructor
}

// -----------------------------------------------------------------------

void AliPMDDiscriminator::Discrimination(TObjArray *pmdcontin, TObjArray *pmdcontout)
{
  // Does Photon/Hadron discrimination

  if(fDiscrim == 0)
    {
      EmpDiscrimination(pmdcontin, pmdcontout);
    }
  else if(fDiscrim == 1)
    {
      NNDiscrimination();
    }
}
// -----------------------------------------------------------------------

void AliPMDDiscriminator::EmpDiscrimination(TObjArray *pmdcontin, TObjArray *pmdcontout)
{
  // Does Photon/Hadron discrimination
  // matching the clusters of CPV and PREshower plane
  //
  const  Int_t kumperdet = 24;
  static Int_t neibx[6]={1,0,-1,-1,0,1}, neiby[6]={0,1,1,0,-1,-1}; 
  Int_t   det,smn;
  Int_t   iprecount[24], icpvcount[24];
  Float_t xpos,ypos;
  Float_t adc, ncell, rad;
  Float_t clusdata[6];

  for(Int_t i = 0; i < kumperdet; i++)
    {
      iprecount[i] = 0;
      icpvcount[i] = 0;
    }
  AliPMDcluster  *pmdcl    = 0;
  AliPMDclupid   *pmdclout = 0;

  Int_t nentries1 = pmdcontin->GetEntries();
  cout << " nentries1 = " << nentries1 << endl;
  for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
    {
      pmdcl = (AliPMDcluster*)pmdcontin->UncheckedAt(ient1);

      det   = pmdcl->GetDetector();
      smn   = pmdcl->GetSMN();
      if(det == 0) iprecount[smn]++;
      if(det == 1) icpvcount[smn]++;
    } // Entries of TObjArray loop

  Int_t   idet, ismn;
  Float_t edepcpv[48][96];
  Int_t   statuscpv[48][96];

  for(Int_t i = 0; i < kumperdet; i++) // unit module
    {
      // Initialisation of the ADC of CPV (1UM)
      for (Int_t ix = 0; ix < 48;ix++)
	{
	  for (Int_t iy = 0; iy < 96;iy++)
	    {
	      edepcpv[ix][iy]   = 0.;
	      statuscpv[ix][iy] = 0;
	    }
	}
      Int_t precounter   = iprecount[i];
      Int_t cpvcounter   = icpvcount[i];

      Float_t *xpadpre   = new Float_t[precounter];
      Float_t *ypadpre   = new Float_t[precounter];
      Float_t *adcpre    = new Float_t[precounter];
      Float_t *ncellpre  = new Float_t[precounter];
      Float_t *radpre    = new Float_t[precounter];
      Int_t   *sortcoord = new Int_t[precounter];
      Int_t   *clupidpre = new Int_t[precounter];

      Float_t *xpadcpv   = new Float_t[cpvcounter];
      Float_t *ypadcpv   = new Float_t[cpvcounter];
      Float_t *adccpv    = new Float_t[cpvcounter];
      Float_t *ncellcpv  = new Float_t[cpvcounter];
      Float_t *radcpv    = new Float_t[cpvcounter];

      Int_t ii = 0;
      Int_t ij = 0;
      for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
	{
	  pmdcl = (AliPMDcluster*)pmdcontin->UncheckedAt(ient1);
	  
	  det   = pmdcl->GetDetector();
	  smn   = pmdcl->GetSMN();
	  xpos  = pmdcl->GetClusX();
	  ypos  = pmdcl->GetClusY();
	  adc   = pmdcl->GetClusADC();
	  ncell = pmdcl->GetClusCells();
	  rad   = pmdcl->GetClusRadius();
	  
	  if(det == 0 && smn == i)
	    {
	      xpadpre[ii]  = xpos;
	      ypadpre[ii]  = ypos;
	      adcpre[ii]   = adc;
	      ncellpre[ii] = ncell;
	      radpre[ii]   = rad;
	      ii++;
	    }
	  if(det == 1 && smn == i)
	    {
	      Int_t ix = (Int_t) (xpos+0.5);
	      Int_t iy = (Int_t) (ypos+0.5);
	      if(ix > 47) ix = 47;
	      if(iy > 95) iy = 95;
	      edepcpv[ix][iy] = adc;
	      xpadcpv[ij]  = xpos;
	      ypadcpv[ij]  = ypos;
	      adccpv[ij]   = adc;
	      ncellcpv[ij] = ncell;
	      radcpv[ij]   = rad;
	      ij++;
	    }
	} // Entries of TObjArray loop
      // sorts from lowest ADC to highest ADC
      // and returns the coordinates
      Bool_t jsort = false;
      TMath::Sort(precounter,adcpre,sortcoord,jsort);

      Int_t jjsort = 0;
      for(Int_t jj=0; jj<precounter; jj++)
	{
	  // PRE information
	  // PIDs for PRE clusters are 0(photon) and 1(hadron)

	  jjsort = sortcoord[jj];

	  Int_t ix = (Int_t) (xpadpre[jjsort]+0.5);
	  Int_t iy = (Int_t) (ypadpre[jjsort]+0.5);
	  if(ix > 47) ix = 47;
	  if(iy > 95) iy = 95;

	  for(Int_t jk=0; jk<6; jk++)
	    {
	      Int_t jd1 = ix + neibx[jk]; 
	      Int_t jd2 = iy + neiby[jk];
	      if(jd1 <0 ) jd1 = 0;
	      if(jd1 >47) jd1 = 47;
	      if(jd2 <0 ) jd2 = 0;
	      if(jd2 >47) jd2 = 47;
	      if(edepcpv[jd1][jd2] > 0.0 && statuscpv[jd1][jd2] == 0)
		{
		  statuscpv[jd1][jd2] = 1;
		  clupidpre[jjsort]   = 1;
		  break;
		}
	    }

	  idet        = 0;
	  ismn        = i;
	  clusdata[0] = xpadpre[jjsort];
	  clusdata[1] = ypadpre[jjsort];
	  clusdata[2] = adcpre[jjsort];
	  clusdata[3] = ncellpre[jjsort];
	  clusdata[4] = radpre[jjsort];
	  //PH	  clusdata[5] = (Float_t) clupidpre[jjsort];

	  // Temporary the cluster PID is set to 1 if the
	  // adc > 3MIP units which will be changed later on.

	  if (adcpre[jjsort] > 4500.)
	    {
	      clusdata[5] = 1.0;
	    }
	  else
	    {
	      clusdata[5] = 0.0;
	    }
	  pmdclout = new AliPMDclupid(idet,ismn,clusdata);
	  pmdcontout->Add(pmdclout);
	} // jj loop

      for(Int_t jj=0; jj<cpvcounter; jj++)
	{
	  // CPV information
	  // PID for CPV clusters is 1

	  idet        = 1;
	  ismn        = i;
	  clusdata[0] = xpadcpv[jj];
	  clusdata[1] = ypadcpv[jj];
	  clusdata[2] = adccpv[jj];
	  clusdata[3] = ncellcpv[jj];
	  clusdata[4] = radcpv[jj];
	  clusdata[5] = 0.;

	  pmdclout = new AliPMDclupid(idet,ismn,clusdata);
	  pmdcontout->Add(pmdclout);
	}

      // delete all the pointers
      delete [] xpadpre;
      delete [] ypadpre;
      delete [] adcpre;
      delete [] ncellpre;
      delete [] radpre;
      delete [] clupidpre;
      delete [] xpadcpv;
      delete [] ypadcpv;
      delete [] adccpv;
      delete [] ncellcpv;
      delete [] radcpv;
    } // i loop

}
// -----------------------------------------------------------------------
void AliPMDDiscriminator::NNDiscrimination()
{
  // This method does discrimination using Neural Network technique

}
// -----------------------------------------------------------------------
void AliPMDDiscriminator::SetDiscrimination(Int_t idiscrim)
{
  fDiscrim = idiscrim;
}
// -----------------------------------------------------------------------
Commit	Line	Data
7dcaf913	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	// Date : August 05 2003 //
	19	// This reads the file PMD.digits.root(TreeD), //
	20	// calls the Clustering algorithm and stores the //
	21	// clustering output in PMD.RecPoints.root(TreeR) //
	22	// //
	23	//-----------------------------------------------------//
	24
	25	#include <Riostream.h>
	26	#include <TMath.h>
	27	#include <TBRIK.h>
	28	#include <TNode.h>
	29	#include <TTree.h>
	30	#include <TGeometry.h>
	31	#include <TObjArray.h>
	32	#include <TClonesArray.h>
	33	#include <TFile.h>
	34	#include <TNtuple.h>
	35	#include <TParticle.h>
	36
	37	#include "AliPMDcluster.h"
	38	#include "AliPMDclupid.h"
	39	#include "AliPMDDiscriminator.h"
	40
	41	ClassImp(AliPMDDiscriminator)
	42
	43	AliPMDDiscriminator::AliPMDDiscriminator() :
	44	fDiscrim(0),
	45	fDebug(0)
	46	{
	47	//
	48	// Default Constructor
	49	//
	50
	51	}
	52	// -----------------------------------------------------------------------
	53	AliPMDDiscriminator::~AliPMDDiscriminator()
	54	{
	55	// Destructor
	56	}
	57
	58	// -----------------------------------------------------------------------
	59
	60	void AliPMDDiscriminator::Discrimination(TObjArray pmdcontin, TObjArray pmdcontout)
	61	{
	62	// Does Photon/Hadron discrimination
	63
	64	if(fDiscrim == 0)
65	{
66	EmpDiscrimination(pmdcontin, pmdcontout);
67	}
68	else if(fDiscrim == 1)
69	{
70	NNDiscrimination();
71	}
72	}
73	// -----------------------------------------------------------------------
74
75	void AliPMDDiscriminator::EmpDiscrimination(TObjArray pmdcontin, TObjArray pmdcontout)
76	{
77	// Does Photon/Hadron discrimination
78	// matching the clusters of CPV and PREshower plane
79	//
80	const Int_t kumperdet = 24;
81	static Int_t neibx[6]={1,0,-1,-1,0,1}, neiby[6]={0,1,1,0,-1,-1};
82	Int_t det,smn;
83	Int_t iprecount[24], icpvcount[24];
84	Float_t xpos,ypos;
85	Float_t adc, ncell, rad;
86	Float_t clusdata[6];
87
88	for(Int_t i = 0; i < kumperdet; i++)
89	{
90	iprecount[i] = 0;
91	icpvcount[i] = 0;
92	}
93	AliPMDcluster *pmdcl = 0;
94	AliPMDclupid *pmdclout = 0;
95
96	Int_t nentries1 = pmdcontin->GetEntries();
97	cout << " nentries1 = " << nentries1 << endl;
98	for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
99	{
100	pmdcl = (AliPMDcluster*)pmdcontin->UncheckedAt(ient1);
101
102	det = pmdcl->GetDetector();
103	smn = pmdcl->GetSMN();
104	if(det == 0) iprecount[smn]++;
105	if(det == 1) icpvcount[smn]++;
106	} // Entries of TObjArray loop
107
108	Int_t idet, ismn;
109	Float_t edepcpv[48][96];
110	Int_t statuscpv[48][96];
111
112	for(Int_t i = 0; i < kumperdet; i++) // unit module
113	{
114	// Initialisation of the ADC of CPV (1UM)
115	for (Int_t ix = 0; ix < 48;ix++)
116	{
117	for (Int_t iy = 0; iy < 96;iy++)
118	{
119	edepcpv[ix][iy] = 0.;
120	statuscpv[ix][iy] = 0;
121	}
122	}
123	Int_t precounter = iprecount[i];
124	Int_t cpvcounter = icpvcount[i];
125
126	Float_t *xpadpre = new Float_t[precounter];
127	Float_t *ypadpre = new Float_t[precounter];
128	Float_t *adcpre = new Float_t[precounter];
129	Float_t *ncellpre = new Float_t[precounter];
130	Float_t *radpre = new Float_t[precounter];
131	Int_t *sortcoord = new Int_t[precounter];
132	Int_t *clupidpre = new Int_t[precounter];
133
134	Float_t *xpadcpv = new Float_t[cpvcounter];
135	Float_t *ypadcpv = new Float_t[cpvcounter];
136	Float_t *adccpv = new Float_t[cpvcounter];
137	Float_t *ncellcpv = new Float_t[cpvcounter];
138	Float_t *radcpv = new Float_t[cpvcounter];
139
140	Int_t ii = 0;
141	Int_t ij = 0;
142	for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
143	{
144	pmdcl = (AliPMDcluster*)pmdcontin->UncheckedAt(ient1);
145
146	det = pmdcl->GetDetector();
147	smn = pmdcl->GetSMN();
148	xpos = pmdcl->GetClusX();
149	ypos = pmdcl->GetClusY();
150	adc = pmdcl->GetClusADC();
151	ncell = pmdcl->GetClusCells();
152	rad = pmdcl->GetClusRadius();
153
154	if(det == 0 && smn == i)
155	{
156	xpadpre[ii] = xpos;
157	ypadpre[ii] = ypos;
158	adcpre[ii] = adc;
159	ncellpre[ii] = ncell;
160	radpre[ii] = rad;
161	ii++;
162	}
163	if(det == 1 && smn == i)
164	{
165	Int_t ix = (Int_t) (xpos+0.5);
166	Int_t iy = (Int_t) (ypos+0.5);
167	if(ix > 47) ix = 47;
168	if(iy > 95) iy = 95;
169	edepcpv[ix][iy] = adc;
170	xpadcpv[ij] = xpos;
171	ypadcpv[ij] = ypos;
172	adccpv[ij] = adc;
173	ncellcpv[ij] = ncell;
174	radcpv[ij] = rad;
175	ij++;
176	}
177	} // Entries of TObjArray loop
178	// sorts from lowest ADC to highest ADC
179	// and returns the coordinates
180	Bool_t jsort = false;
181	TMath::Sort(precounter,adcpre,sortcoord,jsort);
182
183	Int_t jjsort = 0;
184	for(Int_t jj=0; jj<precounter; jj++)
185	{
186	// PRE information
187	// PIDs for PRE clusters are 0(photon) and 1(hadron)
188
189	jjsort = sortcoord[jj];
190
191	Int_t ix = (Int_t) (xpadpre[jjsort]+0.5);
192	Int_t iy = (Int_t) (ypadpre[jjsort]+0.5);
193	if(ix > 47) ix = 47;
194	if(iy > 95) iy = 95;
195
196	for(Int_t jk=0; jk<6; jk++)
197	{
198	Int_t jd1 = ix + neibx[jk];
199	Int_t jd2 = iy + neiby[jk];
200	if(jd1 <0 ) jd1 = 0;
201	if(jd1 >47) jd1 = 47;
202	if(jd2 <0 ) jd2 = 0;
203	if(jd2 >47) jd2 = 47;
204	if(edepcpv[jd1][jd2] > 0.0 && statuscpv[jd1][jd2] == 0)
205	{
206	statuscpv[jd1][jd2] = 1;
207	clupidpre[jjsort] = 1;
208	break;
209	}
210	}
211
212	idet = 0;
213	ismn = i;
214	clusdata[0] = xpadpre[jjsort];
215	clusdata[1] = ypadpre[jjsort];
216	clusdata[2] = adcpre[jjsort];
217	clusdata[3] = ncellpre[jjsort];
218	clusdata[4] = radpre[jjsort];
278747fb	219	//PH clusdata[5] = (Float_t) clupidpre[jjsort];
7dcaf913	220
	221	// Temporary the cluster PID is set to 1 if the
	222	// adc > 3MIP units which will be changed later on.
	223
	224	if (adcpre[jjsort] > 4500.)
	225	{
	226	clusdata[5] = 1.0;
	227	}
	228	else
	229	{
	230	clusdata[5] = 0.0;
	231	}
	232	pmdclout = new AliPMDclupid(idet,ismn,clusdata);
	233	pmdcontout->Add(pmdclout);
	234	} // jj loop
	235
	236	for(Int_t jj=0; jj<cpvcounter; jj++)
	237	{
	238	// CPV information
	239	// PID for CPV clusters is 1
	240
	241	idet = 1;
	242	ismn = i;
	243	clusdata[0] = xpadcpv[jj];
	244	clusdata[1] = ypadcpv[jj];
	245	clusdata[2] = adccpv[jj];
	246	clusdata[3] = ncellcpv[jj];
	247	clusdata[4] = radcpv[jj];
	248	clusdata[5] = 0.;
	249
	250	pmdclout = new AliPMDclupid(idet,ismn,clusdata);
	251	pmdcontout->Add(pmdclout);
	252	}
	253
	254	// delete all the pointers
	255	delete [] xpadpre;
	256	delete [] ypadpre;
	257	delete [] adcpre;
	258	delete [] ncellpre;
	259	delete [] radpre;
	260	delete [] clupidpre;
	261	delete [] xpadcpv;
	262	delete [] ypadcpv;
	263	delete [] adccpv;
	264	delete [] ncellcpv;
	265	delete [] radcpv;
	266	} // i loop
	267
	268	}
	269	// -----------------------------------------------------------------------
	270	void AliPMDDiscriminator::NNDiscrimination()
	271	{
	272	// This method does discrimination using Neural Network technique
	273
	274	}
	275	// -----------------------------------------------------------------------
	276	void AliPMDDiscriminator::SetDiscrimination(Int_t idiscrim)
	277	{
	278	fDiscrim = idiscrim;
	279	}
	280	// -----------------------------------------------------------------------
	281