Skip to content

Commit

Permalink
Merge pull request #106 from astamm/statistical_distributions
Browse files Browse the repository at this point in the history
Statistical distributions
  • Loading branch information
astamm authored Nov 15, 2023
2 parents dd203c5 + c04ba06 commit fc05522
Show file tree
Hide file tree
Showing 62 changed files with 4,426 additions and 3,139 deletions.
2 changes: 1 addition & 1 deletion Anima/diffusion/mcm_tools/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -7,4 +7,4 @@ 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)
add_subdirectory(test_averaging)
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@
#include <animaMCMImage.h>
#include <animaMultiCompartmentModel.h>
#include <animaReadWriteFunctions.h>

#include <itkImageRegionConstIterator.h>
#include <itkImageRegionIterator.h>

Expand All @@ -10,41 +11,47 @@
int main(int argc, char **argv)
{
TCLAP::CmdLine cmd("INRIA / IRISA - VisAGeS/Empenn Team", ' ', ANIMA_VERSION);

TCLAP::ValueArg<std::string> inArg("i", "inputfile", "Input MCM image", true, "", "input mcm image", cmd);
TCLAP::ValueArg<std::string> outArg("o", "outputfile", "Output image with combined anisotropic weights", true, "", "output conbined weight image", cmd);


TCLAP::ValueArg<std::string> inArg(
"i", "inputfile",
"A string specifying the name of the file that stores the input MCM image.",
true, "", "input mcm image", cmd);
TCLAP::ValueArg<std::string> outArg(
"o", "outputfile",
"A string specifying the name of the file that stores the output image with combined anisotropic weights.",
true, "", "output combined weight image", cmd);

try
{
cmd.parse(argc,argv);
cmd.parse(argc, argv);
}
catch (TCLAP::ArgException& e)
catch (TCLAP::ArgException &e)
{
std::cerr << "Error: " << e.error() << "for argument " << e.argId() << std::endl;
return EXIT_FAILURE;
}

using InputImageType = anima::MCMImage <double, 3>;
using OutputImageType = itk::VectorImage <double, 3>;
using InputImageType = anima::MCMImage<double, 3>;
using OutputImageType = itk::VectorImage<double, 3>;
using InputPixelType = InputImageType::PixelType;
using OutputPixelType = OutputImageType::PixelType;
using MCModelType = anima::MultiCompartmentModel;
using MCModelPointer = MCModelType::Pointer;

// Read input sample
anima::MCMFileReader <double, 3> mcmReader;
anima::MCMFileReader<double, 3> mcmReader;
mcmReader.SetFileName(inArg.getValue());
mcmReader.Update();
InputImageType::Pointer mcmImage = mcmReader.GetModelVectorImage();
MCModelPointer mcmPtr = mcmImage->GetDescriptionModel()->Clone();
InputPixelType mcmValue;
unsigned int numCompartments = mcmPtr->GetNumberOfCompartments();
unsigned int numIsoCompartments = mcmPtr->GetNumberOfIsotropicCompartments();
using InputImageIteratorType = itk::ImageRegionConstIterator <InputImageType>;
using OutputImageIteratorType = itk::ImageRegionIterator <OutputImageType>;

using InputImageIteratorType = itk::ImageRegionConstIterator<InputImageType>;
using OutputImageIteratorType = itk::ImageRegionIterator<OutputImageType>;
InputImageIteratorType inItr(mcmImage, mcmImage->GetLargestPossibleRegion());

// Initialize output image
unsigned int nbOutputComponents = numIsoCompartments + 1;
OutputImageType::Pointer outputImage = OutputImageType::New();
Expand All @@ -59,15 +66,15 @@ int main(int argc, char **argv)

std::cout << "- Number of compartments in input image: " << numCompartments << std::endl;
std::cout << "- Number of compartments in output image: " << nbOutputComponents << std::endl;

OutputImageIteratorType outItr(outputImage, outputImage->GetLargestPossibleRegion());

while (!outItr.IsAtEnd())
{
while (!outItr.IsAtEnd())
{
mcmValue = inItr.Get();

bool backgroundVoxel = true;
for (unsigned int i = 0;i < mcmValue.GetNumberOfElements();++i)
for (unsigned int i = 0; i < mcmValue.GetNumberOfElements(); ++i)
{
if (mcmValue[i] != 0.0)
{
Expand All @@ -85,22 +92,22 @@ int main(int argc, char **argv)

mcmPtr->SetModelVector(mcmValue);

for (unsigned int i = 0;i < numIsoCompartments;++i)
for (unsigned int i = 0; i < numIsoCompartments; ++i)
outputValue[i] = mcmPtr->GetCompartmentWeight(i);

double anisoWeight = 0.0;
for (unsigned int i = numIsoCompartments;i < numCompartments;++i)
for (unsigned int i = numIsoCompartments; i < numCompartments; ++i)
anisoWeight += mcmPtr->GetCompartmentWeight(i);

outputValue[numIsoCompartments] = anisoWeight;

outItr.Set(outputValue);

++inItr;
++outItr;
}
anima::writeImage <OutputImageType> (outArg.getValue(), outputImage);
return EXIT_SUCCESS;
}

anima::writeImage<OutputImageType>(outArg.getValue(), outputImage);

return EXIT_SUCCESS;
}
44 changes: 18 additions & 26 deletions Anima/diffusion/mcm_tools/mcm_orientation_priors/CMakeLists.txt
Original file line number Diff line number Diff line change
@@ -1,38 +1,30 @@
if(BUILD_TOOLS)

project(animaMCMOrientationPriors)
project(animaMCMOrientationPriors)

## #############################################################################
## List Sources
## #############################################################################
# ############################################################################
# List Sources
# ############################################################################

list_source_files(${PROJECT_NAME}
${CMAKE_CURRENT_SOURCE_DIR}
)
list_source_files(${PROJECT_NAME} ${CMAKE_CURRENT_SOURCE_DIR})

## #############################################################################
## add executable
## #############################################################################
# ############################################################################
# add executable
# ############################################################################

add_executable(${PROJECT_NAME}
${${PROJECT_NAME}_CFILES}
)
add_executable(${PROJECT_NAME} ${${PROJECT_NAME}_CFILES})

## #############################################################################
## Link
## #############################################################################
# ############################################################################
# Link
# ############################################################################

target_link_libraries(${PROJECT_NAME}
${ITKIO_LIBRARIES}
${TinyXML2_LIBRARY}
AnimaMCM
AnimaStatisticalDistributions
)
target_link_libraries(${PROJECT_NAME} ${ITKIO_LIBRARIES} ${TinyXML2_LIBRARY}
AnimaMCM AnimaStatisticalDistributions)

## #############################################################################
## install
## #############################################################################
# ############################################################################
# install
# ############################################################################

set_exe_install_rules(${PROJECT_NAME})
set_exe_install_rules(${PROJECT_NAME})

endif()
Original file line number Diff line number Diff line change
@@ -1,55 +1,59 @@
#include <animaMCMFileReader.h>
#include <animaReadWriteFunctions.h>
#include <animaMCMOrientationPriorsImageFilter.h>
#include <animaReadWriteFunctions.h>

#include <itkTimeProbe.h>

#include <tclap/CmdLine.h>

int main(int argc, char **argv)
int main(int argc, char **argv)
{
TCLAP::CmdLine cmd("INRIA / IRISA - VisAGeS Team", ' ',ANIMA_VERSION);
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
);
"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
);
"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
);

"w", "output-weights",
"A string specifying the filename for the output vector image that will store fractions prior.",
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
);
"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<std::string> meanOrientationsArg(
"", "mean-orientations",
"A string specifying the basename for the output vector images that will store mean orientations.",
false, "", "output mean orientation images", cmd);
TCLAP::ValueArg<std::string> meanFractionsArg(
"", "mean-fractions",
"A string specifying the filename for the output vector image that will store fractions mean.",
false, "", "output mean fractions image", 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
);
"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);
cmd.parse(argc, argv);
}
catch (TCLAP::ArgException& e)
catch (TCLAP::ArgException &e)
{
std::cerr << "Error: " << e.error() << "for argument " << e.argId() << std::endl;
return EXIT_FAILURE;
}

using MainFilterType = anima::MCMOrientationPriorsImageFilter <double>;
using MainFilterType = anima::MCMOrientationPriorsImageFilter<double>;
MainFilterType::Pointer mainFilter = MainFilterType::New();

using MCMReaderType = anima::MCMFileReader <double,3>;
using MCMReaderType = anima::MCMFileReader<double, 3>;
using MaskImageType = MainFilterType::MaskImageType;
using OutputImageType = MainFilterType::OutputImageType;

Expand All @@ -67,9 +71,9 @@ int main(int argc, char **argv)
while (!inputFile.eof())
{
char tmpStr[2048];
inputFile.getline(tmpStr,2048);
inputFile.getline(tmpStr, 2048);

if (strcmp(tmpStr,"") == 0)
if (strcmp(tmpStr, "") == 0)
continue;

std::cout << "Loading image " << nbOfImages + 1 << ": " << tmpStr << std::endl;
Expand All @@ -92,26 +96,38 @@ int main(int argc, char **argv)
char tmpStr[2048];
while (!masksIn.eof())
{
masksIn.getline(tmpStr,2048);
masksIn.getline(tmpStr, 2048);

if (strcmp(tmpStr,"") == 0)
if (strcmp(tmpStr, "") == 0)
continue;

mainFilter->AddMaskImage(anima::readImage <MaskImageType> (tmpStr));
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)
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>(fileName, mainFilter->GetOutput(i));
}

anima::writeImage <OutputImageType> (outWeightsArg.getValue(), mainFilter->GetOutput(numberOfAnisotropicCompartments));

anima::writeImage<OutputImageType>(outWeightsArg.getValue(), mainFilter->GetOutput(numberOfAnisotropicCompartments));

if (meanOrientationsArg.getValue() != "")
{
for (unsigned int i = 0; i < numberOfAnisotropicCompartments; ++i)
{
std::string fileName = meanOrientationsArg.getValue() + "_" + std::to_string(i) + ".nrrd";
anima::writeImage<OutputImageType>(fileName, mainFilter->GetMeanOrientationImage(i));
}
}

if (meanFractionsArg.getValue() != "")
anima::writeImage<OutputImageType>(meanFractionsArg.getValue(), mainFilter->GetMeanFractionsImage());

return EXIT_SUCCESS;
}
Loading

0 comments on commit fc05522

Please sign in to comment.