--- PamCut/TofCuts/TofDedxHeCut/TofDedxHeCut.cpp	2009/08/14 10:30:49	1.2
+++ PamCut/TofCuts/TofDedxHeCut/TofDedxHeCut.cpp	2009/09/22 10:24:26	1.3
@@ -9,30 +9,81 @@
 
 #include "TofDedxHeCut.h"
 
+const float TofDedxHeCut::_S11_SA_params[] = { 0.27, 0.045, -4.9, 4.7, 2.6, 0.38, 0.045, 1.7, 2., 2.5 };
+const float TofDedxHeCut::_S12_SA_params[] = { 0.27, 0.043, -4.9, 4.7, 2.6, 0.38, 0.042, 1.7, 2., 2.5 };
+const float TofDedxHeCut::_S21_SA_params[] = { 0.29, 0.037, -4.52, 4.4, 2.4, 0.28, 0.031, 1., 2.2, 2. };
+const float TofDedxHeCut::_S22_SA_params[] = { 0.31, 0.033, -4.52, 4.4, 2.4, 0.29, 0.02, 1., 2.4, 1.9 };
+const float TofDedxHeCut::_S31_SA_params[] = { 0.44, 0.029, -2.4, 3.2, 2.5, 0.56, 0.008, 3.2, 1.2, 2.3 };
+const float TofDedxHeCut::_S32_SA_params[] = { 0.48, 0.024, -2.4, 3.4, 2.5, 0.64, -0.0071, 3.2, 1.3, 2.3 };
+
+TofDedxHeCut::TofDedxHeCut(const char *cutName, unsigned int layer, bool standAlone) :
+  PamCut(cutName), _layer(layer), _standAlone(standAlone), _p(NULL) {
+
+  switch (_layer) {
+  case S11:
+    _p = _S11_SA_params;
+    break;
+  case S12:
+    _p = _S12_SA_params;
+    break;
+  case S21:
+    _p = _S21_SA_params;
+    break;
+  case S22:
+    _p = _S22_SA_params;
+    break;
+  case S31:
+    _p = _S31_SA_params;
+    break;
+  case S32:
+    _p = _S32_SA_params;
+    break;
+  default:
+    _p = _S21_SA_params;
+  }
+
+}
+
 int TofDedxHeCut::Check(PamLevel2 *event) {
 
-  int trkSeqNo;
-  if(_standAlone)
+  static int trkSeqNo;
+  if (_standAlone)
     trkSeqNo = 0;
-  else
-    trkSeqNo = event->GetTrack(0)->GetToFTrack()->trkseqno;
+  else {
+    for (int i = 0; i < event->GetToFLevel2()->ntrk(); i++) {
+      if (event->GetToFLevel2()->GetToFTrkVar(i)->trkseqno == event->GetTrack(0)->GetTrkTrack()->GetSeqNo())
+        trkSeqNo = i;
+    }
+  }
 
   float beta = event->GetToFLevel2()->GetToFTrkVar(trkSeqNo)->beta[12];
 
-  //float beta = event->GetTrack(0)->GetToFTrack()->beta[12];
-  //int trkSeqNo = event->GetTrack(0)->GetToFTrack()->trkseqno;
-  float dEdx = event->GetToFLevel2()->GetdEdx(trkSeqNo, 0, 100) + event->GetToFLevel2()->GetdEdx(trkSeqNo, 1, 100)
-      + event->GetToFLevel2()->GetdEdx(trkSeqNo, 2, 100) + event->GetToFLevel2()->GetdEdx(trkSeqNo, 3, 100);
-  dEdx /= 4;
+  static unsigned int auxLayer;
+  static float dEdx;
+  static unsigned int currLayer;
+  static unsigned int nLayers;
+
+  auxLayer = _layer;
+  dEdx = 0;
+  currLayer = 0;
+  nLayers = 0;
+  while (auxLayer > 0) {
+    if ((auxLayer & 0x1) == 1) {
+      //cout << auxLayer << "  " << currLayer << endl;
+      dEdx += event->GetToFLevel2()->GetdEdx(trkSeqNo, currLayer, 100);
+      nLayers++;
+    }
+    currLayer++;
+    auxLayer >>= 1;
+  }
 
-  if (dEdx < 1 / (0.2 * (beta - 0.06)) - 2.)
-    return 0;
+  dEdx /= nLayers;
 
-  if (dEdx > 1 / (0.25 * (beta - 0.12)))
-    return 1;
+  if (dEdx < 1. / (_p[0] * pow(beta, _p[4]) + _p[1]) + _p[2] + _p[3] * beta)
+    return 0;
 
-  if (dEdx < 11. * beta - 7.7)
-    return 2;
+  if (dEdx > 1. / (_p[5] * pow(beta, _p[9]) + _p[6]) + _p[7] + _p[8] * beta)
+    return 0;
 
   return CUTOK;
 }