diff --git a/Framework/MDAlgorithms/inc/MantidMDAlgorithms/MDNormDirectSC.h b/Framework/MDAlgorithms/inc/MantidMDAlgorithms/MDNormDirectSC.h
index 8239e46772954ab534e42e4523f88c232dd0c18a..8c2825335fa316d406499f380aecd6458f662fbb 100644
--- a/Framework/MDAlgorithms/inc/MantidMDAlgorithms/MDNormDirectSC.h
+++ b/Framework/MDAlgorithms/inc/MantidMDAlgorithms/MDNormDirectSC.h
@@ -55,13 +55,13 @@ private:
   DataObjects::MDHistoWorkspace_sptr binInputWS();
   void createNormalizationWS(const DataObjects::MDHistoWorkspace &dataWS);
   std::vector<coord_t>
-  getValuesFromOtherDimensions(bool &skipNormalization) const;
+  getValuesFromOtherDimensions(bool &skipNormalization, uint16_t expInfoIndex=0) const;
   Kernel::Matrix<coord_t>
   findIntergratedDimensions(const std::vector<coord_t> &otherDimValues,
                             bool &skipNormalization);
   void cacheDimensionXValues();
   void calculateNormalization(const std::vector<coord_t> &otherValues,
-                              const Kernel::Matrix<coord_t> &affineTrans);
+                              const Kernel::Matrix<coord_t> &affineTrans, uint16_t expInfoIndex);
 
   void calculateIntersections(std::vector<std::array<double, 4>> &intersections,
                               const double theta, const double phi);
@@ -88,7 +88,10 @@ private:
   Kernel::V3D m_beamDir;
   /// ki-kf for Inelastic convention; kf-ki for Crystallography convention
   std::string convention;
+  /// internal flag to accumulate to an existing workspace
   bool m_accumulate{false};
+  /// number of experiment infos
+  uint16_t m_numExptInfos;
 };
 
 } // namespace MDAlgorithms
diff --git a/Framework/MDAlgorithms/src/MDNormDirectSC.cpp b/Framework/MDAlgorithms/src/MDNormDirectSC.cpp
index 2feb1ff6a126e39b27eea151c0dc5bcb4043bba2..97bd1c18f2c5f4fb703baf4f998d8a1c80f33ed0 100644
--- a/Framework/MDAlgorithms/src/MDNormDirectSC.cpp
+++ b/Framework/MDAlgorithms/src/MDNormDirectSC.cpp
@@ -142,22 +142,31 @@ void MDNormDirectSC::exec() {
   m_normWS->setDisplayNormalization(Mantid::API::NoNormalization);
   setProperty("OutputNormalizationWorkspace", m_normWS);
 
-  // Check for other dimensions if we could measure anything in the original
-  // data
-  bool skipNormalization = false;
-  const std::vector<coord_t> otherValues =
-      getValuesFromOtherDimensions(skipNormalization);
-  const auto affineTrans =
-      findIntergratedDimensions(otherValues, skipNormalization);
-  cacheDimensionXValues();
-
-  if (!skipNormalization) {
-    calculateNormalization(otherValues, affineTrans);
-  } else {
-    g_log.warning("Binning limits are outside the limits of the MDWorkspace. "
-                  "Not applying normalization.");
+  m_numExptInfos = outputWS->getNumExperimentInfo();
+  // loop over all experiment infos
+  for (uint16_t expInfoIndex=0; expInfoIndex < m_numExptInfos; expInfoIndex++)
+  {
+      if (expInfoIndex > 0) {
+          m_accumulate = true;
+      }
+      // Check for other dimensions if we could measure anything in the original
+      // data
+      bool skipNormalization = false;
+      const std::vector<coord_t> otherValues =
+          getValuesFromOtherDimensions(skipNormalization,expInfoIndex);
+      const auto affineTrans =
+          findIntergratedDimensions(otherValues, skipNormalization);
+      cacheDimensionXValues();
+
+      if (!skipNormalization) {
+        calculateNormalization(otherValues, affineTrans,expInfoIndex);
+      } else {
+        g_log.warning("Binning limits are outside the limits of the MDWorkspace. "
+                      "Not applying normalization.");
+      }
   }
 
+
   // Set the display normalization based on the input workspace
   outputWS->setDisplayNormalization(m_inputWS->displayNormalizationHisto());
 }
@@ -299,8 +308,8 @@ void MDNormDirectSC::createNormalizationWS(const MDHistoWorkspace &dataWS) {
  * MD position calculation
  */
 std::vector<coord_t>
-MDNormDirectSC::getValuesFromOtherDimensions(bool &skipNormalization) const {
-  const auto &runZero = m_inputWS->getExperimentInfo(0)->run();
+MDNormDirectSC::getValuesFromOtherDimensions(bool &skipNormalization, uint16_t expInfoIndex) const {
+  const auto &currentRun = m_inputWS->getExperimentInfo(expInfoIndex)->run();
 
   std::vector<coord_t> otherDimValues;
   for (size_t i = 4; i < m_inputWS->getNumDims(); i++) {
@@ -308,7 +317,7 @@ MDNormDirectSC::getValuesFromOtherDimensions(bool &skipNormalization) const {
     float dimMin = static_cast<float>(dimension->getMinimum());
     float dimMax = static_cast<float>(dimension->getMaximum());
     auto *dimProp = dynamic_cast<Kernel::TimeSeriesProperty<double> *>(
-        runZero.getProperty(dimension->getName()));
+        currentRun.getProperty(dimension->getName()));
     if (dimProp) {
       coord_t value = static_cast<coord_t>(dimProp->firstValue());
       otherDimValues.push_back(value);
@@ -396,7 +405,8 @@ Kernel::Matrix<coord_t> MDNormDirectSC::findIntergratedDimensions(
 }
 
 /**
- * Stores the X values from each H,K,L dimension as member variables
+ * Stores the X values from each H,K,L,E dimension as member variables
+ * Energy dimension is transformed to final wavevector.
  */
 void MDNormDirectSC::cacheDimensionXValues() {
   constexpr double energyToK = 8.0 * M_PI * M_PI *
@@ -442,17 +452,17 @@ void MDNormDirectSC::cacheDimensionXValues() {
  * @param otherValues
  * @param affineTrans
  */
-void MDNormDirectSC::calculateNormalization(
-    const std::vector<coord_t> &otherValues,
-    const Kernel::Matrix<coord_t> &affineTrans) {
+void MDNormDirectSC::calculateNormalization(const std::vector<coord_t> &otherValues,
+    const Kernel::Matrix<coord_t> &affineTrans,
+    uint16_t expInfoIndex) {
   constexpr double energyToK = 8.0 * M_PI * M_PI *
                                PhysicalConstants::NeutronMass *
                                PhysicalConstants::meV * 1e-20 /
                                (PhysicalConstants::h * PhysicalConstants::h);
-  const auto &exptInfoZero = *(m_inputWS->getExperimentInfo(0));
+  const auto &currentExptInfo = *(m_inputWS->getExperimentInfo(expInfoIndex));
   using VectorDoubleProperty = Kernel::PropertyWithValue<std::vector<double>>;
   auto *rubwLog =
-      dynamic_cast<VectorDoubleProperty *>(exptInfoZero.getLog("RUBW_MATRIX"));
+      dynamic_cast<VectorDoubleProperty *>(currentExptInfo.getLog("RUBW_MATRIX"));
   if (!rubwLog) {
     throw std::runtime_error(
         "Wokspace does not contain a log entry for the RUBW matrix."
@@ -460,12 +470,12 @@ void MDNormDirectSC::calculateNormalization(
   } else {
     Kernel::DblMatrix rubwValue(
         (*rubwLog)()); // includes the 2*pi factor but not goniometer for now :)
-    m_rubw = exptInfoZero.run().getGoniometerMatrix() * rubwValue;
+    m_rubw = currentExptInfo.run().getGoniometerMatrix() * rubwValue;
     m_rubw.Invert();
   }
-  const double protonCharge = exptInfoZero.run().getProtonCharge();
+  const double protonCharge = currentExptInfo.run().getProtonCharge();
 
-  const auto &spectrumInfo = exptInfoZero.spectrumInfo();
+  const auto &spectrumInfo = currentExptInfo.spectrumInfo();
 
   // Mapping
   const int64_t ndets = static_cast<int64_t>(spectrumInfo.size());
@@ -482,7 +492,8 @@ void MDNormDirectSC::calculateNormalization(
   std::vector<std::atomic<signal_t>> signalArray(m_normWS->getNPoints());
   std::vector<std::array<double, 4>> intersections;
   std::vector<coord_t> pos, posNew;
-  auto prog = make_unique<API::Progress>(this, 0.3, 1.0, ndets);
+  double progStep = 0.7/m_numExptInfos;
+  auto prog = make_unique<API::Progress>(this, 0.3 + progStep * expInfoIndex, 0.3 + progStep * (expInfoIndex + 1.), ndets);
   // cppcheck-suppress syntaxError
 PRAGMA_OMP(parallel for private(intersections, pos, posNew))
 for (int64_t i = 0; i < ndets; i++) {
diff --git a/docs/source/algorithms/MDNormDirectSC-v1.rst b/docs/source/algorithms/MDNormDirectSC-v1.rst
index 353dbdbe8a1ec881188b6b50e6605215dd215b6c..bb20fb152b42add92c41ab003a9c6debd0a8b8aa 100644
--- a/docs/source/algorithms/MDNormDirectSC-v1.rst
+++ b/docs/source/algorithms/MDNormDirectSC-v1.rst
@@ -14,17 +14,15 @@ The algorithm calculates a normalization MD workspace for single crystal direct
 Trajectories of each detector in reciprocal space are calculated, and the flux is integrated between intersections with each
 MDBox. A brief introduction to the multi-dimensional data normalization can be found :ref:`here <MDNorm>`.
 
-.. Note::
-
-    This is an experimental algorithm in Release 3.3. Please check the nightly Mantid build, and the Mantid webpage
-    for better offline help and usage examples.
-
 .. Note::
 
     If the MDEvent input workspace is generated from an event workspace, the algorithm gives the correct normalization
     only if the event workspace is cropped and binned to the same energy transfer range. If the workspace is not cropped, 
     one might have events in places where the normalization is calculated to be 0.
 
+.. Note::
+
+    As of :ref:`Release 3.14.0 <v3.14.0>`, the algorithm can handle merged MD workspaces. Make sure all original MDEvent workspaces have the same dimensions
 
 Usage
 -----
diff --git a/docs/source/release/v3.14.0/direct_inelastic.rst b/docs/source/release/v3.14.0/direct_inelastic.rst
index bd9955f081f110392f1c3585718f87aa85c49086..9188a52307dc45cf94f53aa19778a2670b6aa919 100644
--- a/docs/source/release/v3.14.0/direct_inelastic.rst
+++ b/docs/source/release/v3.14.0/direct_inelastic.rst
@@ -25,11 +25,8 @@ Improvements
 
 - Improved ``Save``-section of the TOFTOF reduction dialog.
 - Behavior of the :ref:`LoadDNSLegacy <algm-LoadDNSLegacy>` for TOF data has been changed: the algorithm does not try to guess elastic channel any more, but asks for the user input.
+- :ref:`LoadDNSSCD <algm-LoadDNSSCD>` has been improved to be able to load TOF data.
+- :ref:`MDNormDirectSC <algm-MDNormDirectSC>` now can handle merged MD workspaces.
 
 :ref:`Release 3.14.0 <v3.14.0>`
 
-
-Improvements
-############
-
-- :ref:`LoadDNSSCD <algm-LoadDNSSCD>` has been improved to be able to load TOF data.