FEDRA emulsion software from the OPERA Collaboration
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
EdbViewCell Class Reference

#include <EdbViewRec.h>

Inheritance diagram for EdbViewCell:
Collaboration diagram for EdbViewCell:

Public Member Functions

int AddCluster (EdbCluster *c)
 
void CalcN ()
 
void CalcStat ()
 
void CleanCell ()
 
void Delete ()
 
 EdbViewCell ()
 
int FillCell (TClonesArray &v)
 
EdbCluster ** GetCell (int j) const
 
void Init ()
 
void InitMem ()
 
int IXcell (float x) const
 
int IYcell (float y) const
 
int Jcell (float x, float y, float z) const
 
int Jcell (float x, float y, int ifr) const
 
int Jcell (int ix, int iy, int ifr) const
 
int Jcell (int ixy, int ifr) const
 
int JcellXY (float x, float y) const
 
int Jneib (int i) const
 
void Print ()
 
void SetBin (float binx, float biny, float binz=-1)
 
void SetCellLimits (int ncell, int ncl)
 
void SetLimits (float xmin, float xmax, float ymin, float ymax, float zmin, float zmax)
 
void SetNfr (int nfr, float zmin, float zmax, int ifz=0)
 
virtual ~EdbViewCell ()
 

Public Attributes

Float_t eBinX
 
Float_t eBinY
 
Float_t eBinZ
 
Int_t eIFZ
 
Int_t * eNC
 
Int_t eNcell
 
Int_t eNcellXY
 
Int_t eNcl
 
Int_t eNfr
 
Int_t eNx
 
Int_t eNy
 
Float_t eXmax
 
Float_t eXmin
 
Float_t eYmax
 
Float_t eYmin
 
Float_t eZmax
 
Float_t eZmin
 

Private Attributes

Int_t eCellLim
 
Int_t * eFrame
 pointers to cells [eNcellsLim] More...
 
Int_t eNcellsLim
 [eNcell] number of clusters/cell More...
 
Int_t eNeib [9]
 index of the first cell of the frame in epCell array More...
 
EdbCluster ** epC
 
EdbCluster *** epCell
 pointers to clusters [eNcellsLim*eCellLim] More...
 

Constructor & Destructor Documentation

◆ EdbViewCell()

EdbViewCell::EdbViewCell ( )
28{
29 epC=0;
30 epCell=0;
31 eNcellsLim=0;
32 eCellLim=0;
33 eFrame=0;
34 eNC=0;
35
36 eNcell=eNfr=0;
40}
Int_t * eNC
Definition: EdbViewRec.h:40
Int_t eNx
Definition: EdbViewRec.h:31
Float_t eBinY
Definition: EdbViewRec.h:28
Float_t eZmin
Definition: EdbViewRec.h:25
Float_t eBinX
Definition: EdbViewRec.h:27
Float_t eBinZ
Definition: EdbViewRec.h:26
Int_t eNfr
Definition: EdbViewRec.h:29
Int_t eNcellsLim
[eNcell] number of clusters/cell
Definition: EdbViewRec.h:43
Float_t eYmax
Definition: EdbViewRec.h:24
Float_t eYmin
Definition: EdbViewRec.h:24
Float_t eZmax
Definition: EdbViewRec.h:25
Int_t eNcell
Definition: EdbViewRec.h:34
EdbCluster ** epC
Definition: EdbViewRec.h:46
Float_t eXmax
Definition: EdbViewRec.h:23
Int_t eNy
Definition: EdbViewRec.h:32
Int_t * eFrame
pointers to cells [eNcellsLim]
Definition: EdbViewRec.h:48
Float_t eXmin
Definition: EdbViewRec.h:23
Int_t eCellLim
Definition: EdbViewRec.h:44
Int_t eNcl
Definition: EdbViewRec.h:35
EdbCluster *** epCell
pointers to clusters [eNcellsLim*eCellLim]
Definition: EdbViewRec.h:47
Int_t eNcellXY
Definition: EdbViewRec.h:33

◆ ~EdbViewCell()

EdbViewCell::~EdbViewCell ( )
virtual
44{
45 //printf("EdbViewCell::~EdbViewCell() \n");
46 Delete();
47}
void Delete()
Definition: EdbViewRec.cxx:50

Member Function Documentation

◆ AddCluster()

int EdbViewCell::AddCluster ( EdbCluster c)
inline
68 {
69 int j;
70 if(!eIFZ) j = Jcell(c->eX,c->eY,c->eFrame);
71 else j = Jcell(c->eX,c->eY,c->eZ);
72 printf("j=%d eNC[j]=%d\n", j, eNC[j]);
73 if( j<0 || j>eNcellsLim ) return 0;
74 if( eNC[j] >= eCellLim ) return 0;
75 epCell[j][eNC[j]] = c;
76 eNC[j]++;
77 return 1;
78 }
Float_t eZ
cluster coordinates in pixels(?)
Definition: EdbCluster.h:25
Int_t eFrame
frame index
Definition: EdbCluster.h:28
Float_t eY
cluster coordinates in pixels(?)
Definition: EdbCluster.h:24
Float_t eX
cluster coordinates in pixels(?)
Definition: EdbCluster.h:23
int Jcell(float x, float y, int ifr) const
Definition: EdbViewRec.h:89
Int_t eIFZ
Definition: EdbViewRec.h:37

◆ CalcN()

void EdbViewCell::CalcN ( )
60{
61 // just for eNx, eNy calculation
62 eNx = (int)((eXmax-eXmin)/eBinX)+1;
63 eNy = (int)((eYmax-eYmin)/eBinY)+1;
66
67 eNeib[0]=0;
68 eNeib[1]=1; eNeib[2]=-1; eNeib[3]=eNx; eNeib[4]=-eNx;
69 eNeib[5]=eNx+1; eNeib[6]=eNx-1; eNeib[7]=-eNx+1; eNeib[8]=-eNx-1;
70}
Int_t eNeib[9]
index of the first cell of the frame in epCell array
Definition: EdbViewRec.h:50

◆ CalcStat()

void EdbViewCell::CalcStat ( )
146{
147 int nmax=0;
148 for(int i=0; i<eNcell; i++) if( eNC[i]>nmax ) nmax=eNC[i];
149 printf("max cell: %d \n",nmax);
150}

◆ CleanCell()

void EdbViewCell::CleanCell ( )
116{
117 memset(eNC,'\0',eNcellsLim*sizeof(Int_t));
118}

◆ Delete()

void EdbViewCell::Delete ( )
51{
52 if(eFrame) { delete [] eFrame; eFrame=0; }
53 if(epCell) { delete [] epCell; epCell=0; }
54 if(epC) { delete [] epC; epC = 0; }
55 if(eNC) { delete [] eNC; eNC = 0; }
56}

◆ FillCell()

int EdbViewCell::FillCell ( TClonesArray &  v)
122{
123 int cnt=0;
124 eNcl=v.GetLast()+1;
125 if(eNcl<1) return 0;
126
127 CleanCell();
128
129 printf("Fill %d (%d) Cells with %d clusters... ",eNcell,eNcellsLim,eNcl);
130 EdbCluster *c=0;
131 for(int i=0; i<eNcl; i++) {
132 //printf("\n 1 ");
133 c = (EdbCluster*)v.UncheckedAt(i);
134 //printf("\n 2 ");
135 //c->Print();
136 if( c->eX<eXmin || c->eX>eXmax || c->eY<eYmin || c->eY>eYmax) continue;
137 AddCluster(c);
138 cnt++;
139 }
140 printf(" %d are accepted\n",cnt);
141 return cnt;
142}
Definition: EdbCluster.h:19
int AddCluster(EdbCluster *c)
Definition: EdbViewRec.h:67
void CleanCell()
Definition: EdbViewRec.cxx:115

◆ GetCell()

EdbCluster ** EdbViewCell::GetCell ( int  j) const
inline
95{return epCell[j];}

◆ Init()

void EdbViewCell::Init ( void  )
106{
107 // initialization and memory allocation
108 CalcN();
109 InitMem();
110 printf("EdbViewCell::Init: \n");
111 Print();
112}
void Print()
Definition: EdbViewRec.cxx:153
void InitMem()
Definition: EdbViewRec.cxx:91
void CalcN()
Definition: EdbViewRec.cxx:59

◆ InitMem()

void EdbViewCell::InitMem ( )
92{
93 // memory allocation
94
95 Delete();
96 epC = new EdbCluster*[eNcellsLim*eCellLim]; //- pointers to clusters
97 epCell = new EdbCluster**[eNcellsLim]; //- pointers to cells
98 EdbCluster **pc = epC;
99 epCell[0]=pc;
100 for(int i=1; i<eNcellsLim; i++) {pc+=eCellLim; epCell[i]=pc;}
101 eNC = new Int_t[eNcellsLim];
102}

◆ IXcell()

int EdbViewCell::IXcell ( float  x) const
inline
86{return (int)((x-eXmin)/eBinX); }

◆ IYcell()

int EdbViewCell::IYcell ( float  y) const
inline
87{return (int)((y-eYmin)/eBinY); }

◆ Jcell() [1/4]

int EdbViewCell::Jcell ( float  x,
float  y,
float  z 
) const
inline
90{return TMath::Nint((z-eZmin)/eBinZ)*eNcellXY + IYcell(y)*eNx+IXcell(x);}
int IYcell(float y) const
Definition: EdbViewRec.h:87
int IXcell(float x) const
Definition: EdbViewRec.h:86

◆ Jcell() [2/4]

int EdbViewCell::Jcell ( float  x,
float  y,
int  ifr 
) const
inline
89{return eFrame[ifr] + IYcell(y)*eNx+IXcell(x); }

◆ Jcell() [3/4]

int EdbViewCell::Jcell ( int  ix,
int  iy,
int  ifr 
) const
inline
92{return eFrame[ifr]+iy*eNx + ix;}

◆ Jcell() [4/4]

int EdbViewCell::Jcell ( int  ixy,
int  ifr 
) const
inline
91{return eFrame[ifr]+ixy;}

◆ JcellXY()

int EdbViewCell::JcellXY ( float  x,
float  y 
) const
inline
88{return IYcell(y)*eNx+IXcell(x); }

◆ Jneib()

int EdbViewCell::Jneib ( int  i) const
inline
93{return eNeib[i];}

◆ Print()

void EdbViewCell::Print ( )
154{
155 printf("\n");
156 printf("eNfr = %d\n",eNfr);
157 printf("eNcl = %d\n",eNcl);
158 printf("eNx: (%f - %f)/ %f = %d \n", eXmax,eXmin,eBinX, eNx);
159 printf("eNy: (%f - %f)/ %f = %d \n", eYmax,eYmin,eBinY, eNy);
160 printf("eNz: (%f - %f)/ %f = %d \n", eZmax,eZmin,eBinZ, (int)((eZmax-eZmin)/eBinZ));
161 printf("eNcell = eNfr*eNcellXY: %d = %d * %d \n",eNcell,eNfr,eNcellXY);
162 printf("\n");
163}

◆ SetBin()

void EdbViewCell::SetBin ( float  binx,
float  biny,
float  binz = -1 
)
inline
59 {eBinX=binx; eBinY=biny; eBinZ=binz;}

◆ SetCellLimits()

void EdbViewCell::SetCellLimits ( int  ncell,
int  ncl 
)
inline
60{ eNcellsLim=ncell; eCellLim=ncl; }

◆ SetLimits()

void EdbViewCell::SetLimits ( float  xmin,
float  xmax,
float  ymin,
float  ymax,
float  zmin,
float  zmax 
)
inline
57 { eXmin=xmin; eXmax=xmax; eYmin=ymin; eYmax=ymax; eZmin=zmin; eZmax=zmax; }

◆ SetNfr()

void EdbViewCell::SetNfr ( int  nfr,
float  zmin,
float  zmax,
int  ifz = 0 
)
74{
75 eNfr=nfr;
76 eZmin=zmin;
77 eZmax=zmax;
78 eIFZ=ifz;
79 if(ifz) eBinZ = (zmax-zmin)/(nfr-1);
80 else eBinZ=0;
82
83 if(eFrame) delete [] eFrame;
84 eFrame = new Int_t[eNfr];
85 for(int i=0; i<eNfr; i++) eFrame[i]=i*eNcellXY;
86
87 printf("SetNfr: eNfr=%d; eZmin=%f; eZmax=%f; eBinZ=%f eNcell=%d\n", eNfr, eZmin, eZmax, eBinZ, eNcell);
88}

Member Data Documentation

◆ eBinX

Float_t EdbViewCell::eBinX

◆ eBinY

Float_t EdbViewCell::eBinY

◆ eBinZ

Float_t EdbViewCell::eBinZ

◆ eCellLim

Int_t EdbViewCell::eCellLim
private

◆ eFrame

Int_t* EdbViewCell::eFrame
private

pointers to cells [eNcellsLim]

◆ eIFZ

Int_t EdbViewCell::eIFZ

◆ eNC

Int_t* EdbViewCell::eNC

◆ eNcell

Int_t EdbViewCell::eNcell

◆ eNcellsLim

Int_t EdbViewCell::eNcellsLim
private

[eNcell] number of clusters/cell

◆ eNcellXY

Int_t EdbViewCell::eNcellXY

◆ eNcl

Int_t EdbViewCell::eNcl

◆ eNeib

Int_t EdbViewCell::eNeib[9]
private

index of the first cell of the frame in epCell array

◆ eNfr

Int_t EdbViewCell::eNfr

◆ eNx

Int_t EdbViewCell::eNx

◆ eNy

Int_t EdbViewCell::eNy

◆ epC

EdbCluster** EdbViewCell::epC
private

◆ epCell

EdbCluster*** EdbViewCell::epCell
private

pointers to clusters [eNcellsLim*eCellLim]

◆ eXmax

Float_t EdbViewCell::eXmax

◆ eXmin

Float_t EdbViewCell::eXmin

◆ eYmax

Float_t EdbViewCell::eYmax

◆ eYmin

Float_t EdbViewCell::eYmin

◆ eZmax

Float_t EdbViewCell::eZmax

◆ eZmin

Float_t EdbViewCell::eZmin

The documentation for this class was generated from the following files: