diff --git a/v2.5/toolkit/fastnlotoolkit/fastNLOCoeffAddBase.cc b/v2.5/toolkit/fastnlotoolkit/fastNLOCoeffAddBase.cc
index 1bca96f556c5c113305317fa62a3ba2ed918381c..cfa0ed2068dd60b815a66cbb9484f053fea6878a 100644
--- a/v2.5/toolkit/fastnlotoolkit/fastNLOCoeffAddBase.cc
+++ b/v2.5/toolkit/fastnlotoolkit/fastNLOCoeffAddBase.cc
@@ -868,9 +868,23 @@ void fastNLOCoeffAddBase::ExtendX(int ObsBin, std::vector<fastNLOInterpolBase*>&
 }
 
 
+//________________________________________________________________________________________________________________ //
+void fastNLOCoeffAddBase::ExtendMuS(int ObsBin, std::vector<std::vector<fastNLOInterpolBase*>>& KernMuS) {
+   error[__func__] << "Method not implemented for fastNLOCoeffAddBase subclass. Exiting." << endl;
+   exit(1);
+}
+
+
+//________________________________________________________________________________________________________________ //
+void fastNLOCoeffAddBase::ExtendMu12(int ObsBin, std::vector<fastNLOInterpolBase*>& KernMu1, std::vector<fastNLOInterpolBase*>& KernMu2) {
+   error[__func__] << "Method not implemented for fastNLOCoeffAddBase subclass. Exiting." << endl;
+   exit(1);
+}
+
+
 //________________________________________________________________________________________________________________ //
 void fastNLOCoeffAddBase::ExtendSigmaTildeX(int ObsBin, unsigned int OldXSize1, unsigned int OldXSize2) {
-   error["ExtendSigmaTildeX"] << "Method not implemented for fastNLOCoeffAddBase subclass. Exiting." << endl;
+   error[__func__] << "Method not implemented for fastNLOCoeffAddBase subclass. Exiting." << endl;
    exit(1);
 }
 
@@ -878,6 +892,6 @@ void fastNLOCoeffAddBase::ExtendSigmaTildeX(int ObsBin, unsigned int OldXSize1,
 //________________________________________________________________________________________________________________ //
 void fastNLOCoeffAddBase::Fill(fnloEvent& Event, int ObsBin, int X, int scalevar,
       const std::vector<std::pair<int, double>>& nmu1, const std::vector<std::pair<int, double>>& nmu2, int SubProcess, double w) {
-   error["Fill"] << "Method not implemented for fastNLOCoeffAddBase subclass. Exiting." << endl;
+   error[__func__] << "Method not implemented for fastNLOCoeffAddBase subclass. Exiting." << endl;
    exit(1);
 }
diff --git a/v2.5/toolkit/fastnlotoolkit/fastNLOCreate.cc b/v2.5/toolkit/fastnlotoolkit/fastNLOCreate.cc
index da8b5dadf03be4ecbbfc772e0b3d139165f80220..8f61330e1beb05cda2a40890f577c0e6d23ee8bb 100644
--- a/v2.5/toolkit/fastnlotoolkit/fastNLOCreate.cc
+++ b/v2.5/toolkit/fastnlotoolkit/fastNLOCreate.cc
@@ -2969,8 +2969,10 @@ void fastNLOCreate::FillContribution(int scalevar, const double wgtfac) {
    if (fIsFlexibleScale) {
       nmu1 = fKernMu1[ObsBin]->GetNodeValues(fScenario._m1);
       nmu2 = fKernMu2[ObsBin]->GetNodeValues(fScenario._m2);
+      c->ExtendMu12(ObsBin, fKernMu1, fKernMu2);
    } else {
       nmu1  = fKernMuS[ObsBin][scalevar]->GetNodeValues(fScenario._m1 * fScaleFac[scalevar]);
+      c->ExtendMuS(ObsBin, fKernMuS);
    }
 
    if (fApplyPDFReweight) {
@@ -4633,8 +4635,8 @@ void  fastNLOCreate::InitInterpolationKernels() {
          }
       } else {
          // error
-         logger.error["InitInterpolationKernels"]<<"Cannot understand node counting: "<<fScenConsts.X_NNodeCounting<<"."<<endl;
-         logger.error["InitInterpolationKernels"]<<"Supported options are: 'NodesPerMagnitude', 'NodesPerBin', 'NodesMax' and 'NodeDensity'."<<endl;
+         logger.error["InitInterpolationKernels"]<<"Unknown value for X_NNodeCounting: "<<fScenConsts.X_NNodeCounting<<"."<<endl;
+         logger.error["InitInterpolationKernels"]<<"Supported values are: 'NodesPerMagnitude', 'NodesPerBin', 'NodesMax' and 'NodeDensity'."<<endl;
       }
 
       // Remove last node at x = 1; is multiplied by PDFs equalling zero anyway.
@@ -4647,22 +4649,69 @@ void  fastNLOCreate::InitInterpolationKernels() {
       // ------------------------------------------------
       int nqtot1 = fScenConsts.Mu1_NNodes;
       if (fIsFlexibleScale) {
-         logger.debug["InitInterpolationKernels"]<<"Make Mu1 grid for obsbin="<<i<<endl;
-         fKernMu1[i] = MakeInterpolationKernels(fScenConsts.Mu1_Kernel,wrmMu1Dn[i],wrmMu1Up[i],fScenConsts.Mu1_DistanceMeasure);
-         fKernMu1[i]->MakeGrids(nqtot1);
+         if (fScenConsts.Mu1_NNodeCounting == "NodesPerBin" ||
+             fScenConsts.Mu1_NNodeCounting == "NodesPerDim1" ||
+             fScenConsts.Mu1_NNodeCounting == "NodesPerDim2" ||
+             fScenConsts.Mu1_NNodeCounting == "NodesPerDim3") {
+            logger.debug["InitInterpolationKernels"]<<"Make Mu1 grid for obsbin="<<i<<endl;
+            fKernMu1[i] = MakeInterpolationKernels(fScenConsts.Mu1_Kernel, wrmMu1Dn[i], wrmMu1Up[i],
+                                                   fScenConsts.Mu1_DistanceMeasure);
+            fKernMu1[i]->MakeGrids(nqtot1);
+         } else if (fScenConsts.Mu1_NNodeCounting == "NodeDensity") {
+            fKernMu1[i] = MakeInterpolationKernels(fScenConsts.Mu1_Kernel, nqtot1, fScenConsts.Mu1_DistanceMeasure);
+            fKernMu1[i]->MakeGrids(4);
+         } else if (fScenConsts.Mu1_NNodeCounting != "NodesPerDim1" &&
+                    fScenConsts.Mu1_NNodeCounting != "NodesPerDim2" &&
+                    fScenConsts.Mu1_NNodeCounting != "NodesPerDim3"){
+            // error
+            logger.error["InitInterpolationKernels"]<<"Unknown value for Mu1_NNodeCounting: "<<fScenConsts.Mu1_NNodeCounting<<"."<<endl;
+            logger.error["InitInterpolationKernels"]<<"Supported values are: 'NodesPerBin', 'NodesPerDim1', 'NodesPerDim2', 'NodesPerDim3', and 'NodeDensity'."<<endl;
+         }
+
          // ------------------------------------------------
          // init scale2-interpolation kernels
          // ------------------------------------------------
          int nqtot2 = fScenConsts.Mu2_NNodes;
-         logger.debug["InitInterpolationKernels"]<<"Make Mu2 grid for obsbin="<<i<<endl;
-         fKernMu2[i] = MakeInterpolationKernels(fScenConsts.Mu2_Kernel,wrmMu2Dn[i],wrmMu2Up[i],fScenConsts.Mu2_DistanceMeasure);
-         fKernMu2[i]->MakeGrids(nqtot2);
+         if (fScenConsts.Mu2_NNodeCounting == "NodesPerBin" ||
+             fScenConsts.Mu2_NNodeCounting == "NodesPerDim1" ||
+             fScenConsts.Mu2_NNodeCounting == "NodesPerDim2" ||
+             fScenConsts.Mu2_NNodeCounting == "NodesPerDim3") {
+            logger.debug["InitInterpolationKernels"]<<"Make Mu2 grid for obsbin="<<i<<endl;
+            fKernMu2[i] = MakeInterpolationKernels(fScenConsts.Mu2_Kernel, wrmMu2Dn[i], wrmMu2Up[i],
+                                                   fScenConsts.Mu2_DistanceMeasure);
+            fKernMu2[i]->MakeGrids(nqtot2);
+         } else if (fScenConsts.Mu2_NNodeCounting == "NodeDensity") {
+            fKernMu2[i] = MakeInterpolationKernels(fScenConsts.Mu2_Kernel, nqtot2, fScenConsts.Mu2_DistanceMeasure);
+            fKernMu2[i]->MakeGrids(4);
+         } else if (fScenConsts.Mu2_NNodeCounting != "NodesPerDim1" &&
+                    fScenConsts.Mu2_NNodeCounting != "NodesPerDim2" &&
+                    fScenConsts.Mu2_NNodeCounting != "NodesPerDim3"){
+            // error
+            logger.error["InitInterpolationKernels"]<<"Unknown value for Mu2_NNodeCounting: "<<fScenConsts.Mu2_NNodeCounting<<"."<<endl;
+            logger.error["InitInterpolationKernels"]<<"Supported values are: 'NodesPerBin', 'NodesPerDim1', 'NodesPerDim2', 'NodesPerDim3', and 'NodeDensity'."<<endl;
+         }
       } else {
-         for (unsigned int k = 0 ; k<fScaleFac.size() ; k++) {
-            double muDn = fScaleFac[k]*wrmMu1Dn[i];
-            double muUp = fScaleFac[k]*wrmMu1Up[i];
-            fKernMuS[i][k] = MakeInterpolationKernels(fScenConsts.Mu1_Kernel,muDn,muUp,fScenConsts.Mu1_DistanceMeasure);
-            fKernMuS[i][k]->MakeGrids(nqtot1);
+         if (fScenConsts.Mu1_NNodeCounting == "NodesPerBin" ||
+             fScenConsts.Mu1_NNodeCounting == "NodesPerDim1" ||
+             fScenConsts.Mu1_NNodeCounting == "NodesPerDim2" ||
+             fScenConsts.Mu1_NNodeCounting == "NodesPerDim3") {
+            for (unsigned int k = 0 ; k<fScaleFac.size() ; k++) {
+               double muDn = fScaleFac[k]*wrmMu1Dn[i];
+               double muUp = fScaleFac[k]*wrmMu1Up[i];
+               fKernMuS[i][k] = MakeInterpolationKernels(fScenConsts.Mu1_Kernel,muDn,muUp,fScenConsts.Mu1_DistanceMeasure);
+               fKernMuS[i][k]->MakeGrids(nqtot1);
+            }
+         } else if (fScenConsts.Mu1_NNodeCounting == "NodeDensity") {
+            for (unsigned int k = 0 ; k<fScaleFac.size() ; k++) {
+               fKernMuS[i][k] = MakeInterpolationKernels(fScenConsts.Mu1_Kernel, nqtot1, fScenConsts.Mu1_DistanceMeasure);
+               fKernMuS[i][k]->MakeGrids(4);
+            }
+         } else if (fScenConsts.Mu1_NNodeCounting != "NodesPerDim1" &&
+                    fScenConsts.Mu1_NNodeCounting != "NodesPerDim2" &&
+                    fScenConsts.Mu1_NNodeCounting != "NodesPerDim3"){
+            // error
+            logger.error["InitInterpolationKernels"]<<"Unknown value for Mu1_NNodeCounting: "<<fScenConsts.Mu1_NNodeCounting<<"."<<endl;
+            logger.error["InitInterpolationKernels"]<<"Supported values are: 'NodesPerBin', 'NodesPerDim1', 'NodesPerDim2', 'NodesPerDim3', and 'NodeDensity'."<<endl;
          }
       }
    }
diff --git a/v2.5/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOCoeffAddBase.h b/v2.5/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOCoeffAddBase.h
index 9023a4089c290d42e6f4ecab8407d669aebecd3f..592886897a7b10a7a8b86650affb003a6ee60605 100644
--- a/v2.5/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOCoeffAddBase.h
+++ b/v2.5/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOCoeffAddBase.h
@@ -180,6 +180,8 @@ public:
 
    void ExtendX(int ObsBin, std::vector<fastNLOInterpolBase*>& KernX1, std::vector<fastNLOInterpolBase*>& KernX2);
    virtual void ExtendSigmaTildeX(int ObsBin, unsigned int OldXSize1, unsigned int OldXSize2);
+   virtual void ExtendMuS(int ObsBin, std::vector<std::vector<fastNLOInterpolBase*>>& KernMuS);
+   virtual void ExtendMu12(int ObsBin, std::vector<fastNLOInterpolBase*>& KernMu1, std::vector<fastNLOInterpolBase*>& KernMu2);
    virtual void Fill(fnloEvent& Event, int ObsBin, int X, int scalevar, const std::vector<std::pair<int, double>>& nmu1,
       const std::vector<std::pair<int, double>>& nmu2, int SubProcess, double w);