Index: modules/CscClusterId.cc
===================================================================
--- modules/CscClusterId.cc	(revision 157fe1366e381630b70ab906a5dd4234b8a8f487)
+++ modules/CscClusterId.cc	(revision 4aec3830c8c09116e39d3a66ee5ef9d06db4e0b4)
@@ -21,4 +21,5 @@
   *  This module is specific to the CMS paper searching for neutral LLPs in the CMS endcap muon detectors: https://arxiv.org/abs/2107.04838
   *  It is implemented based on the cut_based_id.py function provided in the HEPData entry of the paper: https://www.hepdata.net/record/104408
+  *  to reproduce the cut-based ID efficiency of the CMS paper.
   *
   *  \author Christina Wang
@@ -48,5 +49,5 @@
 #include <sstream>
 #include <stdexcept>
-
+#include "assert.h"
 using namespace std;
 
@@ -54,7 +55,8 @@
 
 CscClusterId::CscClusterId() :
-  fFormula(0), fItInputArray(0)
+  fFormula(0), fEtaFormula(0), fItInputArray(0)
 {
   fFormula = new DelphesCscClusterFormula;
+  fEtaFormula = new DelphesCscClusterFormula;
 }
 
@@ -64,4 +66,5 @@
 {
   if(fFormula) delete fFormula;
+  if(fEtaFormula) delete fEtaFormula;
 }
 
@@ -73,4 +76,6 @@
 
   fFormula->Compile(GetString("EfficiencyFormula", "1.0"));
+  fEtaFormula->Compile(GetString("EtaCutFormula", "1.0"));
+  fEtaCutMax = GetDouble("EtaCutMax", 999.0);
 
   // import input array
@@ -114,19 +119,10 @@
     eta = (cosTheta == 1.0 ? signPz * 999.9 : momentum.Eta());
 
+    // calculate the NStation > 1 efficiency, implemented according to Additional Figure 8 in HEPData
     NStationEff = fFormula->Eval(decayR, decayZ, Ehad);
 
-    // assign average station for the cluster
-    if (decayZ < 6320) avgStation = 1;
-    else if (decayZ < 7240 && decayR > 2750)avgStation = 1;
-    else if (decayZ < 8500) avgStation = 2;
-    else if (decayZ < 9700) avgStation = 3;
-    else avgStation = 4;
-
-    // if NStation == 1, different eta cut is applied
-    if (avgStation == 1) eta_cut = 1.8;
-    else if (avgStation == 2) eta_cut = 1.6;
-    else if (avgStation == 3) eta_cut = 1.6;
-    else if (avgStation == 4) eta_cut = 1.8;
-    if(gRandom->Uniform() > NStationEff*(abs(eta)<1.9)+(1.0-NStationEff)*(abs(eta)<eta_cut)) continue;
+    // depending on the decay region (station Number), different eta cut is applied, implemented based on cut_based_id.py in HEPData
+    float eta_cut = fEtaFormula->Eval(decayR, decayZ);
+    if(gRandom->Uniform() > NStationEff*(abs(eta)<fEtaCutMax)+(1.0-NStationEff)*(abs(eta)<eta_cut)) continue;
 
     fOutputArray->Add(candidate);
Index: modules/CscClusterId.h
===================================================================
--- modules/CscClusterId.h	(revision 157fe1366e381630b70ab906a5dd4234b8a8f487)
+++ modules/CscClusterId.h	(revision 4aec3830c8c09116e39d3a66ee5ef9d06db4e0b4)
@@ -46,4 +46,6 @@
 private:
   DelphesCscClusterFormula *fFormula; //!
+  DelphesCscClusterFormula *fEtaFormula; //!
+  Double_t fEtaCutMax;
 
   TIterator *fItInputArray; //!
Index: modules/LLPFilter.cc
===================================================================
--- modules/LLPFilter.cc	(revision 157fe1366e381630b70ab906a5dd4234b8a8f487)
+++ modules/LLPFilter.cc	(revision 4aec3830c8c09116e39d3a66ee5ef9d06db4e0b4)
@@ -81,10 +81,10 @@
   fRequireDecayRegion = GetBool("RequireDecayRegion", 0);
 
-  fDecayRegionRMax = GetInt("DecayRegionRMax", 0); //mm
-  fDecayRegionRMin = GetInt("DecayRegionRMin", 0); //mm
-  fDecayRegionZMax = GetInt("DecayRegionZMax", 0); //mm
-  fDecayRegionZMin = GetInt("DecayRegionZMin", 0); //mm
-  fDecayRegionEtaMax = GetInt("DecayRegionEtaMax", 0); // requirement on abs(eta)
-  fDecayRegionEtaMin = GetInt("DecayRegionEtaMin", 0); //requirement on abs(eta)
+  fDecayRegionRMax = GetDouble("DecayRegionRMax", 0.0); //mm
+  fDecayRegionRMin = GetDouble("DecayRegionRMin", 0.0); //mm
+  fDecayRegionZMax = GetDouble("DecayRegionZMax", 0.0); //mm
+  fDecayRegionZMin = GetDouble("DecayRegionZMin", 0.0); //mm
+  fDecayRegionEtaMax = GetDouble("DecayRegionEtaMax", 0.0); // requirement on abs(eta)
+  fDecayRegionEtaMin = GetDouble("DecayRegionEtaMin", 0.0); //requirement on abs(eta)
 
 
Index: modules/LLPFilter.h
===================================================================
--- modules/LLPFilter.h	(revision 157fe1366e381630b70ab906a5dd4234b8a8f487)
+++ modules/LLPFilter.h	(revision 4aec3830c8c09116e39d3a66ee5ef9d06db4e0b4)
@@ -52,10 +52,10 @@
   Double_t fPTMin; //!
   Bool_t fRequireDecayRegion;
-  Int_t fDecayRegionRMax;
-  Int_t fDecayRegionRMin;
-  Int_t fDecayRegionZMax;
-  Int_t fDecayRegionZMin;
-  Int_t fDecayRegionEtaMax;
-  Int_t fDecayRegionEtaMin;
+  Double_t fDecayRegionRMax;
+  Double_t fDecayRegionRMin;
+  Double_t fDecayRegionZMax;
+  Double_t fDecayRegionZMin;
+  Double_t fDecayRegionEtaMax;
+  Double_t fDecayRegionEtaMin;
   Int_t fDaughterNumber;
   Bool_t fInvert; //!
