diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/CMakeLists.txt b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/CMakeLists.txt
index b17bb2488639e2ae39cc482550dbd234eca4afd0..bcd6cbe7aa18dd8e0fd7c3e20a1881e2de1d9397 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/CMakeLists.txt
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/CMakeLists.txt
@@ -36,7 +36,10 @@ atlas_add_root_dictionary( ISF_FastCaloSimEvent _dictSource
                            ISF_FastCaloSimEvent/TFCS2DFunctionHistogram.h 
                            ISF_FastCaloSimEvent/TFCSParametrizationBase.h 
                            ISF_FastCaloSimEvent/TFCSParametrization.h 
+                           ISF_FastCaloSimEvent/TFCSInvisibleParametrization.h
                            ISF_FastCaloSimEvent/TFCSInitWithEkin.h
+                           ISF_FastCaloSimEvent/TFCSEnergyInterpolationLinear.h
+                           ISF_FastCaloSimEvent/TFCSEnergyInterpolationSpline.h
                            ISF_FastCaloSimEvent/TFCSParametrizationChain.h 
                            ISF_FastCaloSimEvent/TFCSParametrizationBinnedChain.h
                            ISF_FastCaloSimEvent/TFCSParametrizationFloatSelectChain.h
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/LinkDef.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/LinkDef.h
index 828419585279297148a719dbe1f1be6e8eb14373..3bbf13ed340d9e995a1c626b5ab3c8bab74a1b31 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/LinkDef.h
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/LinkDef.h
@@ -14,7 +14,10 @@
 
 #include "ISF_FastCaloSimEvent/TFCSParametrizationBase.h"
 #include "ISF_FastCaloSimEvent/TFCSParametrization.h"
+#include "ISF_FastCaloSimEvent/TFCSInvisibleParametrization.h"
 #include "ISF_FastCaloSimEvent/TFCSInitWithEkin.h"
+#include "ISF_FastCaloSimEvent/TFCSEnergyInterpolationLinear.h"
+#include "ISF_FastCaloSimEvent/TFCSEnergyInterpolationSpline.h"
 #include "ISF_FastCaloSimEvent/TFCSParametrizationChain.h"
 #include "ISF_FastCaloSimEvent/TFCSParametrizationBinnedChain.h"
 #include "ISF_FastCaloSimEvent/TFCSParametrizationFloatSelectChain.h"
@@ -50,10 +53,12 @@
 #pragma link C++ class TFCS1DFunctionRegressionTF+;
 #pragma link C++ class TFCS2DFunction+;
 #pragma link C++ class TFCS2DFunctionHistogram+;
-
 #pragma link C++ class TFCSParametrizationBase+;
 #pragma link C++ class TFCSParametrization+;
+#pragma link C++ class TFCSInvisibleParametrization+;
 #pragma link C++ class TFCSInitWithEkin+;
+#pragma link C++ class TFCSEnergyInterpolationLinear+;
+#pragma link C++ class TFCSEnergyInterpolationSpline+;
 #pragma link C++ class TFCSParametrizationChain+;
 #pragma link C++ class TFCSParametrizationBinnedChain+;
 #pragma link C++ class TFCSParametrizationFloatSelectChain+;
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSEnergyInterpolationLinear.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSEnergyInterpolationLinear.h
new file mode 100644
index 0000000000000000000000000000000000000000..863f1fef0dc8c2700696c143c751f1c8aeb4822c
--- /dev/null
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSEnergyInterpolationLinear.h
@@ -0,0 +1,37 @@
+/*
+  Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef ISF_FASTCALOSIMEVENT_TFCSEnergyInterpolationLinear_h
+#define ISF_FASTCALOSIMEVENT_TFCSEnergyInterpolationLinear_h
+
+#include "ISF_FastCaloSimEvent/TFCSParametrization.h"
+
+class TFCSEnergyInterpolationLinear:public TFCSParametrization {
+public:
+  TFCSEnergyInterpolationLinear(const char* name=nullptr, const char* title=nullptr);
+
+  virtual bool is_match_Ekin_bin(int /*Ekin_bin*/) const {return true;};
+  virtual bool is_match_calosample(int /*calosample*/) const {return true;};
+  
+  void set_slope(float slope) {m_slope=slope;};
+  void set_offset(float offset) {m_offset=offset;};
+
+  // Initialize simulstate with the mean reconstructed energy in the calorimater expeted from the true kinetic energy
+  virtual void simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol);
+
+  void Print(Option_t *option="") const;
+
+  static void unit_test(TFCSSimulationState* simulstate=nullptr,TFCSTruthState* truth=nullptr, const TFCSExtrapolationState* extrapol=nullptr);
+private:
+  float m_slope;
+  float m_offset;
+
+  ClassDef(TFCSEnergyInterpolationLinear,1)  //TFCSEnergyInterpolationLinear
+};
+
+#if defined(__ROOTCLING__) && defined(__FastCaloSimStandAlone__)
+#pragma link C++ class TFCSEnergyInterpolationLinear+;
+#endif
+
+#endif
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSEnergyInterpolationSpline.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSEnergyInterpolationSpline.h
new file mode 100644
index 0000000000000000000000000000000000000000..568be359c0c047d3157b6be5ae5ce02f182c7bac
--- /dev/null
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSEnergyInterpolationSpline.h
@@ -0,0 +1,50 @@
+/*
+  Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef ISF_FASTCALOSIMEVENT_TFCSEnergyInterpolationSpline_h
+#define ISF_FASTCALOSIMEVENT_TFCSEnergyInterpolationSpline_h
+
+#include "ISF_FastCaloSimEvent/TFCSParametrization.h"
+#include "TSpline.h"
+
+class TGraph;
+
+class TFCSEnergyInterpolationSpline:public TFCSParametrization {
+public:
+  TFCSEnergyInterpolationSpline(const char* name=nullptr, const char* title=nullptr);
+
+  virtual bool is_match_Ekin_bin(int /*Ekin_bin*/) const {return true;};
+  virtual bool is_match_calosample(int /*calosample*/) const {return true;};
+  
+  ///Initialize interpolation from spline
+  ///x values should be log(Ekin), y values should <E(reco)/Ekin(true)>
+  void InitFromSpline(const TSpline3& spline) {m_spline=spline;};
+
+  ///Initialize spline interpolation from arrays in log(Ekin) and response=<E(reco)/Ekin(true)>
+  ///opt, valbeg and valend as defined for TSpline3
+  void InitFromArrayInLogEkin(Int_t np, Double_t logEkin[], Double_t response[], const char *opt=nullptr,Double_t valbeg=0, Double_t valend=0);
+
+  ///Initialize spline interpolation from arrays in Ekin and response=<E(reco)/Ekin(true)>
+  ///opt, valbeg and valend as defined for TSpline3
+  void InitFromArrayInEkin(Int_t np, Double_t Ekin[], Double_t response[], const char *opt=nullptr,Double_t valbeg=0, Double_t valend=0);
+
+  const TSpline3& spline() const {return m_spline;};
+
+  ///Initialize simulstate with the mean reconstructed energy in the calorimater expeted from the true kinetic energy
+  virtual void simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol);
+
+  void Print(Option_t *option="") const;
+
+  static void unit_test(TFCSSimulationState* simulstate=nullptr,TFCSTruthState* truth=nullptr, const TFCSExtrapolationState* extrapol=nullptr);
+private:
+  TSpline3 m_spline;
+
+  ClassDef(TFCSEnergyInterpolationSpline,1)  //TFCSEnergyInterpolationSpline
+};
+
+#if defined(__ROOTCLING__) && defined(__FastCaloSimStandAlone__)
+#pragma link C++ class TFCSEnergyInterpolationSpline+;
+#endif
+
+#endif
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSInitWithEkin.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSInitWithEkin.h
index 5e73b4abb3546e2a7e9c6b136bf5f87f14c3dfc0..b782795dbd60cb62021625a0a146bf6143a01439 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSInitWithEkin.h
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSInitWithEkin.h
@@ -5,20 +5,14 @@
 #ifndef ISF_FASTCALOSIMEVENT_TFCSInitWithEkin_h
 #define ISF_FASTCALOSIMEVENT_TFCSInitWithEkin_h
 
-#include "ISF_FastCaloSimEvent/TFCSParametrizationBase.h"
+#include "ISF_FastCaloSimEvent/TFCSParametrization.h"
 
-class TFCSInitWithEkin:public TFCSParametrizationBase {
+class TFCSInitWithEkin:public TFCSParametrization {
 public:
   TFCSInitWithEkin(const char* name=nullptr, const char* title=nullptr);
 
-  virtual bool is_match_Ekin(float /*Ekin*/) const {return true;};
-  virtual bool is_match_eta(float /*eta*/) const {return true;};
-
   virtual bool is_match_Ekin_bin(int /*Ekin_bin*/) const {return true;};
   virtual bool is_match_calosample(int /*calosample*/) const {return true;};
-
-  virtual bool is_match_all_Ekin() const {return true;};
-  virtual bool is_match_all_eta() const {return true;};
   virtual bool is_match_all_Ekin_bin() const {return true;};
   virtual bool is_match_all_calosample() const {return true;};
 
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSInvisibleParametrization.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSInvisibleParametrization.h
new file mode 100644
index 0000000000000000000000000000000000000000..cac15460c1802a7a9c9a5dbc6aa282115c44a6ef
--- /dev/null
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSInvisibleParametrization.h
@@ -0,0 +1,27 @@
+/*
+  Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef ISF_FASTCALOSIMEVENT_TFCSInvisibleParametrization_h
+#define ISF_FASTCALOSIMEVENT_TFCSInvisibleParametrization_h
+
+#include "ISF_FastCaloSimEvent/TFCSParametrization.h"
+
+class TFCSInvisibleParametrization:public TFCSParametrization {
+public:
+  TFCSInvisibleParametrization(const char* name=nullptr, const char* title=nullptr):TFCSParametrization(name,title) {};
+
+  virtual bool is_match_Ekin_bin(int /*Ekin_bin*/) const {return true;};
+  virtual bool is_match_calosample(int /*calosample*/) const {return true;};
+
+  void simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol);
+private:
+
+  ClassDef(TFCSInvisibleParametrization,1)  //TFCSInvisibleParametrization
+};
+
+#if defined(__ROOTCLING__) && defined(__FastCaloSimStandAlone__)
+#pragma link C++ class TFCSInvisibleParametrization+;
+#endif
+
+#endif
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSParametrizationEkinSelectChain.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSParametrizationEkinSelectChain.h
index 74a8a809ca1e6d677c862e3e7d2a18a738ade86e..04cefb369498113173a7c0e645c8a4d551992970 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSParametrizationEkinSelectChain.h
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSParametrizationEkinSelectChain.h
@@ -10,8 +10,17 @@
 
 class TFCSParametrizationEkinSelectChain:public TFCSParametrizationFloatSelectChain {
 public:
-  TFCSParametrizationEkinSelectChain(const char* name=nullptr, const char* title=nullptr):TFCSParametrizationFloatSelectChain(name,title) {};
-  TFCSParametrizationEkinSelectChain(const TFCSParametrizationEkinSelectChain& ref):TFCSParametrizationFloatSelectChain(ref) {};
+  TFCSParametrizationEkinSelectChain(const char* name=nullptr, const char* title=nullptr):TFCSParametrizationFloatSelectChain(name,title) {reset_DoRandomInterpolation();};
+  TFCSParametrizationEkinSelectChain(const TFCSParametrizationEkinSelectChain& ref):TFCSParametrizationFloatSelectChain(ref) {reset_DoRandomInterpolation();};
+
+  ///Status bit for Ekin Selection
+  enum FCSEkinStatusBits {
+     kDoRandomInterpolation = BIT(15) ///< Set this bit in the TObject bit field if a random selection between neighbouring Ekin bins should be done
+  };
+
+  bool DoRandomInterpolation() const {return TestBit(kDoRandomInterpolation);};
+  void set_DoRandomInterpolation() {SetBit(kDoRandomInterpolation);};
+  void reset_DoRandomInterpolation() {ResetBit(kDoRandomInterpolation);};
 
   using TFCSParametrizationFloatSelectChain::push_back_in_bin;
   virtual void push_back_in_bin(TFCSParametrizationBase* param);
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSParametrizationPDGIDSelectChain.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSParametrizationPDGIDSelectChain.h
index 88492fe86b60781bd91a5785aa18b8bcb52c99ba..56d8d978ef656925376f7fb6e21cdad5241819d0 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSParametrizationPDGIDSelectChain.h
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSParametrizationPDGIDSelectChain.h
@@ -9,8 +9,17 @@
 
 class TFCSParametrizationPDGIDSelectChain:public TFCSParametrizationChain {
 public:
-  TFCSParametrizationPDGIDSelectChain(const char* name=nullptr, const char* title=nullptr):TFCSParametrizationChain(name,title) {};
-  TFCSParametrizationPDGIDSelectChain(const TFCSParametrizationPDGIDSelectChain& ref):TFCSParametrizationChain(ref) {};
+  TFCSParametrizationPDGIDSelectChain(const char* name=nullptr, const char* title=nullptr):TFCSParametrizationChain(name,title) {reset_SimulateOnlyOnePDGID();};
+  TFCSParametrizationPDGIDSelectChain(const TFCSParametrizationPDGIDSelectChain& ref):TFCSParametrizationChain(ref) {reset_SimulateOnlyOnePDGID();};
+
+  ///Status bit for PDGID Selection
+  enum FCSPDGIDStatusBits {
+     kSimulateOnlyOnePDGID = BIT(15) ///< Set this bit in the TObject bit field if the PDGID selection loop should be aborted after the first successful match
+  };
+
+  bool SimulateOnlyOnePDGID() const {return TestBit(kSimulateOnlyOnePDGID);};
+  void set_SimulateOnlyOnePDGID() {SetBit(kSimulateOnlyOnePDGID);};
+  void reset_SimulateOnlyOnePDGID() {ResetBit(kSimulateOnlyOnePDGID);};
 
   virtual void simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol);
 
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCS2DFunctionHistogram.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCS2DFunctionHistogram.cxx
index 9b10dbbf738cd8647b29431826cd1a1a7a2b0429..de1726138cef6fae1f9b765c134d0bb47e11a205 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCS2DFunctionHistogram.cxx
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCS2DFunctionHistogram.cxx
@@ -31,7 +31,10 @@ void TFCS2DFunctionHistogram::Initialize(TH2* hist)
       float binval=hist->GetBinContent(ix,iy);
       if(binval<0) {
         //Can't work if a bin is negative, forcing bins to 0 in this case
-        std::cout<<"ERROR: bin content is negative in histogram "<<hist->GetName()<<" : "<<hist->GetTitle()<<std::endl;
+        double fraction=binval/hist->Integral();
+        if(TMath::Abs(fraction)>1e-5) {
+          std::cout<<"WARNING: bin content is negative in histogram "<<hist->GetName()<<" : "<<hist->GetTitle()<<" binval="<<binval<<" "<<fraction*100<<"% of integral="<<hist->Integral()<<". Forcing bin to 0."<<std::endl;
+        }  
         binval=0;
       }
       integral+=binval;
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSEnergyInterpolationLinear.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSEnergyInterpolationLinear.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..9744d1a9eb82df21e4502af90bb43b6e3c6404eb
--- /dev/null
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSEnergyInterpolationLinear.cxx
@@ -0,0 +1,82 @@
+/*
+  Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "ISF_FastCaloSimEvent/TFCSEnergyInterpolationLinear.h"
+#include "ISF_FastCaloSimEvent/TFCSSimulationState.h"
+#include "ISF_FastCaloSimEvent/TFCSTruthState.h"
+#include "ISF_FastCaloSimEvent/TFCSExtrapolationState.h"
+#include "TCanvas.h"
+#include "TGraph.h"
+#include "TAxis.h"
+#include <iostream>
+
+//=============================================
+//======= TFCSEnergyInterpolationLinear =========
+//=============================================
+
+TFCSEnergyInterpolationLinear::TFCSEnergyInterpolationLinear(const char* name, const char* title):TFCSParametrization(name,title),m_slope(1),m_offset(0)
+{
+}
+
+void TFCSEnergyInterpolationLinear::simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState*)
+{
+  float Emean=m_slope*truth->Ekin()+m_offset;
+
+  ATH_MSG_DEBUG("set E="<<Emean<<" for true Ekin="<<truth->Ekin());
+  simulstate.set_E(Emean);
+}
+
+void TFCSEnergyInterpolationLinear::Print(Option_t *option) const
+{
+  TString opt(option);
+  bool shortprint=opt.Index("short")>=0;
+  bool longprint=msgLvl(MSG::DEBUG) || (msgLvl(MSG::INFO) && !shortprint);
+  TString optprint=opt;optprint.ReplaceAll("short","");
+  TFCSParametrization::Print(option);
+
+  if(longprint) ATH_MSG_INFO(optprint <<"  Emean="<<m_slope<<"*Ekin(true) + "<<m_offset);
+}
+
+void TFCSEnergyInterpolationLinear::unit_test(TFCSSimulationState* simulstate,TFCSTruthState* truth, const TFCSExtrapolationState* extrapol)
+{
+  if(!simulstate) simulstate=new TFCSSimulationState();
+  if(!truth) truth=new TFCSTruthState();
+  if(!extrapol) extrapol=new TFCSExtrapolationState();
+  
+  TFCSEnergyInterpolationLinear test("testTFCSEnergyInterpolation","test linear TFCSEnergyInterpolation");
+  test.set_pdgid(22);
+  test.set_Ekin_nominal(1000);
+  test.set_Ekin_min(1000);
+  test.set_Ekin_max(100000);
+  test.set_eta_nominal(0.225);
+  test.set_eta_min(0.2);
+  test.set_eta_max(0.25);
+  test.set_slope(0.95);
+  test.set_offset(-50);
+  test.Print();
+  
+  truth->set_pdgid(22);
+  
+  TGraph* gr=new TGraph();
+  gr->SetNameTitle("testTFCSEnergyInterpolation","test linear TFCSEnergyInterpolation");
+  gr->GetXaxis()->SetTitle("Ekin [MeV]");
+  gr->GetYaxis()->SetTitle("<E(reco)>/Ekin(true)");
+  
+  int ip=0;
+  for(float Ekin=1000;Ekin<=100000;Ekin*=1.1) {
+    //Init LorentzVector for truth. For photon Ekin=E
+    truth->SetPxPyPzE(Ekin,0,0,Ekin);
+    test.simulate(*simulstate,truth,extrapol);
+    gr->SetPoint(ip,Ekin,simulstate->E()/Ekin);
+    ++ip;
+  }  
+
+  //Drawing doesn't make sense inside athena and necessary libraries not linked by default
+  #if defined(__FastCaloSimStandAlone__)
+  TCanvas* c=new TCanvas("testTFCSEnergyInterpolation","test linear TFCSEnergyInterpolation");
+  gr->Draw("APL");
+  c->SetLogx();
+  #endif
+}
+
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSEnergyInterpolationSpline.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSEnergyInterpolationSpline.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..fbab697946ca43d94748f9e1d4829edac9d70c52
--- /dev/null
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSEnergyInterpolationSpline.cxx
@@ -0,0 +1,139 @@
+/*
+  Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "ISF_FastCaloSimEvent/TFCSEnergyInterpolationSpline.h"
+#include "ISF_FastCaloSimEvent/TFCSSimulationState.h"
+#include "ISF_FastCaloSimEvent/TFCSTruthState.h"
+#include "ISF_FastCaloSimEvent/TFCSExtrapolationState.h"
+#include "TFile.h"
+#include "TCanvas.h"
+#include "TGraph.h"
+#include "TAxis.h"
+#include <iostream>
+#include <vector>
+
+//=============================================
+//======= TFCSEnergyInterpolation =========
+//=============================================
+
+TFCSEnergyInterpolationSpline::TFCSEnergyInterpolationSpline(const char* name, const char* title):TFCSParametrization(name,title)
+{
+}
+
+void TFCSEnergyInterpolationSpline::InitFromArrayInLogEkin(Int_t np, Double_t logEkin[], Double_t response[], const char *opt,Double_t valbeg, Double_t valend)
+{
+  TSpline3 initspline(GetName(),logEkin,response,np,opt,valbeg,valend);
+  m_spline=initspline;
+}
+
+void TFCSEnergyInterpolationSpline::InitFromArrayInEkin(Int_t np, Double_t Ekin[], Double_t response[], const char *opt,Double_t valbeg, Double_t valend)
+{
+  std::vector<Double_t> logEkin(np);
+  for(int i=0;i<np;++i) logEkin[i]=TMath::Log(Ekin[i]);
+  InitFromArrayInLogEkin(np,logEkin.data(),response,opt,valbeg,valend);
+}
+
+void TFCSEnergyInterpolationSpline::simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState*)
+{
+  float Emean;
+  float logEkin=TMath::Log(truth->Ekin());
+  if(logEkin<m_spline.GetXmin()) {
+    Emean=m_spline.Eval(m_spline.GetXmin())*truth->Ekin();
+  } else {
+    Emean=m_spline.Eval(logEkin)*truth->Ekin();
+  }  
+
+  ATH_MSG_DEBUG("set E="<<Emean<<" for true Ekin="<<truth->Ekin());
+  simulstate.set_E(Emean);
+}
+
+void TFCSEnergyInterpolationSpline::Print(Option_t *option) const
+{
+  TString opt(option);
+  bool shortprint=opt.Index("short")>=0;
+  bool longprint=msgLvl(MSG::DEBUG) || (msgLvl(MSG::INFO) && !shortprint);
+  TString optprint=opt;optprint.ReplaceAll("short","");
+  TFCSParametrization::Print(option);
+
+  if(longprint) ATH_MSG_INFO(optprint <<"  Spline N="<<m_spline.GetNp()
+                           <<" "<<m_spline.GetXmin()<<"<=log(Ekin)<="<<m_spline.GetXmax()
+                           <<" "<<TMath::Exp(m_spline.GetXmin())<<"<=Ekin<="<<TMath::Exp(m_spline.GetXmax()));
+}
+
+void TFCSEnergyInterpolationSpline::unit_test(TFCSSimulationState* simulstate,TFCSTruthState* truth, const TFCSExtrapolationState* extrapol)
+{
+  if(!simulstate) simulstate=new TFCSSimulationState();
+  if(!truth) truth=new TFCSTruthState();
+  if(!extrapol) extrapol=new TFCSExtrapolationState();
+
+  const int Graph0_n=9;
+  Double_t Graph0_fx1001[Graph0_n] = {
+  1.024,
+  2.048,
+  4.094,
+  8.192,
+  16.384,
+  32.768,
+  65.536,
+  131.072,
+  262.144};
+  for(int i=0;i<Graph0_n;++i) Graph0_fx1001[i]*=1000;
+
+  Double_t Graph0_fy1001[Graph0_n] = {
+  0.6535402,
+  0.6571529,
+  0.6843001,
+  0.7172835,
+  0.7708416,
+  0.798819,
+  0.8187628,
+  0.8332745,
+  0.8443931};
+  TGraph* grspline = new TGraph(Graph0_n,Graph0_fx1001,Graph0_fy1001);
+  
+  /*
+  TFile* file=TFile::Open("Example.root");
+  TGraph* grspline=(TGraph*)file->Get("Graph");
+  file->Close();
+  */
+  TGraph* grdraw=(TGraph*)grspline->Clone();
+  grdraw->SetMarkerColor(46);
+  grdraw->SetMarkerStyle(8);
+  
+  TFCSEnergyInterpolationSpline test("testTFCSEnergyInterpolationSpline","test TFCSEnergyInterpolationSpline");
+  test.set_pdgid(22);
+  test.set_Ekin_nominal(0.5*(grdraw->GetX()[0]+grdraw->GetX()[grdraw->GetN()-1]));
+  test.set_Ekin_min(grdraw->GetX()[0]);
+  test.set_Ekin_max(grdraw->GetX()[grdraw->GetN()-1]);
+  test.set_eta_nominal(0.225);
+  test.set_eta_min(0.2);
+  test.set_eta_max(0.25);
+  test.InitFromArrayInEkin(Graph0_n,grspline->GetX(),grspline->GetY(),"b2e2",0,0);
+  test.Print();
+  
+  truth->set_pdgid(22);
+  
+  TGraph* gr=new TGraph();
+  gr->SetNameTitle("testTFCSEnergyInterpolationSplineLogX","test TFCSEnergyInterpolationSpline log x-axis");
+  gr->GetXaxis()->SetTitle("Ekin [MeV]");
+  gr->GetYaxis()->SetTitle("<E(reco)>/Ekin(true)");
+  
+  int ip=0;
+  for(float Ekin=test.Ekin_min()*0.25;Ekin<=test.Ekin_max()*4;Ekin*=1.05) {
+    //Init LorentzVector for truth. For photon Ekin=E
+    truth->SetPxPyPzE(Ekin,0,0,Ekin);
+    test.simulate(*simulstate,truth,extrapol);
+    gr->SetPoint(ip,Ekin,simulstate->E()/Ekin);
+    ++ip;
+  }  
+
+  //Drawing doesn't make sense inside athena and necessary libraries not linked by default
+  #if defined(__FastCaloSimStandAlone__)
+  TCanvas* c=new TCanvas(gr->GetName(),gr->GetTitle());
+  gr->Draw("APL");
+  grdraw->Draw("Psame");
+  c->SetLogx();
+  #endif
+}
+
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSInitWithEkin.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSInitWithEkin.cxx
index 9145226cbbcf605e131898759bcc710e921c5fe6..5eeed2cd5cb05c30f0d38083bc7875921dfff2d6 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSInitWithEkin.cxx
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSInitWithEkin.cxx
@@ -10,7 +10,7 @@
 //======= TFCSInitWithEkin =========
 //=============================================
 
-TFCSInitWithEkin::TFCSInitWithEkin(const char* name, const char* title):TFCSParametrizationBase(name,title)
+TFCSInitWithEkin::TFCSInitWithEkin(const char* name, const char* title):TFCSParametrization(name,title)
 {
   set_match_all_pdgid();
 }
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSInvisibleParametrization.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSInvisibleParametrization.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..2e0bd47d8ce829c0540dafdaaca2e64d9b745a6d
--- /dev/null
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSInvisibleParametrization.cxx
@@ -0,0 +1,15 @@
+/*
+  Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "ISF_FastCaloSimEvent/TFCSInvisibleParametrization.h"
+
+//=============================================
+//======= TFCSInvisibleParametrization =========
+//=============================================
+
+void TFCSInvisibleParametrization::simulate(TFCSSimulationState& /*simulstate*/,const TFCSTruthState* /*truth*/, const TFCSExtrapolationState* /*extrapol*/)
+{
+  ATH_MSG_VERBOSE("now in TFCSInvisibleParametrization::simulate(). Don't do anything for invisible");
+}
+
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx
index de3d05ce6d052932224d064cebd324d2448d2cd6..aeb539de5f3a51492495d0896126f0ee46f87b3f 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx
@@ -83,7 +83,9 @@ void TFCSLateralShapeParametrizationHitChain::Print(Option_t *option) const
     m_number_of_hits_simul->Print(opt+"#:");
   }
   if(longprint) ATH_MSG_INFO(optprint <<"- Simulation chain:");
+  char count='A';
   for(TFCSLateralShapeParametrizationHitBase* hitsim : m_chain) {
-    hitsim->Print(opt+"- ");
+    hitsim->Print(opt+count+" ");
+    count++;
   } 
 }
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationBinnedChain.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationBinnedChain.cxx
index 44964fcfc2e0cad559e8c7dca8f76a69038addd6..5bbcafad86342b1a9e84dbfe0b68b9bbcfffca38 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationBinnedChain.cxx
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationBinnedChain.cxx
@@ -92,18 +92,19 @@ void TFCSParametrizationBinnedChain::Print(Option_t *option) const
 
   TString prefix="- ";
   for(unsigned int ichain=0;ichain<size();++ichain) {
-    if(ichain==0) {
-      prefix="- ";
+    if(ichain==0 && ichain!=m_bin_start.front()) {
+      prefix="> ";
       if(longprint) ATH_MSG_INFO(optprint<<prefix<<"Run for all bins");
     }  
     for(unsigned int ibin=0;ibin<get_number_of_bins();++ibin) {
       if(ichain==m_bin_start[ibin]) {
+        if(ibin<get_number_of_bins()-1) if(ichain==m_bin_start[ibin+1]) continue;
         prefix=Form("%-2d",ibin);
         if(longprint) ATH_MSG_INFO(optprint<<prefix<<"Run for "<<get_bin_text(ibin));
       }  
     }  
     if(ichain==m_bin_start.back()) {
-      prefix="- ";
+      prefix="< ";
       if(longprint) ATH_MSG_INFO(optprint<<prefix<<"Run for all bins");
     }  
     chain()[ichain]->Print(opt+prefix);
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationChain.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationChain.cxx
index a680ddf8f78101185a5295fad766f5387e7d99b4..34d26f82a518d83580dafde06feca76d9120aad4 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationChain.cxx
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationChain.cxx
@@ -138,9 +138,10 @@ void TFCSParametrizationChain::Print(Option_t *option) const
   TString opt(option);
   //bool shortprint=opt.Index("short")>=0;
   //bool longprint=msgLvl(MSG::DEBUG) || (msgLvl(MSG::INFO) && !shortprint);
-  opt=opt+"- ";
 
+  char count='A';
   for(auto param: m_chain) {
-    param->Print(opt);
+    param->Print(opt+count+' ');
+    count++;
   }
 }
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationEkinSelectChain.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationEkinSelectChain.cxx
index 0c52b16a0f65c9a4d750c4ae4aaa8e99a7304a16..08f1758f05a9aa2e5a37d242aa4b9513926fe8cf 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationEkinSelectChain.cxx
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationEkinSelectChain.cxx
@@ -3,10 +3,13 @@
 */
 
 #include "ISF_FastCaloSimEvent/TFCSParametrizationEkinSelectChain.h"
+#include "ISF_FastCaloSimEvent/TFCSInvisibleParametrization.h"
 #include "ISF_FastCaloSimEvent/TFCSSimulationState.h"
 #include "ISF_FastCaloSimEvent/TFCSTruthState.h"
 #include "ISF_FastCaloSimEvent/TFCSExtrapolationState.h"
+
 #include <iostream>
+#include "TRandom.h"
 
 //=============================================
 //======= TFCSParametrizationEkinSelectChain =========
@@ -31,7 +34,61 @@ void TFCSParametrizationEkinSelectChain::push_back_in_bin(TFCSParametrizationBas
 
 int TFCSParametrizationEkinSelectChain::get_bin(TFCSSimulationState&,const TFCSTruthState* truth, const TFCSExtrapolationState*) const
 {
-  return val_to_bin(truth->Ekin());
+  float Ekin=truth->Ekin();
+  int bin=val_to_bin(Ekin);
+  
+  if(!DoRandomInterpolation()) return bin;
+
+  if(bin<0) return bin;
+  if(bin>=(int)get_number_of_bins()) return bin;
+
+  //if no parametrizations for this bin, return
+  if(m_bin_start[bin+1]==m_bin_start[bin]) return bin;
+
+  TFCSParametrizationBase* first_in_bin=chain()[m_bin_start[bin]];
+  if(!first_in_bin) return bin;
+  
+  if(Ekin<first_in_bin->Ekin_nominal()) {
+    if(bin==0) return bin;
+    int prevbin=bin-1;
+    //if no parametrizations for previous bin, return
+    if(m_bin_start[prevbin+1]==m_bin_start[prevbin]) return bin;
+    
+    TFCSParametrizationBase* first_in_prevbin=chain()[m_bin_start[prevbin]];
+    if(!first_in_prevbin) return bin;
+    
+    float logEkin=TMath::Log(Ekin);
+    float logEkin_nominal=TMath::Log(first_in_bin->Ekin_nominal());
+    float logEkin_previous=TMath::Log(first_in_prevbin->Ekin_nominal());
+    float numerator=logEkin-logEkin_previous;
+    float denominator=logEkin_nominal-logEkin_previous;
+    if(denominator<=0) return bin;
+
+    float rnd=gRandom->Rndm();
+    if(numerator/denominator<rnd) bin=prevbin;
+    ATH_MSG_DEBUG("logEkin="<<logEkin<<" logEkin_previous="<<logEkin_previous<<" logEkin_nominal="<<logEkin_nominal<<" (rnd="<<1-rnd<<" < p(previous)="<<(1-numerator/denominator)<<")? => orgbin="<<prevbin+1<<" selbin="<<bin);
+  } else {
+    if(bin==(int)get_number_of_bins()-1) return bin;
+    int nextbin=bin+1;
+    //if no parametrizations for previous bin, return
+    if(m_bin_start[nextbin+1]==m_bin_start[nextbin]) return bin;
+    
+    TFCSParametrizationBase* first_in_nextbin=chain()[m_bin_start[nextbin]];
+    if(!first_in_nextbin) return bin;
+    
+    float logEkin=TMath::Log(Ekin);
+    float logEkin_nominal=TMath::Log(first_in_bin->Ekin_nominal());
+    float logEkin_next=TMath::Log(first_in_nextbin->Ekin_nominal());
+    float numerator=logEkin-logEkin_nominal;
+    float denominator=logEkin_next-logEkin_nominal;
+    if(denominator<=0) return bin;
+
+    float rnd=gRandom->Rndm();
+    if(rnd<numerator/denominator) bin=nextbin;
+    ATH_MSG_DEBUG("logEkin="<<logEkin<<" logEkin_nominal="<<logEkin_nominal<<" logEkin_next="<<logEkin_next<<" (rnd="<<rnd<<" < p(next)="<<numerator/denominator<<")? => orgbin="<<nextbin-1<<" selbin="<<bin);
+  }
+
+  return bin;
 }
 
 const std::string TFCSParametrizationEkinSelectChain::get_variable_text(TFCSSimulationState&,const TFCSTruthState* truth, const TFCSExtrapolationState*) const
@@ -57,13 +114,13 @@ void TFCSParametrizationEkinSelectChain::unit_test(TFCSSimulationState* simulsta
   chain.setLevel(MSG::DEBUG);
 
   TFCSParametrization* param;
-  param=new TFCSParametrization("A begin all","A begin all");
+  param=new TFCSInvisibleParametrization("A begin all","A begin all");
   param->setLevel(MSG::DEBUG);
   param->set_Ekin_nominal(2);
   param->set_Ekin_min(2);
   param->set_Ekin_max(5);
   chain.push_before_first_bin(param);
-  param=new TFCSParametrization("A end all","A end all");
+  param=new TFCSInvisibleParametrization("A end all","A end all");
   param->setLevel(MSG::DEBUG);
   param->set_Ekin_nominal(2);
   param->set_Ekin_min(2);
@@ -72,29 +129,29 @@ void TFCSParametrizationEkinSelectChain::unit_test(TFCSSimulationState* simulsta
 
   const int n_params=5;
   for(int i=2;i<n_params;++i) {
-    param=new TFCSParametrization(Form("A%d",i),Form("A %d",i));
+    param=new TFCSInvisibleParametrization(Form("A%d",i),Form("A %d",i));
     param->setLevel(MSG::DEBUG);
-    param->set_Ekin_nominal(i*i+0.1);
-    param->set_Ekin_min(i*i);
-    param->set_Ekin_max((i+1)*(i+1));
+    param->set_Ekin_nominal(TMath::Power(2.0,i));
+    param->set_Ekin_min(TMath::Power(2.0,i-0.5));
+    param->set_Ekin_max(TMath::Power(2.0,i+0.5));
     chain.push_back_in_bin(param);
   }
   for(int i=n_params;i>=1;--i) {
-    param=new TFCSParametrization(Form("B%d",i),Form("B %d",i));
+    param=new TFCSInvisibleParametrization(Form("B%d",i),Form("B %d",i));
     param->setLevel(MSG::DEBUG);
-    param->set_Ekin_nominal(i*i+0.1);
-    param->set_Ekin_min(i*i);
-    param->set_Ekin_max((i+1)*(i+1));
+    param->set_Ekin_nominal(TMath::Power(2.0,i));
+    param->set_Ekin_min(TMath::Power(2.0,i-0.5));
+    param->set_Ekin_max(TMath::Power(2.0,i+0.5));
     chain.push_back_in_bin(param);
   }
 
   std::cout<<"====         Chain setup       ===="<<std::endl;
   chain.Print();
 
-  param=new TFCSParametrization("B end all","B end all");
+  param=new TFCSInvisibleParametrization("B end all","B end all");
   param->setLevel(MSG::DEBUG);
   chain.push_back(param);
-  param=new TFCSParametrization("B begin all","B begin all");
+  param=new TFCSInvisibleParametrization("B begin all","B begin all");
   param->setLevel(MSG::DEBUG);
   chain.push_before_first_bin(param);
   
@@ -103,19 +160,23 @@ void TFCSParametrizationEkinSelectChain::unit_test(TFCSSimulationState* simulsta
   std::cout<<"==== Simulate with E=0.3      ===="<<std::endl;
   truth->SetPtEtaPhiM(0.3,0,0,0);
   chain.simulate(*simulstate,truth,extrapol);
-  std::cout<<"==== Simulate with E=1.0      ===="<<std::endl;
-  truth->SetPtEtaPhiM(1,0,0,0);
-  chain.simulate(*simulstate,truth,extrapol);
-  std::cout<<"==== Simulate with E=1.3      ===="<<std::endl;
-  truth->SetPtEtaPhiM(1.3,0,0,0);
-  chain.simulate(*simulstate,truth,extrapol);
-  std::cout<<"==== Simulate with E=4.3      ===="<<std::endl;
-  truth->SetPtEtaPhiM(4.3,0,0,0);
-  chain.simulate(*simulstate,truth,extrapol);
+  for(double E=1;E<10.1;E+=1) {
+    std::cout<<"==== Simulate with E="<<E<<"      ===="<<std::endl;
+    truth->SetPtEtaPhiM(E,0,0,0);
+    chain.simulate(*simulstate,truth,extrapol);
+  }  
   std::cout<<"==== Simulate with E=100      ===="<<std::endl;
   truth->SetPtEtaPhiM(100,0,0,0);
   chain.simulate(*simulstate,truth,extrapol);
   std::cout<<"==================================="<<std::endl<<std::endl;
+  std::cout<<"====== now with random bin ========"<<std::endl<<std::endl;
+  chain.set_DoRandomInterpolation();
+  for(double E=15;E<35.1;E+=4) {
+    std::cout<<"==== Simulate with E="<<E<<"      ===="<<std::endl;
+    truth->SetPtEtaPhiM(E,0,0,0);
+    for(int i=0;i<10;++i) chain.simulate(*simulstate,truth,extrapol);
+  }  
+  std::cout<<"==================================="<<std::endl<<std::endl;
 }
 
 
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationPDGIDSelectChain.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationPDGIDSelectChain.cxx
index c1b3bb66875bfb68db8e371a7942b97c53184568..92f95764f7d00e1fd53d77873ef06edc66d44ebf 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationPDGIDSelectChain.cxx
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSParametrizationPDGIDSelectChain.cxx
@@ -3,6 +3,7 @@
 */
 
 #include "ISF_FastCaloSimEvent/TFCSParametrizationPDGIDSelectChain.h"
+#include "ISF_FastCaloSimEvent/TFCSInvisibleParametrization.h"
 #include "ISF_FastCaloSimEvent/TFCSSimulationState.h"
 #include "ISF_FastCaloSimEvent/TFCSTruthState.h"
 #include "ISF_FastCaloSimEvent/TFCSExtrapolationState.h"
@@ -27,8 +28,9 @@ void TFCSParametrizationPDGIDSelectChain::simulate(TFCSSimulationState& simulsta
 {
   for(auto param: chain()) {
     if(param->is_match_pdgid(truth->pdgid())) {
-      ATH_MSG_DEBUG("pdgid="<<truth->pdgid()<<", now run: "<<param->GetName());
+      ATH_MSG_DEBUG("pdgid="<<truth->pdgid()<<", now run: "<<param->GetName()<< ((SimulateOnlyOnePDGID()==true) ? ", abort PDGID loop afterwards" : ", continue PDGID loop afterwards"));
       param->simulate(simulstate,truth,extrapol);
+      if(SimulateOnlyOnePDGID()) break;
     }  
   }
 }
@@ -43,33 +45,33 @@ void TFCSParametrizationPDGIDSelectChain::unit_test(TFCSSimulationState* simulst
   chain.setLevel(MSG::DEBUG);
 
   TFCSParametrization* param;
-  param=new TFCSParametrization("A begin all","A begin all");
+  param=new TFCSInvisibleParametrization("A begin all","A begin all");
   param->setLevel(MSG::DEBUG);
   param->set_pdgid(0);
   chain.push_back(param);
-  param=new TFCSParametrization("A end all","A end all");
+  param=new TFCSInvisibleParametrization("A end all","A end all");
   param->setLevel(MSG::DEBUG);
   param->set_pdgid(0);
   chain.push_back(param);
 
   for(int i=0;i<3;++i) {
-    param=new TFCSParametrization(Form("A%d",i),Form("A %d",i));
+    param=new TFCSInvisibleParametrization(Form("A%d",i),Form("A %d",i));
     param->setLevel(MSG::DEBUG);
     param->set_pdgid(i);
     chain.push_back(param);
   }
 
   for(int i=3;i>0;--i) {
-    param=new TFCSParametrization(Form("B%d",i),Form("B %d",i));
+    param=new TFCSInvisibleParametrization(Form("B%d",i),Form("B %d",i));
     param->setLevel(MSG::DEBUG);
     param->set_pdgid(i);
     chain.push_back(param);
   }
-  param=new TFCSParametrization("B end all","B end all");
+  param=new TFCSInvisibleParametrization("B end all","B end all");
   param->setLevel(MSG::DEBUG);
-  param->set_pdgid(1);
+  param->set_match_all_pdgid();
   chain.push_back(param);
-  param=new TFCSParametrization("B begin all","B begin all");
+  param=new TFCSInvisibleParametrization("B begin all","B begin all");
   param->setLevel(MSG::DEBUG);
   param->set_pdgid(1);
   chain.push_back(param);
@@ -85,6 +87,21 @@ void TFCSParametrizationPDGIDSelectChain::unit_test(TFCSSimulationState* simulst
   std::cout<<"==== Simulate with pdgid=2      ===="<<std::endl;
   truth->set_pdgid(2);
   chain.simulate(*simulstate,truth,extrapol);
-  std::cout<<"==================================="<<std::endl<<std::endl;
+  std::cout<<"====================================="<<std::endl<<std::endl;
+
+  std::cout<<"====================================="<<std::endl;
+  std::cout<<"= Now only one simul for each PDGID ="<<std::endl;
+  std::cout<<"====================================="<<std::endl;
+  chain.set_SimulateOnlyOnePDGID();
+  std::cout<<"==== Simulate with pdgid=0      ===="<<std::endl;
+  truth->set_pdgid(0);
+  chain.simulate(*simulstate,truth,extrapol);
+  std::cout<<"==== Simulate with pdgid=1      ===="<<std::endl;
+  truth->set_pdgid(1);
+  chain.simulate(*simulstate,truth,extrapol);
+  std::cout<<"==== Simulate with pdgid=2      ===="<<std::endl;
+  truth->set_pdgid(2);
+  chain.simulate(*simulstate,truth,extrapol);
+  
 }
 
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/CaloGeometryFromFile.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/CaloGeometryFromFile.cxx
index c7cabfe0a0a05addd307f0c84f96d3e550ee8d60..cbe656e1982d4cc0583ae743b99a58da731c18bd 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/CaloGeometryFromFile.cxx
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/CaloGeometryFromFile.cxx
@@ -17,7 +17,15 @@ map<unsigned long long, unsigned long long> g_cellId_vs_cellHashId_map;
 
 CaloGeometryFromFile::CaloGeometryFromFile() : CaloGeometry()
 {
-	ifstream textfile("/afs/cern.ch/atlas/groups/Simulation/FastCaloSimV2/cellId_vs_cellHashId_map.txt");
+}
+
+CaloGeometryFromFile::~CaloGeometryFromFile()
+{
+}
+
+bool CaloGeometryFromFile::LoadGeometryFromFile(TString filename,TString treename,TString hashfile)
+{
+	ifstream textfile(hashfile);
 	unsigned long long id, hash_id; 
 	cout << "Loading cellId_vs_cellHashId_map" << endl;
 	int i=0;
@@ -35,14 +43,7 @@ CaloGeometryFromFile::CaloGeometryFromFile() : CaloGeometry()
 	}
 	cout << "Done." << endl;
 	
-}
-
-CaloGeometryFromFile::~CaloGeometryFromFile()
-{
-}
 
-bool CaloGeometryFromFile::LoadGeometryFromFile(TString filename,TString treename)
-{
   TTree *tree;
   TFile *f = TFile::Open(filename);
   if(!f) return false;
@@ -103,7 +104,10 @@ bool CaloGeometryFromFile::LoadGeometryFromFile(TString filename,TString treenam
     if (ientry < 0) break;
     fChain->GetEntry(jentry);
     
-    if (g_cellId_vs_cellHashId_map.find(cell.m_identify)!=g_cellId_vs_cellHashId_map.end()) cell.m_hash_id=g_cellId_vs_cellHashId_map[cell.m_identify];
+    if (g_cellId_vs_cellHashId_map.find(cell.m_identify)!=g_cellId_vs_cellHashId_map.end()) {
+      cell.m_hash_id=g_cellId_vs_cellHashId_map[cell.m_identify];
+      if(cell.m_hash_id!=jentry) cout<<jentry<<" : ERROR hash="<<cell.m_hash_id<<endl;
+    }  
     else cout << endl << "ERROR: Cell id not found in the cellId_vs_cellHashId_map!!!" << endl << endl;
     
 
diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/CaloGeometryFromFile.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/CaloGeometryFromFile.h
index 254ebf0562f7802f99e4aae440eafdb43413b74b..75f10de4630bdafcdaacb65615080dca10506a3c 100644
--- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/CaloGeometryFromFile.h
+++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/CaloGeometryFromFile.h
@@ -6,13 +6,14 @@
 #define CaloGeometryFromFile_h
 
 #include "ISF_FastCaloSimParametrization/CaloGeometry.h"
+#include "TString.h"
 
 class CaloGeometryFromFile:public CaloGeometry {
 public :
    CaloGeometryFromFile();
    virtual ~CaloGeometryFromFile();
    
-   virtual bool LoadGeometryFromFile(TString filename,TString treename);
+   virtual bool LoadGeometryFromFile(TString filename,TString treename,TString hashfile="/afs/cern.ch/atlas/groups/Simulation/FastCaloSimV2/cellId_vs_cellHashId_map.txt");
    virtual void LoadFCalGeometryFromFiles(TString filename1,TString filename2,TString filename3);
    void DrawFCalGraph(int isam,int color);
 };