TrkCuts/TrkNucleiZCut/TrkNucleiZCut.cpp

/*
 * TrkNucleiZCut.cpp
 *
 *  Created on: 10-mag-2009
 *      Author: Nicola Mori
 */

/*! @file TrkNucleiZCut.cpp The TrkNucleiZCut class implementation file */
#ifndef NO_TRKNUCLEI

#include "TrkNucleiZCut.h"

const float TrkNucleiZCut::_meanRigMean[] = {.992, 1.99, 3., 4., 5., 6.}; // the charge distribution peak for Z=1..6
const float TrkNucleiZCut::_sigmaRigMean[] = {0.06, 0.1, 0.17, 0.34, 0.44, 0.59}; // the charge-width for Z=1..6

const float TrkNucleiZCut::_meanBetaMean[] = {1., 2., 3., 4., 5., 6.}; // the charge distribution peak for Z=1..6
const float TrkNucleiZCut::_sigmaBetaMean[] = {0.09, 0.1, 0.17, 0.34, 0.44, 0.59}; // the charge-width for Z=1..6


int TrkNucleiZCut::Check(PamLevel2 *event) {

  if (_lowZ < 1 || _lowZ > 6) {
    return TRKNUCLEIZ_ILLEGALLOWZ;
  }
  if ((_highZ < _lowZ && _highZ != -1) || _highZ > 6) {
    return TRKNUCLEIZ_ILLEGALLOWZ;
  }

  /* Check if we must perform a range selection */
  bool range = false;
  if (_highZ 1= -1)
  range = true;

  // ******** Check if TrkNuclei has already been initialized ********
  if (_trkNuclei == NULL)
  _trkNuclei = new TrkNuclei();
  _trkNuclei->Set(event->GetTrack(0)->GetTrkTrack());

  // ******** First of all, check if the minimum number of planes is matched. ********
  // NOTE: in current implementation a plane is considered valid if it is tagged as good in Tracker level 2
  //       routines. This should correspond to a valid charge.
  // TODO Check the above statement.

  // Prepare the valid charges counter
  unsigned int nValid = 0;
  for (unsigned int ip = 0; ip < 6; ip++) {
    if (_trkNuclei->GetTrkTrack()->XGood(ip) && _trkNuclei->GetTrkTrack()->XGood(ip))
    nValid++;
  }
  if (nValid < _minLayers)
  return TRKNUCLEIZ_TOOFEWLAYERS;

  // ******** Check the mean charge value ********
  if (_howToCheck == TRKNUCLEIZ_CHECKMEAN) {

    Float_t charge = 0., meanLowZ = 0., sigmaLowZ = 0., meanHighZ = 0., sigmaHighZ = 0.;;
    switch (_method) {
      case TRKNUCLEIZ_RIG:
      charge = _trkNuclei->GetZ_Rigidity();
      meanLowZ = _meanRigMean[_lowZ-1];
      sigmaLowZ = _sigmaRigMean[_lowZ-1];
      if(range) {
        meanHighZ = _meanRigMean[_highZ-1];
        sigmaHighZ = _sigmaRigMean[_highZ-1];
      }
      break;
      case TRKNUCLEIZ_BETA:
      charge = _trkNuclei->GetZ_Beta(event->GetTrack(0)->GetToFTrack()->beta[12]);
      meanLowZ = _meanBetaMean[_Z-1];
      sigmaLowZ = _sigmaBetaMean[_Z-1];
      if(range) {
        meanHighZ = _meanBetaMean[_highZ-1];
        sigmaHighZ = _sigmaBetaMean[_highZ-1];
      }

      break;
    }

    if(!range) {
      if ((charge > meanLowZ + _upperLimit * sigmaLowZ) || (charge < meanLowZ - _lowerLimit * sigmaLowZ)) {
        return TRKNUCLEIZ_OUTOFBOUNDS;
      }
    }
    else {
      if ((charge > meanHighZ + _upperLimit * sigmaHighZ) || (charge < meanLowZ - _lowerLimit * sigmaLowZ)) {
        return TRKNUCLEIZ_OUTOFBOUNDS;
      }
    }
  }
  // ******** Check the charge value from each plane ********
  // NB. In current implementation only 2 check methods (mean and single values) are used, so an "else" statement
  //     is sufficient. Implement an "if" cascade or a "switch" if another check method is necessary.

  else {
    Float_t charge[6] = {0., 0., 0., 0., 0., 0.};
    nValid = 0;
    switch (_method) {
      case TRKNUCLEIZ_RIG:
      for (int ip = 0; ip < 6; ip++) {
        charge[ip] = _trkNuclei->GetZ_Rigidity(ip);
        if (charge[ip] > 0.) {
          nValid++;
        }
      }
      break;
      case TRKNUCLEIZ_BETA:
      for (int ip = 0; ip < 6; ip++) {
        charge[ip] = _trkNuclei->GetZ_Beta(ip, event->GetTrack(0)->GetToFTrack()->beta[12]);
        if (charge[ip] > 0.) {
          nValid++;
        }
      }
      break;
    }

    for (int ip = 0; ip < 6; ip++) {
      if (charge[ip] > 0.) // Check only good layers
      // Currently, no calibration for single layers is available. Z is used as a mean value, while for sigma
      // we use the standard values use for mean dE/dx Vs. beta.
      if(!range) {
        if ((charge[ip] > _lowZ + _upperLimit * _sigmaBetaMean[_lowZ - 1]) || (charge[ip] < _lowZ - _lowerLimit * _sigmaBetaMean[_lowZ - 1])) {
          return TRKNUCLEIZ_OUTOFBOUNDS;
        }
      }
      else {
        if ((charge[ip] > _highZ + _upperLimit * _sigmaBetaMean[_highZ - 1]) || (charge[ip] < _lowZ - _lowerLimit * _sigmaBetaMean[_lowZ - 1])) {
          return TRKNUCLEIZ_OUTOFBOUNDS;
        }
      }
    }
  }

  return CUTOK;

}

#endif /* NO_TRKNUCLEI */
1	/*
2	* TrkNucleiZCut.cpp
3	*
4	* Created on: 10-mag-2009
5	* Author: Nicola Mori
6	*/
7
8	/! @file TrkNucleiZCut.cpp The TrkNucleiZCut class implementation file /
9	#ifndef NO_TRKNUCLEI
10
11	#include "TrkNucleiZCut.h"
12
13	const float TrkNucleiZCut::_meanRigMean[] = {.992, 1.99, 3., 4., 5., 6.}; // the charge distribution peak for Z=1..6
14	const float TrkNucleiZCut::_sigmaRigMean[] = {0.06, 0.1, 0.17, 0.34, 0.44, 0.59}; // the charge-width for Z=1..6
15
16	const float TrkNucleiZCut::_meanBetaMean[] = {1., 2., 3., 4., 5., 6.}; // the charge distribution peak for Z=1..6
17	const float TrkNucleiZCut::_sigmaBetaMean[] = {0.09, 0.1, 0.17, 0.34, 0.44, 0.59}; // the charge-width for Z=1..6
18
19
20	int TrkNucleiZCut::Check(PamLevel2 *event) {
21
22	if (_lowZ < 1 \|\| _lowZ > 6) {
23	return TRKNUCLEIZ_ILLEGALLOWZ;
24	}
25	if ((_highZ < _lowZ && _highZ != -1) \|\| _highZ > 6) {
26	return TRKNUCLEIZ_ILLEGALLOWZ;
27	}
28
29	/* Check if we must perform a range selection */
30	bool range = false;
31	if (_highZ 1= -1)
32	range = true;
33
34	// ****** Check if TrkNuclei has already been initialized ******
35	if (_trkNuclei == NULL)
36	_trkNuclei = new TrkNuclei();
37	_trkNuclei->Set(event->GetTrack(0)->GetTrkTrack());
38
39	// ****** First of all, check if the minimum number of planes is matched. ******
40	// NOTE: in current implementation a plane is considered valid if it is tagged as good in Tracker level 2
41	// routines. This should correspond to a valid charge.
42	// TODO Check the above statement.
43
44	// Prepare the valid charges counter
45	unsigned int nValid = 0;
46	for (unsigned int ip = 0; ip < 6; ip++) {
47	if (_trkNuclei->GetTrkTrack()->XGood(ip) && _trkNuclei->GetTrkTrack()->XGood(ip))
48	nValid++;
49	}
50	if (nValid < _minLayers)
51	return TRKNUCLEIZ_TOOFEWLAYERS;
52
53	// ****** Check the mean charge value ******
54	if (_howToCheck == TRKNUCLEIZ_CHECKMEAN) {
55
56	Float_t charge = 0., meanLowZ = 0., sigmaLowZ = 0., meanHighZ = 0., sigmaHighZ = 0.;;
57	switch (_method) {
58	case TRKNUCLEIZ_RIG:
59	charge = _trkNuclei->GetZ_Rigidity();
60	meanLowZ = _meanRigMean[_lowZ-1];
61	sigmaLowZ = _sigmaRigMean[_lowZ-1];
62	if(range) {
63	meanHighZ = _meanRigMean[_highZ-1];
64	sigmaHighZ = _sigmaRigMean[_highZ-1];
65	}
66	break;
67	case TRKNUCLEIZ_BETA:
68	charge = _trkNuclei->GetZ_Beta(event->GetTrack(0)->GetToFTrack()->beta[12]);
69	meanLowZ = _meanBetaMean[_Z-1];
70	sigmaLowZ = _sigmaBetaMean[_Z-1];
71	if(range) {
72	meanHighZ = _meanBetaMean[_highZ-1];
73	sigmaHighZ = _sigmaBetaMean[_highZ-1];
74	}
75
76	break;
77	}
78
79	if(!range) {
80	if ((charge > meanLowZ + _upperLimit * sigmaLowZ) \|\| (charge < meanLowZ - _lowerLimit * sigmaLowZ)) {
81	return TRKNUCLEIZ_OUTOFBOUNDS;
82	}
83	}
84	else {
85	if ((charge > meanHighZ + _upperLimit * sigmaHighZ) \|\| (charge < meanLowZ - _lowerLimit * sigmaLowZ)) {
86	return TRKNUCLEIZ_OUTOFBOUNDS;
87	}
88	}
89	}
90	// ****** Check the charge value from each plane ******
91	// NB. In current implementation only 2 check methods (mean and single values) are used, so an "else" statement
92	// is sufficient. Implement an "if" cascade or a "switch" if another check method is necessary.
93
94	else {
95	Float_t charge[6] = {0., 0., 0., 0., 0., 0.};
96	nValid = 0;
97	switch (_method) {
98	case TRKNUCLEIZ_RIG:
99	for (int ip = 0; ip < 6; ip++) {
100	charge[ip] = _trkNuclei->GetZ_Rigidity(ip);
101	if (charge[ip] > 0.) {
102	nValid++;
103	}
104	}
105	break;
106	case TRKNUCLEIZ_BETA:
107	for (int ip = 0; ip < 6; ip++) {
108	charge[ip] = _trkNuclei->GetZ_Beta(ip, event->GetTrack(0)->GetToFTrack()->beta[12]);
109	if (charge[ip] > 0.) {
110	nValid++;
111	}
112	}
113	break;
114	}
115
116	for (int ip = 0; ip < 6; ip++) {
117	if (charge[ip] > 0.) // Check only good layers
118	// Currently, no calibration for single layers is available. Z is used as a mean value, while for sigma
119	// we use the standard values use for mean dE/dx Vs. beta.
120	if(!range) {
121	if ((charge[ip] > _lowZ + _upperLimit * _sigmaBetaMean[_lowZ - 1]) \|\| (charge[ip] < _lowZ - _lowerLimit * _sigmaBetaMean[_lowZ - 1])) {
122	return TRKNUCLEIZ_OUTOFBOUNDS;
123	}
124	}
125	else {
126	if ((charge[ip] > _highZ + _upperLimit * _sigmaBetaMean[_highZ - 1]) \|\| (charge[ip] < _lowZ - _lowerLimit * _sigmaBetaMean[_lowZ - 1])) {
127	return TRKNUCLEIZ_OUTOFBOUNDS;
128	}
129	}
130	}
131	}
132
133	return CUTOK;
134
135	}
136
137	#endif /* NO_TRKNUCLEI */