diff --git a/src/modelchecker/DFTAnalyser.h b/src/modelchecker/DFTAnalyser.h
new file mode 100644
index 000000000..93bc7c7b1
--- /dev/null
+++ b/src/modelchecker/DFTAnalyser.h
@@ -0,0 +1,158 @@
+#pragma once
+
+#include "logic/Formula.h"
+#include "parser/DFTGalileoParser.h"
+#include "builder/ExplicitDFTModelBuilder.h"
+#include "modelchecker/results/CheckResult.h"
+#include "utility/storm.h"
+#include "storage/dft/DFTIsomorphism.h"
+
+
+#include <chrono>
+
+template<typename ValueType>
+class DFTAnalyser {
+    
+    std::chrono::duration<double> buildingTime = std::chrono::duration<double>::zero();
+    std::chrono::duration<double> explorationTime = std::chrono::duration<double>::zero();
+    std::chrono::duration<double> bisimulationTime = std::chrono::duration<double>::zero();
+    std::chrono::duration<double> modelCheckingTime = std::chrono::duration<double>::zero();
+    std::chrono::duration<double> totalTime = std::chrono::duration<double>::zero();
+    ValueType checkResult = storm::utility::zero<ValueType>();
+public:
+    void check(storm::storage::DFT<ValueType> dft , std::shared_ptr<const storm::logic::Formula> const& formula, bool symred = true, bool allowModularisation = true) {
+        
+        // Building DFT
+        std::chrono::high_resolution_clock::time_point totalStart = std::chrono::high_resolution_clock::now();
+        
+        
+        // Optimizing DFT
+        dft = dft.optimize();
+        std::vector<storm::storage::DFT<ValueType>> dfts;
+        checkResult = checkHelper(dft, formula, symred, allowModularisation);
+        totalTime = std::chrono::high_resolution_clock::now() - totalStart;
+            
+        
+            
+    }
+    
+private:
+    ValueType checkHelper(storm::storage::DFT<ValueType> const& dft , std::shared_ptr<const storm::logic::Formula> const& formula, bool symred = true, bool allowModularisation = true)  {
+        std::cout << "check helper called" << std::endl;
+        bool invResults = false;
+        std::vector<storm::storage::DFT<ValueType>> dfts = {dft};
+        std::vector<ValueType> res;
+        
+        if(allowModularisation) {
+            std::cout << dft.topLevelType() << std::endl;
+            if(dft.topLevelType() == storm::storage::DFTElementType::AND) {
+                std::cout << "top modularisation called AND" << std::endl;
+                dfts = dft.topModularisation();
+                if(dfts.size() > 1) {
+                    for(auto const ft : dfts) {
+                        res.push_back(checkHelper(ft, formula, symred, allowModularisation));
+                    }
+                }
+            }
+            if(dft.topLevelType() == storm::storage::DFTElementType::OR) {
+                std::cout << "top modularisation called OR" << std::endl;
+                dfts = dft.topModularisation();
+                if(dfts.size() > 1) {
+                    for(auto const ft : dfts) {
+                        res.push_back(checkHelper(ft, formula, symred, allowModularisation));
+                    }
+                }
+                invResults = true;
+            }
+        } 
+        if(res.empty()) {
+            // Model based modularisation.
+            auto const& models = buildMarkovModels(dfts, formula, symred);
+
+
+            for (auto const& model : models) {
+                // Model checking
+                std::cout << "Model checking..." << std::endl;
+                std::chrono::high_resolution_clock::time_point modelCheckingStart = std::chrono::high_resolution_clock::now();
+                std::unique_ptr<storm::modelchecker::CheckResult> result(storm::verifySparseModel(model, formula));
+                std::cout << "done." << std::endl;
+                assert(result);
+                result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates()));
+                modelCheckingTime += std::chrono::high_resolution_clock::now() -  modelCheckingStart;
+                res.push_back(result->asExplicitQuantitativeCheckResult<ValueType>().getValueMap().begin()->second);
+            }
+        }
+        
+        ValueType result = storm::utility::one<ValueType>();
+        for(auto const& r : res) {
+            if(invResults) {
+                result *= storm::utility::one<ValueType>() - r;
+            } else {
+                result *= r;
+            }
+        }
+        if(invResults) {
+            return storm::utility::one<ValueType>() - result;
+        }
+        return result;
+    }
+    
+    
+    std::vector<std::shared_ptr<storm::models::sparse::Model<ValueType>>> buildMarkovModels(std::vector<storm::storage::DFT<ValueType>> const&  dfts,   std::shared_ptr<const storm::logic::Formula> const& formula, bool symred) {
+        std::vector<std::shared_ptr<storm::models::sparse::Model<ValueType>>> models;
+        for(auto& dft : dfts) {
+            std::chrono::high_resolution_clock::time_point buildingStart = std::chrono::high_resolution_clock::now();
+
+            std::map<size_t, std::vector<std::vector<size_t>>> emptySymmetry;
+            storm::storage::DFTIndependentSymmetries symmetries(emptySymmetry);
+            if(symred) {
+                auto colouring = dft.colourDFT();
+                symmetries = dft.findSymmetries(colouring);
+                std::cout << "Found " << symmetries.groups.size() << " symmetries." << std::endl;
+                STORM_LOG_TRACE("Symmetries: " << std::endl << symmetries);
+            }
+            std::chrono::high_resolution_clock::time_point buildingEnd = std::chrono::high_resolution_clock::now();
+            buildingTime += buildingEnd - buildingStart;
+            
+            // Building Markov Automaton
+            std::cout << "Building Model..." << std::endl;
+            storm::builder::ExplicitDFTModelBuilder<ValueType> builder(dft, symmetries);
+            typename storm::builder::ExplicitDFTModelBuilder<ValueType>::LabelOptions labeloptions; // TODO initialize this with the formula
+            std::shared_ptr<storm::models::sparse::Model<ValueType>> model = builder.buildModel(labeloptions);
+            //model->printModelInformationToStream(std::cout);
+            std::cout << "No. states (Explored): " << model->getNumberOfStates() << std::endl;
+            std::cout << "No. transitions (Explored): " << model->getNumberOfTransitions() << std::endl;
+            std::chrono::high_resolution_clock::time_point explorationEnd = std::chrono::high_resolution_clock::now();
+            explorationTime += explorationEnd -buildingEnd;
+
+            // Bisimulation
+            if (model->isOfType(storm::models::ModelType::Ctmc)) {
+                std::cout << "Bisimulation..." << std::endl;
+                model =  storm::performDeterministicSparseBisimulationMinimization<storm::models::sparse::Ctmc<ValueType>>(model->template as<storm::models::sparse::Ctmc<ValueType>>(), {formula}, storm::storage::BisimulationType::Weak)->template as<storm::models::sparse::Ctmc<ValueType>>();
+                //model->printModelInformationToStream(std::cout);
+            }
+            std::cout << "No. states (Bisimulation): " << model->getNumberOfStates() << std::endl;
+            std::cout << "No. transitions (Bisimulation): " << model->getNumberOfTransitions() << std::endl;
+            bisimulationTime += std::chrono::high_resolution_clock::now() - explorationEnd;
+            models.push_back(model);
+        }
+        return models;
+    }
+    
+public:
+    void printTimings(std::ostream& os = std::cout) {
+        os << "Times:" << std::endl;
+        os << "Building:\t" << buildingTime.count() << std::endl;
+        os << "Exploration:\t" << explorationTime.count() << std::endl;
+        os << "Bisimulation:\t" << bisimulationTime.count() << std::endl;
+        os << "Modelchecking:\t" << modelCheckingTime.count() << std::endl;
+        os << "Total:\t\t" << totalTime.count() << std::endl;
+    }
+    
+    void printResult(std::ostream& os = std::cout) {
+        
+        os << "Result: [";
+        os << checkResult << "]" << std::endl;
+    }
+    
+};
diff --git a/src/storage/dft/DFT.cpp b/src/storage/dft/DFT.cpp
index dd771f230..482e8203c 100644
--- a/src/storage/dft/DFT.cpp
+++ b/src/storage/dft/DFT.cpp
@@ -247,6 +247,46 @@ namespace storm {
             return stateIndex;
         }
         
+        template<typename ValueType>
+        std::vector<DFT<ValueType>>  DFT<ValueType>::topModularisation() const {
+            assert(isGate(mTopLevelIndex));
+            auto const& children = getGate(mTopLevelIndex)->children();
+            std::map<size_t, std::vector<size_t>> subdfts;
+            for(auto const& child : children) {
+                std::vector<size_t> isubdft;
+                if(child->nrParents() > 1 || child->hasOutgoingDependencies()) {
+                    return {*this};
+                }
+                if (isGate(child->id())) {
+                    isubdft = getGate(child->id())->independentSubDft(false);
+                } else {
+                    assert(isBasicElement(child->id()));
+                    if(!getBasicElement(child->id())->hasIngoingDependencies()) {
+                        isubdft = {child->id()};
+                    }
+                    
+                }
+                if(isubdft.empty()) {
+                    return {*this};
+                } else {
+                    subdfts[child->id()] = isubdft;
+                }
+            }
+            
+            std::vector<DFT<ValueType>> res;
+            for(auto const& subdft : subdfts) {
+                DFTBuilder<ValueType> builder;
+            
+                for(size_t id : subdft.second) {
+                    builder.copyElement(mElements[id]);
+                }
+                builder.setTopLevel(mElements[subdft.first]->name());
+                res.push_back(builder.build());
+            }
+            return res;
+            
+        }
+        
         template<typename ValueType>
         DFT<ValueType> DFT<ValueType>::optimize() const {
             std::vector<size_t> modIdea = findModularisationRewrite();
diff --git a/src/storage/dft/DFT.h b/src/storage/dft/DFT.h
index f28d244f0..0f4832fb5 100644
--- a/src/storage/dft/DFT.h
+++ b/src/storage/dft/DFT.h
@@ -6,6 +6,7 @@
 #include <unordered_map>
 #include <list>
 #include <map>
+#include <vector>
 
 #include <boost/iterator/counting_iterator.hpp>
 
@@ -95,6 +96,10 @@ namespace storm {
                 return mTopLevelIndex;
             }
             
+            DFTElementType topLevelType() const {
+                return mElements[getTopLevelIndex()]->type();
+            }
+            
             size_t getMaxSpareChildCount() const {
                 return mMaxSpareChildCount;
             }
@@ -187,6 +192,8 @@ namespace storm {
                 return elements;
             }
             
+            std::vector<DFT<ValueType>> topModularisation() const;
+            
             bool isRepresentative(size_t id) const {
                 for (auto const& parent : getElement(id)->parents()) {
                     if (parent->isSpareGate()) {
diff --git a/src/storm-dyftee.cpp b/src/storm-dyftee.cpp
index 1f170072b..f280bbebf 100644
--- a/src/storm-dyftee.cpp
+++ b/src/storm-dyftee.cpp
@@ -1,11 +1,7 @@
 #include "logic/Formula.h"
-#include "parser/DFTGalileoParser.h"
 #include "utility/initialize.h"
-#include "builder/ExplicitDFTModelBuilder.h"
-#include "modelchecker/results/CheckResult.h"
 #include "utility/storm.h"
-#include "storage/dft/DFTIsomorphism.h"
-
+#include "modelchecker/DFTAnalyser.h"
 #include <boost/lexical_cast.hpp>
 
 /*!
@@ -15,71 +11,20 @@
  * @param property PCTC formula capturing the property to check.
  */
 template <typename ValueType>
-void analyzeDFT(std::string filename, std::string property, bool symred = false) {
+void analyzeDFT(std::string filename, std::string property, bool symred = false, bool allowModularisation = false) {
     storm::settings::SettingsManager& manager = storm::settings::mutableManager();
     manager.setFromString("");
 
-    // Building DFT
-    std::chrono::high_resolution_clock::time_point buildingStart = std::chrono::high_resolution_clock::now();
     storm::parser::DFTGalileoParser<ValueType> parser;
     storm::storage::DFT<ValueType> dft = parser.parseDFT(filename);
     std::vector<std::shared_ptr<storm::logic::Formula>> parsedFormulas = storm::parseFormulasForExplicit(property);
     std::vector<std::shared_ptr<const storm::logic::Formula>> formulas(parsedFormulas.begin(), parsedFormulas.end());
     assert(formulas.size() == 1);
     
-    // Optimizing DFT
-    dft = dft.optimize();
-    std::map<size_t, std::vector<std::vector<size_t>>> emptySymmetry;
-    storm::storage::DFTIndependentSymmetries symmetries(emptySymmetry);
-    if(symred) {
-        auto colouring = dft.colourDFT();
-        symmetries = dft.findSymmetries(colouring);
-        std::cout << "Found " << symmetries.groups.size() << " symmetries." << std::endl;
-        STORM_LOG_TRACE("Symmetries: " << std::endl << symmetries);
-    }
-    std::chrono::high_resolution_clock::time_point buildingEnd = std::chrono::high_resolution_clock::now();
-    
-    // Building Markov Automaton
-    std::cout << "Building Model..." << std::endl;
-    storm::builder::ExplicitDFTModelBuilder<ValueType> builder(dft, symmetries);
-    typename storm::builder::ExplicitDFTModelBuilder<ValueType>::LabelOptions labeloptions; // TODO initialize this with the formula
-    std::shared_ptr<storm::models::sparse::Model<ValueType>> model = builder.buildModel(labeloptions);
-    //model->printModelInformationToStream(std::cout);
-    std::cout << "No. states (Explored): " << model->getNumberOfStates() << std::endl;
-    std::cout << "No. transitions (Explored): " << model->getNumberOfTransitions() << std::endl;
-    std::chrono::high_resolution_clock::time_point explorationEnd = std::chrono::high_resolution_clock::now();
-    
-    // Bisimulation
-    if (model->isOfType(storm::models::ModelType::Ctmc)) {
-        std::cout << "Bisimulation..." << std::endl;
-        model =  storm::performDeterministicSparseBisimulationMinimization<storm::models::sparse::Ctmc<ValueType>>(model->template as<storm::models::sparse::Ctmc<ValueType>>(), formulas, storm::storage::BisimulationType::Weak)->template as<storm::models::sparse::Ctmc<ValueType>>();
-        //model->printModelInformationToStream(std::cout);
-    }
-    std::cout << "No. states (Bisimulation): " << model->getNumberOfStates() << std::endl;
-    std::cout << "No. transitions (Bisimulation): " << model->getNumberOfTransitions() << std::endl;
-    std::chrono::high_resolution_clock::time_point bisimulationEnd = std::chrono::high_resolution_clock::now();
-    
-    // Model checking
-    std::cout << "Model checking..." << std::endl;
-    std::unique_ptr<storm::modelchecker::CheckResult> result(storm::verifySparseModel(model, formulas[0]));
-    assert(result);
-    result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates()));
-    
-    // Times
-    std::chrono::high_resolution_clock::time_point modelCheckingEnd = std::chrono::high_resolution_clock::now();
-    std::chrono::duration<double> buildingTime = buildingEnd - buildingStart;
-    std::chrono::duration<double> explorationTime = explorationEnd - buildingEnd;
-    std::chrono::duration<double> bisimulationTime = bisimulationEnd - explorationEnd;
-    std::chrono::duration<double> modelCheckingTime = modelCheckingEnd - bisimulationEnd;
-    std::chrono::duration<double> totalTime = modelCheckingEnd - buildingStart;
-    std::cout << "Times:" << std::endl;
-    std::cout << "Building:\t" << buildingTime.count() << std::endl;
-    std::cout << "Exploration:\t" << explorationTime.count() << std::endl;
-    std::cout << "Bisimulation:\t" << bisimulationTime.count() << std::endl;
-    std::cout << "Modelchecking:\t" << modelCheckingTime.count() << std::endl;
-    std::cout << "Total:\t\t" << totalTime.count() << std::endl;
-    std::cout << "Result: ";
-    std::cout << *result << std::endl;
+    DFTAnalyser<ValueType> analyser;
+    analyser.check(dft, formulas[0], symred, allowModularisation);
+    analyser.printTimings();
+    analyser.printResult();
 }
 
 /*!
@@ -100,6 +45,8 @@ int main(int argc, char** argv) {
     bool parametric = false;
     bool symred = false;
     bool minimal = true;
+    bool allowModular = true;
+    bool enableModularisation = false;
     std::string filename = argv[1];
     std::string operatorType = "";
     std::string targetFormula = "";
@@ -112,6 +59,7 @@ int main(int argc, char** argv) {
             assert(targetFormula.empty());
             operatorType = "ET";
             targetFormula = "F \"failed\"";
+            allowModular = false;
         } else if (option == "--probability") {
             assert(targetFormula.empty());
             operatorType = "P";;
@@ -142,6 +90,8 @@ int main(int argc, char** argv) {
             pctlFormula = argv[i];
         } else if (option == "--symred") {
             symred = true;
+        } else if (option == "--modularisation") {
+            enableModularisation = true;
         } else if (option == "--min") {
             minimal = true;
         } else if (option == "--max") {
@@ -164,8 +114,8 @@ int main(int argc, char** argv) {
     std::cout << "Running " << (parametric ? "parametric " : "") << "DFT analysis on file " << filename << " with property " << pctlFormula << std::endl;
 
     if (parametric) {
-        analyzeDFT<storm::RationalFunction>(filename, pctlFormula, symred);
+        analyzeDFT<storm::RationalFunction>(filename, pctlFormula, symred, allowModular && enableModularisation );
     } else {
-        analyzeDFT<double>(filename, pctlFormula, symred);
+        analyzeDFT<double>(filename, pctlFormula, symred, allowModular && enableModularisation);
     }
 }