diff --git a/Anima/diffusion/mcm_tools/CMakeLists.txt b/Anima/diffusion/mcm_tools/CMakeLists.txt
index 9c086625c..09a7309a8 100644
--- a/Anima/diffusion/mcm_tools/CMakeLists.txt
+++ b/Anima/diffusion/mcm_tools/CMakeLists.txt
@@ -4,6 +4,7 @@ add_subdirectory(get_scalar_map_from_ddi)
 add_subdirectory(mcm_average_images)
 add_subdirectory(mcm_merge_anisotropic_weights)
 add_subdirectory(mcm_merge_block_images)
+add_subdirectory(mcm_orientation_priors)
 add_subdirectory(mcm_scalar_maps)
 add_subdirectory(mt_estimation_validation)
 add_subdirectory(test_averaging)
\ No newline at end of file
diff --git a/Anima/diffusion/mcm_tools/mcm_orientation_priors/CMakeLists.txt b/Anima/diffusion/mcm_tools/mcm_orientation_priors/CMakeLists.txt
new file mode 100644
index 000000000..8db408adc
--- /dev/null
+++ b/Anima/diffusion/mcm_tools/mcm_orientation_priors/CMakeLists.txt
@@ -0,0 +1,38 @@
+if(BUILD_TOOLS)
+
+project(animaMCMOrientationPriors)
+
+## #############################################################################
+## List Sources
+## #############################################################################
+
+list_source_files(${PROJECT_NAME}
+  ${CMAKE_CURRENT_SOURCE_DIR}
+  )
+
+## #############################################################################
+## add executable
+## #############################################################################
+
+add_executable(${PROJECT_NAME}
+  ${${PROJECT_NAME}_CFILES}
+  )
+
+## #############################################################################
+## Link
+## #############################################################################
+
+target_link_libraries(${PROJECT_NAME}
+  ${ITKIO_LIBRARIES}
+  ${TinyXML2_LIBRARY}
+  AnimaMCM
+  AnimaStatisticalDistributions
+  )
+
+## #############################################################################
+## install
+## #############################################################################
+
+set_exe_install_rules(${PROJECT_NAME})
+
+endif()
diff --git a/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriors.cxx b/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriors.cxx
new file mode 100644
index 000000000..8f31243f8
--- /dev/null
+++ b/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriors.cxx
@@ -0,0 +1,117 @@
+#include <animaMCMFileReader.h>
+#include <animaReadWriteFunctions.h>
+#include <animaMCMOrientationPriorsImageFilter.h>
+
+#include <itkTimeProbe.h>
+#include <tclap/CmdLine.h>
+
+int main(int argc,  char **argv)
+{
+    TCLAP::CmdLine cmd("INRIA / IRISA - VisAGeS Team", ' ',ANIMA_VERSION);
+
+    TCLAP::ValueArg<std::string> inArg(
+        "i", "input-mcms", 
+        "A text file specifying a list of MCM images in a common geometry.", 
+        true, "", "input MCM images", cmd
+    );
+    TCLAP::ValueArg<std::string> outOrientationArg(
+        "o", "output-orientation", 
+        "A string specifying the basename for the output vector images that will store priors on the orientations.", 
+        true, "", "output orientation priors",cmd
+    );
+    TCLAP::ValueArg<std::string> outWeightsArg(
+        "w", "output-weights", 
+        "A string specifying the filename for the output vector image that will store priors on the weights.", 
+        true, "", "output weights priors", cmd
+    );
+    
+    TCLAP::ValueArg<std::string> maskArg(
+        "m", "input-masks", 
+        "A text file specifying a list of mask images in the same common geometry as the input MCM images (default: none).", 
+        false, "", "input mask images", cmd
+    );
+    TCLAP::ValueArg<unsigned int> nbThreadsArg(
+        "T", "nthreads", 
+        "An integer value specifying the number of threads to run on (default: all cores).", 
+        false, itk::MultiThreaderBase::GetGlobalDefaultNumberOfThreads(), "number of threads", cmd
+    );
+
+    try
+    {
+        cmd.parse(argc,argv);
+    }
+    catch (TCLAP::ArgException& e)
+    {
+        std::cerr << "Error: " << e.error() << "for argument " << e.argId() << std::endl;
+        return EXIT_FAILURE;
+    }
+
+    using MainFilterType = anima::MCMOrientationPriorsImageFilter <double>;
+    MainFilterType::Pointer mainFilter = MainFilterType::New();
+
+    using MCMReaderType = anima::MCMFileReader <double,3>;
+    using MaskImageType = MainFilterType::MaskImageType;
+    using OutputImageType = MainFilterType::OutputImageType;
+
+    // Load MCM images
+    std::ifstream inputFile(inArg.getValue().c_str());
+
+    if (!inputFile.is_open())
+    {
+        std::cerr << "Please provide usable file with input MCMs" << std::endl;
+        return EXIT_FAILURE;
+    }
+
+    unsigned int nbOfImages = 0;
+
+    while (!inputFile.eof())
+    {
+        char tmpStr[2048];
+        inputFile.getline(tmpStr,2048);
+
+        if (strcmp(tmpStr,"") == 0)
+            continue;
+
+        std::cout << "Loading image " << nbOfImages + 1 << ": " << tmpStr << std::endl;
+
+        MCMReaderType mcmReader;
+        mcmReader.SetFileName(tmpStr);
+        mcmReader.Update();
+
+        mainFilter->SetInput(nbOfImages, mcmReader.GetModelVectorImage());
+
+        nbOfImages++;
+    }
+
+    std::ifstream masksIn;
+    if (maskArg.getValue() != "")
+        masksIn.open(maskArg.getValue());
+
+    if (masksIn.is_open())
+    {
+        char tmpStr[2048];
+        while (!masksIn.eof())
+        {
+            masksIn.getline(tmpStr,2048);
+
+            if (strcmp(tmpStr,"") == 0)
+                continue;
+
+            mainFilter->AddMaskImage(anima::readImage <MaskImageType> (tmpStr));
+        }
+    }
+
+    mainFilter->SetNumberOfWorkUnits(nbThreadsArg.getValue());
+    mainFilter->Update();
+
+    unsigned int numberOfAnisotropicCompartments = mainFilter->GetNumberOfAnisotropicCompartments();
+    for (unsigned int i = 0;i < numberOfAnisotropicCompartments;++i)
+    {
+        std::string fileName = outOrientationArg.getValue() + "_" + std::to_string(i) + ".nrrd";
+        anima::writeImage <OutputImageType> (fileName, mainFilter->GetOutput(i));
+    }
+    
+    anima::writeImage <OutputImageType> (outWeightsArg.getValue(), mainFilter->GetOutput(numberOfAnisotropicCompartments));
+
+    return EXIT_SUCCESS;
+}
diff --git a/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriorsImageFilter.h b/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriorsImageFilter.h
new file mode 100644
index 000000000..168982e07
--- /dev/null
+++ b/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriorsImageFilter.h
@@ -0,0 +1,83 @@
+#pragma once
+
+#include <animaMCMImage.h>
+#include <animaMultiCompartmentModel.h>
+
+#include <itkImageToImageFilter.h>
+#include <itkVectorImage.h>
+
+namespace anima
+{
+
+template <class TPixelType>
+class MCMOrientationPriorsImageFilter :
+public itk::ImageToImageFilter< anima::MCMImage <TPixelType, 3>, itk::VectorImage<TPixelType, 3> >
+{
+public:
+    /** Standard class type def */
+
+    using Self = MCMOrientationPriorsImageFilter;
+    using InputImageType = anima::MCMImage <TPixelType, 3>;
+    using OutputImageType = itk::VectorImage <TPixelType, 3>;
+    using Superclass = itk::ImageToImageFilter <InputImageType, OutputImageType >;
+    using Pointer = itk::SmartPointer<Self>;
+    using ConstPointer = itk::SmartPointer<const Self>;
+
+    /** Method for creation through the object factory. */
+    itkNewMacro(Self)
+
+    /** Run-time type information (and related methods) */
+    itkTypeMacro(MCMOrientationPriorsImageFilter, ImageToImageFilter)
+
+    using InputImagePointer = typename InputImageType::ConstPointer;
+    using InputRegionType = typename InputImageType::RegionType;
+
+    using OutputImagePointer = typename OutputImageType::Pointer;
+    using OutputPixelType = typename OutputImageType::PixelType;
+
+    using MaskImageType = itk::Image <unsigned char, 3>;
+    using MaskImagePointer = MaskImageType::Pointer;
+
+    // Multi-compartment models typedefs
+    using MCModelType = anima::MultiCompartmentModel;
+    using MCModelPointer = MCModelType::Pointer;
+
+    void AddMaskImage(MaskImageType *maskImage) {m_MaskImages.push_back(maskImage);}
+    unsigned int GetNumberOfAnisotropicCompartments() {return m_NumberOfAnisotropicCompartments;}
+
+protected:
+    MCMOrientationPriorsImageFilter()
+    {
+        m_ReferenceInputModels.clear();
+        m_MaskImages.clear();
+        m_NumberOfAnisotropicCompartments = 0;
+    }
+
+    virtual ~MCMOrientationPriorsImageFilter() {}
+
+    bool isZero(const itk::VariableLengthVector <double> &value) const
+    {
+        for (unsigned int i = 0;i < value.GetNumberOfElements();++i)
+        {
+            if (value[i] != 0)
+                return false;
+        }
+
+        return true;
+    }
+
+    void GenerateOutputInformation() ITK_OVERRIDE;
+    void BeforeThreadedGenerateData() ITK_OVERRIDE;
+    void DynamicThreadedGenerateData(const InputRegionType &region) ITK_OVERRIDE;
+
+private:
+    ITK_DISALLOW_COPY_AND_ASSIGN(MCMOrientationPriorsImageFilter);
+
+    std::vector <MCModelPointer> m_ReferenceInputModels;
+    std::vector <MaskImagePointer> m_MaskImages;
+    unsigned int m_NumberOfAnisotropicCompartments;
+};
+
+} // end namespace anima
+
+#include "animaMCMOrientationPriorsImageFilter.hxx"
diff --git a/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriorsImageFilter.hxx b/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriorsImageFilter.hxx
new file mode 100644
index 000000000..ad94a0de0
--- /dev/null
+++ b/Anima/diffusion/mcm_tools/mcm_orientation_priors/animaMCMOrientationPriorsImageFilter.hxx
@@ -0,0 +1,330 @@
+#pragma once
+#include "animaMCMOrientationPriorsImageFilter.h"
+#include <animaSpectralClusteringFilter.h>
+#include <animaDirichletDistribution.h>
+#include <animaWatsonDistribution.h>
+
+#include <itkImageRegionConstIterator.h>
+#include <itkImageRegionIterator.h>
+
+namespace anima
+{
+
+template <class TPixelType>
+void
+MCMOrientationPriorsImageFilter <TPixelType>
+::GenerateOutputInformation()
+{
+    //-------------------------------
+    // Creates the additional outputs
+    //-------------------------------
+    unsigned int prevNum = this->GetNumberOfOutputs();
+
+    InputImageType *input = const_cast <InputImageType *> (this->GetInput(0));
+    MCModelPointer tmpMCM = input->GetDescriptionModel();
+    unsigned int numIsoCompartments = tmpMCM->GetNumberOfIsotropicCompartments();
+    unsigned int numCompartments = tmpMCM->GetNumberOfCompartments();
+    m_NumberOfAnisotropicCompartments = numCompartments - numIsoCompartments;
+    unsigned int numOutputs = m_NumberOfAnisotropicCompartments + 1;
+    
+    this->SetNumberOfIndexedOutputs(numOutputs);
+    
+    for (unsigned int i = prevNum;i < numOutputs;++i)
+        this->SetNthOutput(i,this->MakeOutput(i).GetPointer());
+    
+    //---------------------------------------------------
+    // Call the superclass' implementation of this method
+    //---------------------------------------------------
+    Superclass::GenerateOutputInformation();
+
+    //--------------------------------------------------
+    // Set length of vectors in the output vector images
+    //--------------------------------------------------
+    for (unsigned int i = 0;i < m_NumberOfAnisotropicCompartments;++i)
+    {
+        OutputImageType *output = this->GetOutput(i);
+        output->SetVectorLength(4); // AST: now fixed for Watson
+    }
+
+    OutputImageType *output = this->GetOutput(m_NumberOfAnisotropicCompartments);
+    output->SetVectorLength(m_NumberOfAnisotropicCompartments);
+}
+
+template <class TPixelType>
+void
+MCMOrientationPriorsImageFilter <TPixelType>
+::BeforeThreadedGenerateData()
+{
+    this->Superclass::BeforeThreadedGenerateData();
+
+    m_ReferenceInputModels.resize(this->GetNumberOfIndexedInputs());
+
+    InputImageType *input = const_cast <InputImageType *> (this->GetInput(0));
+    MCModelPointer tmpMCM = input->GetDescriptionModel();
+    unsigned int numIsoCompartments = tmpMCM->GetNumberOfIsotropicCompartments();
+    unsigned int numCompartments = tmpMCM->GetNumberOfCompartments();
+
+    for (unsigned int i = 0;i < this->GetNumberOfIndexedInputs();++i)
+    {
+        InputImageType *input = const_cast <InputImageType *> (this->GetInput(i));
+        MCModelPointer model = input->GetDescriptionModel();
+
+        if (model->GetNumberOfIsotropicCompartments() != numIsoCompartments || model->GetNumberOfCompartments() != numCompartments)
+            itkExceptionMacro("All input MCM images should have the same compartment structure.");
+
+        m_ReferenceInputModels[i] = model;
+    }
+}
+
+template <class TPixelType>
+void
+MCMOrientationPriorsImageFilter <TPixelType>
+::DynamicThreadedGenerateData(const InputRegionType &region)
+{
+    unsigned int numInputs = this->GetNumberOfIndexedInputs();
+    unsigned int numOutputs = this->GetNumberOfIndexedOutputs();
+
+    using InputIteratorType = itk::ImageRegionConstIterator <InputImageType>;
+    std::vector<InputIteratorType> inputIterators(numInputs);
+    for (unsigned int i = 0; i < numInputs; ++i)
+        inputIterators[i] = InputIteratorType(this->GetInput(i), region);
+
+    using OutputIteratorType = itk::ImageRegionIterator <OutputImageType>;
+    std::vector<OutputIteratorType> outputIterators(numOutputs);
+    for (unsigned int i = 0; i < numOutputs; ++i)
+        outputIterators[i] = OutputIteratorType(this->GetOutput(i), region);
+
+    using MaskIteratorType = itk::ImageRegionConstIterator <MaskImageType>;
+    std::vector<MaskIteratorType> maskIterators;
+    if (m_MaskImages.size() == numInputs)
+    {
+        for (unsigned int i = 0; i < numInputs; ++i)
+            maskIterators.push_back(MaskIteratorType(m_MaskImages[i], region));
+    }
+    unsigned int numMasks = maskIterators.size();
+
+    std::vector<MCModelPointer> workInputModels(numInputs);
+    for (unsigned int i = 0;i < numInputs;++i)
+        workInputModels[i] = m_ReferenceInputModels[i]->Clone();
+    
+    using ClusterFilterType = anima::SpectralClusteringFilter<double>;
+    using OrientationVectorType = itk::Vector<double,3>;
+
+    ClusterFilterType clusterFilter;
+    ClusterFilterType::MatrixType squaredDistanceMatrix;
+    std::vector<OrientationVectorType> workAllOrientations(numInputs * m_NumberOfAnisotropicCompartments);
+    std::vector<OrientationVectorType> workInputOrientations;
+    std::vector<double> workAllWeights(numInputs * m_NumberOfAnisotropicCompartments);
+    vnl_matrix<double> workInputWeights;
+    std::vector<double> workConcentrationParameters;
+    std::vector<bool> usefulInputsForWeights(numInputs, false);
+    std::vector<unsigned int> workAllMemberships(numInputs * m_NumberOfAnisotropicCompartments);
+    OrientationVectorType workSphericalOrientation, workCartesianOrientation;
+    workSphericalOrientation[2] = 1.0;
+    std::vector<unsigned int> clusterMembers;
+    anima::WatsonDistribution watsonDistribution;
+    anima::DirichletDistribution dirichletDistribution;
+
+    std::vector<OutputPixelType> outputOrientationValues(m_NumberOfAnisotropicCompartments);
+    for (unsigned int i = 0;i < m_NumberOfAnisotropicCompartments;++i)
+        outputOrientationValues[i].SetSize(4);
+    OutputPixelType outputWeightValues;
+    outputWeightValues.SetSize(m_NumberOfAnisotropicCompartments);
+
+    while (!inputIterators[0].IsAtEnd())
+    {
+        unsigned int numUsefulInputs = 0;
+        for (unsigned int i = 0;i < numInputs;++i)
+        {
+            if (isZero(inputIterators[i].Get()))
+                continue;
+
+            if (numMasks == numInputs)
+            {
+                if (maskIterators[i].Get() == 0)
+                    continue;
+            }
+
+            workInputModels[i]->SetModelVector(inputIterators[i].Get());
+
+            ++numUsefulInputs;
+        }
+
+        if (numUsefulInputs == 0)
+        {
+            for (unsigned int i = 0;i < m_NumberOfAnisotropicCompartments;++i)
+            {
+                outputOrientationValues[i].Fill(0.0);
+                outputIterators[i].Set(outputOrientationValues[i]);
+                ++outputIterators[i];
+            }
+
+            outputWeightValues.Fill(0.0);
+            outputIterators[m_NumberOfAnisotropicCompartments].Set(outputWeightValues);
+            ++outputIterators[m_NumberOfAnisotropicCompartments];
+
+            for (unsigned int i = 0;i < numInputs;++i)
+                ++inputIterators[i];
+            
+            if (numMasks == numInputs)
+            {
+                for (unsigned int i = 0;i < numMasks;++i)
+                    ++maskIterators[i];
+            }
+
+            continue;
+        }
+
+        // Grab list of all orientations, memberships and weights from all MCM models
+        for (unsigned int i = 0;i < numInputs;++i)
+        {
+            unsigned int numIsoCompartments = workInputModels[i]->GetNumberOfIsotropicCompartments();
+            
+            double sumOfAnisotropicWeights = 0.0;
+            for (unsigned int j = 0;j < m_NumberOfAnisotropicCompartments;++j)
+            {
+                workSphericalOrientation[0] = workInputModels[i]->GetCompartment(j + numIsoCompartments)->GetOrientationTheta();
+                workSphericalOrientation[1] = workInputModels[i]->GetCompartment(j + numIsoCompartments)->GetOrientationPhi();
+                anima::TransformSphericalToCartesianCoordinates(workSphericalOrientation, workCartesianOrientation);
+                workAllOrientations[i * m_NumberOfAnisotropicCompartments + j] = workCartesianOrientation;
+                workAllMemberships[i * m_NumberOfAnisotropicCompartments + j] = i;
+                sumOfAnisotropicWeights += workInputModels[i]->GetCompartmentWeight(j + numIsoCompartments);
+            }
+
+            for (unsigned int j = 0;j < m_NumberOfAnisotropicCompartments;++j)
+                workAllWeights[i * m_NumberOfAnisotropicCompartments + j] = workInputModels[i]->GetCompartmentWeight(j + numIsoCompartments) / sumOfAnisotropicWeights;
+        }
+
+        // Compute squared distance matrix between orientations
+        unsigned int numOrientations = workAllOrientations.size();
+        squaredDistanceMatrix.set_size(numOrientations, numOrientations);
+        squaredDistanceMatrix.fill_diagonal(0.0);
+        for (unsigned int i = 0;i < numOrientations - 1;++i)
+        {
+            workCartesianOrientation = workAllOrientations[i];
+            for (unsigned int j = i + 1;j < numOrientations;++j)
+            {
+                double cosValue = anima::ComputeScalarProduct(workAllOrientations[j], workCartesianOrientation);
+                cosValue = std::abs(cosValue);
+                if (cosValue > 1.0)
+                    cosValue = 1.0;
+                double distanceValue = std::acos(cosValue) / M_PI;
+                squaredDistanceMatrix.put(i, j, distanceValue * distanceValue);
+                if (j != i)
+                    squaredDistanceMatrix.put(j, i, distanceValue * distanceValue);
+            }
+        }
+
+        // Perform clustering based on orientations
+        clusterFilter.SetInputData(squaredDistanceMatrix);
+        clusterFilter.SetNbClass(m_NumberOfAnisotropicCompartments);
+        clusterFilter.SetVerbose(false);
+        clusterFilter.Update();
+
+        // Characterize orientations in each cluster by a Watson distribution
+        for (unsigned int i = 0;i < m_NumberOfAnisotropicCompartments;++i)
+        {
+            clusterMembers = clusterFilter.GetClassMembers(i);
+            unsigned int clusterSize = clusterMembers.size();
+            workInputOrientations.resize(clusterSize);
+            for (unsigned int j = 0;j < clusterSize;++j)
+                workInputOrientations[j] = workAllOrientations[clusterMembers[j]];
+                
+            watsonDistribution.Fit(workInputOrientations, "mle");
+            workCartesianOrientation = watsonDistribution.GetMeanAxis();
+            for (unsigned int j = 0;j < 3;++j)
+                outputOrientationValues[i][j] = workCartesianOrientation[j];
+            outputOrientationValues[i][3] = watsonDistribution.GetConcentrationParameter();
+            outputIterators[i].Set(outputOrientationValues[i]);
+            ++outputIterators[i];
+        }
+
+        // Characterize anisotropic weights in each cluster by a Dirichlet distribution
+
+        // First, find MCMs whose anisotropic compartments did spread in all clusters
+        std::fill(usefulInputsForWeights.begin(), usefulInputsForWeights.end(), false);
+        unsigned int nbUsefulInputsForWeights = 0;
+        for (unsigned int i = 0;i < numInputs;++i)
+        {
+            unsigned int numCompartments = 0;
+            for (unsigned int j = 0;j < m_NumberOfAnisotropicCompartments;++j)
+            {
+                clusterMembers = clusterFilter.GetClassMembers(j);
+                for (unsigned int k = 0;k < clusterMembers.size();++k)
+                {
+                    if (workAllMemberships[clusterMembers[k]] == i)
+                    {
+                        numCompartments++;
+                        break;
+                    }
+                }
+            }
+
+            if (numCompartments == m_NumberOfAnisotropicCompartments)
+            {
+                usefulInputsForWeights[i] = true;
+                nbUsefulInputsForWeights++;
+            }
+        }
+
+        // Next, compute Dirichlet priors in each cluster
+        if (nbUsefulInputsForWeights < 2)
+        {
+            outputWeightValues.Fill(1.0);
+            outputIterators[m_NumberOfAnisotropicCompartments].Set(outputWeightValues);
+            ++outputIterators[m_NumberOfAnisotropicCompartments];
+
+            for (unsigned int i = 0;i < numInputs;++i)
+                ++inputIterators[i];
+            
+            if (numMasks == numInputs)
+            {
+                for (unsigned int i = 0;i < numMasks;++i)
+                    ++maskIterators[i];
+            }
+
+            continue;
+        }
+        
+        workInputWeights.set_size(nbUsefulInputsForWeights, m_NumberOfAnisotropicCompartments);
+        unsigned int pos = 0;
+        for (unsigned int i = 0;i < numInputs;++i)
+        {
+            if (!usefulInputsForWeights[i])
+                continue;
+            
+            for (unsigned int j = 0;j < m_NumberOfAnisotropicCompartments;++j)
+            {
+                clusterMembers = clusterFilter.GetClassMembers(j);
+                for (unsigned int k = 0;k < clusterMembers.size();++k)
+                {
+                    if (workAllMemberships[clusterMembers[k]] == i)
+                    {
+                        workInputWeights.put(pos, j, workAllWeights[clusterMembers[k]]);
+                        break;
+                    }
+                }
+            }
+            
+            pos++;
+        }
+
+        dirichletDistribution.Fit(workInputWeights, "mle");
+        workConcentrationParameters = dirichletDistribution.GetConcentrationParameters();
+        for (unsigned int i = 0;i < m_NumberOfAnisotropicCompartments;++i)
+            outputWeightValues[i] = workConcentrationParameters[i];
+        outputIterators[m_NumberOfAnisotropicCompartments].Set(outputWeightValues);
+        ++outputIterators[m_NumberOfAnisotropicCompartments];
+
+        for (unsigned int i = 0;i < numInputs;++i)
+            ++inputIterators[i];
+        
+        if (numMasks == numInputs)
+        {
+            for (unsigned int i = 0;i < numMasks;++i)
+                ++maskIterators[i];
+        }
+    }
+}
+
+} // end namespace anima