--- /dev/null
+/*=========================================================================
+ Program: vv http://www.creatis.insa-lyon.fr/rio/vv
+
+ Authors belong to:
+ - University of LYON http://www.universite-lyon.fr/
+ - Léon Bérard cancer center http://www.centreleonberard.fr
+ - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the copyright notices for more information.
+
+ It is distributed under dual licence
+
+ - BSD See included LICENSE.txt file
+ - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+===========================================================================**/
+#ifndef CLITKSCINTIVOLSTATS_CXX
+#define CLITKSCINTIVOLSTATS_CXX
+
+// clitk include
+#include "clitkScintivolStats_ggo.h"
+#include "clitkIO.h"
+#include "clitkImageCommon.h"
+#include "clitkCommon.h"
+
+#include "itkImageFileReader.h"
+#include "itkBinaryImageToLabelMapFilter.h"
+#include "itkLabelMapToLabelImageFilter.h"
+#include "itkLabelStatisticsImageFilter.h"
+#include "itkStatisticsImageFilter.h"
+#include "itkExtractImageFilter.h"
+
+#include <iostream>
+
+//-------------------------------------------------------------------=
+int main(int argc, char * argv[])
+{
+
+ // init command line
+ GGO(clitkScintivolStats, args_info);
+ CLITK_INIT;
+
+ typedef itk::Image<double, 3> Input3DType;
+ typedef itk::Image<double, 2> Input2DType;
+ typedef itk::Image<char, 3> Input3DMaskType;
+ typedef itk::Image<char, 2> Input2DMaskType;
+ typedef itk::ImageFileReader<Input3DType> InputReader3DType;
+ typedef itk::ImageFileReader<Input3DMaskType> InputReader3DMaskType;
+ typedef itk::ImageFileReader<Input2DMaskType> InputReader2DMaskType;
+ typedef itk::LabelStatisticsImageFilter<Input3DType, Input3DMaskType> LabelStatistics3DImageFilterType;
+ typedef itk::LabelStatisticsImageFilter<Input2DType, Input2DMaskType> LabelStatistics2DImageFilterType;
+ typedef itk::StatisticsImageFilter<Input2DType> Statistics2DImageFilterType;
+ typedef itk::ExtractImageFilter<Input3DType, Input2DType> ExtractImageFilter;
+
+ //Determine nbAcquisitionDynamic1 and nbAcquisitionDynamic2 (ie. the number of slice in dynamic1 and dynamic2)
+ //Start to open the dynamic1 and dynamic2 and get the 3rd dimension size
+ InputReader3DType::Pointer readerDynamic1 = InputReader3DType::New();
+ readerDynamic1->SetFileName(args_info.dynamic1_arg);
+ readerDynamic1->Update();
+ Input3DType::Pointer dynamic1 = readerDynamic1->GetOutput();
+ int nbAcquisitionDynamic1 = dynamic1->GetLargestPossibleRegion().GetSize()[2];
+
+ InputReader3DType::Pointer readerDynamic2 = InputReader3DType::New();
+ readerDynamic2->SetFileName(args_info.dynamic2_arg);
+ readerDynamic2->Update();
+ Input3DType::Pointer dynamic2 = readerDynamic2->GetOutput();
+ int nbAcquisitionDynamic2 = dynamic2->GetLargestPossibleRegion().GetSize()[2];
+
+ //Open the csv file
+ //If it's scatter or attenuation correction, just append results to the csv file
+ //If not open it normally and write the number of dynamic acquisition for the first and the second images (given by the ggo and from the dicom tag)
+ std::ofstream csvFile;
+ if (args_info.append_flag)
+ csvFile.open (args_info.output_arg, std::ios::app);
+ else {
+ csvFile.open (args_info.output_arg);
+ csvFile << "1;" << nbAcquisitionDynamic1 << ";" << nbAcquisitionDynamic2 << "\n";
+ }
+
+ //Read Tomo image, total and remnant masks
+ InputReader3DType::Pointer readerTomo = InputReader3DType::New();
+ readerTomo->SetFileName(args_info.tomo_arg);
+ readerTomo->Update();
+ Input3DType::Pointer tomo = readerTomo->GetOutput();
+
+ InputReader3DMaskType::Pointer readerTotalLiver = InputReader3DMaskType::New();
+ readerTotalLiver->SetFileName(args_info.totalLiverMask_arg);
+ readerTotalLiver->Update();
+
+ InputReader3DMaskType::Pointer readerRemnantLiver = InputReader3DMaskType::New();
+ readerRemnantLiver->SetFileName(args_info.remnantLiverMask_arg);
+ readerRemnantLiver->Update();
+
+ //Find number of counts in tomo for total Liver and remnant Liver
+ LabelStatistics3DImageFilterType::Pointer labelStatisticsImageFilterTotalLiver = LabelStatistics3DImageFilterType::New();
+ labelStatisticsImageFilterTotalLiver->SetLabelInput(readerTotalLiver->GetOutput());
+ labelStatisticsImageFilterTotalLiver->SetInput(tomo);
+ labelStatisticsImageFilterTotalLiver->Update();
+
+ LabelStatistics3DImageFilterType::Pointer labelStatisticsImageFilterRemnantLiver = LabelStatistics3DImageFilterType::New();
+ labelStatisticsImageFilterRemnantLiver->SetLabelInput(readerRemnantLiver->GetOutput());
+ labelStatisticsImageFilterRemnantLiver->SetInput(tomo);
+ labelStatisticsImageFilterRemnantLiver->Update();
+
+ //Write them in the csv file
+ csvFile << args_info.acquisitionTimeTomo_arg << ";" << labelStatisticsImageFilterTotalLiver->GetSum(1) << ";" << labelStatisticsImageFilterRemnantLiver->GetSum(1) << "\n";
+
+ //Read liver and heart masks for dynamic1
+ InputReader2DMaskType::Pointer readerLiver = InputReader2DMaskType::New();
+ readerLiver->SetFileName(args_info.liverMask_arg);
+ readerLiver->Update();
+
+ InputReader2DMaskType::Pointer readerHeart = InputReader2DMaskType::New();
+ readerHeart->SetFileName(args_info.heartMask_arg);
+ readerHeart->Update();
+
+ for (unsigned int i=0; i<nbAcquisitionDynamic1; ++i) {
+ //Compute time frame in s
+ double timeFrame = args_info.frameDurationDynamic1_arg * i + args_info.frameDurationDynamic1_arg;
+
+ //Extract the Region of interest of the dynamic1
+ ExtractImageFilter::Pointer extractSlice = ExtractImageFilter::New();
+ Input3DType::IndexType start;
+ start[0] = 0;
+ start[1] = 0;
+ start[2] = i;
+ Input3DType::SizeType size;
+ size[0] = dynamic1->GetLargestPossibleRegion().GetSize()[0];
+ size[1] = dynamic1->GetLargestPossibleRegion().GetSize()[1];
+ size[2] = 0;
+ Input3DType::RegionType desiredRegion;
+ desiredRegion.SetSize(size);
+ desiredRegion.SetIndex(start);
+ extractSlice->SetExtractionRegion(desiredRegion);
+ extractSlice->SetInput(dynamic1);
+#if ITK_VERSION_MAJOR >= 4
+ extractSlice->SetDirectionCollapseToIdentity();
+#endif
+ extractSlice->Update();
+
+ //Find number of counts in dynamic1 slice for Liver and heart
+ LabelStatistics2DImageFilterType::Pointer labelStatisticsImageFilterLiver = LabelStatistics2DImageFilterType::New();
+ labelStatisticsImageFilterLiver->SetLabelInput(readerLiver->GetOutput());
+ labelStatisticsImageFilterLiver->SetInput(extractSlice->GetOutput());
+ labelStatisticsImageFilterLiver->SetCoordinateTolerance(0.001);
+ labelStatisticsImageFilterLiver->SetDirectionTolerance(0.001);
+ labelStatisticsImageFilterLiver->Update();
+
+ LabelStatistics2DImageFilterType::Pointer labelStatisticsImageFilterHeart = LabelStatistics2DImageFilterType::New();
+ labelStatisticsImageFilterHeart->SetLabelInput(readerHeart->GetOutput());
+ labelStatisticsImageFilterHeart->SetInput(extractSlice->GetOutput());
+ labelStatisticsImageFilterHeart->SetCoordinateTolerance(0.001);
+ labelStatisticsImageFilterHeart->SetDirectionTolerance(0.001);
+ labelStatisticsImageFilterHeart->Update();
+
+ //Find number of counts in dynamic1 slice
+ Statistics2DImageFilterType::Pointer statisticsImageFilterDynamic1 = Statistics2DImageFilterType::New();
+ statisticsImageFilterDynamic1->SetInput(extractSlice->GetOutput());
+ statisticsImageFilterDynamic1->Update();
+
+ //Write them in the csv file
+ csvFile << timeFrame << ";" << labelStatisticsImageFilterLiver->GetSum(1) << ";" << labelStatisticsImageFilterHeart->GetSum(1) << ";" << statisticsImageFilterDynamic1->GetSum() << "\n";
+ }
+
+ //Read parenchyma mask for dynamic2
+ InputReader2DMaskType::Pointer readerParenchyma = InputReader2DMaskType::New();
+ readerParenchyma->SetFileName(args_info.parenchymaMask_arg);
+ readerParenchyma->Update();
+
+ for (unsigned int i=0; i<nbAcquisitionDynamic2; ++i) {
+ //Compute time frame in s
+ double timeFrame = args_info.frameDurationDynamic2_arg * i + args_info.acquisitionTimeDynamic2_arg + args_info.frameDurationDynamic2_arg;
+
+ //Extract the Region of interest of the dynamic1
+ ExtractImageFilter::Pointer extractSlice = ExtractImageFilter::New();
+ Input3DType::IndexType start;
+ start[0] = 0;
+ start[1] = 0;
+ start[2] = i;
+ Input3DType::SizeType size;
+ size[0] = dynamic2->GetLargestPossibleRegion().GetSize()[0];
+ size[1] = dynamic2->GetLargestPossibleRegion().GetSize()[1];
+ size[2] = 0;
+ Input3DType::RegionType desiredRegion;
+ desiredRegion.SetSize(size);
+ desiredRegion.SetIndex(start);
+ extractSlice->SetExtractionRegion(desiredRegion);
+ extractSlice->SetInput(dynamic2);
+#if ITK_VERSION_MAJOR >= 4
+ extractSlice->SetDirectionCollapseToIdentity();
+#endif
+ extractSlice->Update();
+
+ //Find number of counts in dynamic2 slice for parenchyma
+ LabelStatistics2DImageFilterType::Pointer labelStatisticsImageFilterParenchyma = LabelStatistics2DImageFilterType::New();
+ labelStatisticsImageFilterParenchyma->SetLabelInput(readerParenchyma->GetOutput());
+ labelStatisticsImageFilterParenchyma->SetInput(extractSlice->GetOutput());
+ labelStatisticsImageFilterParenchyma->SetCoordinateTolerance(0.001);
+ labelStatisticsImageFilterParenchyma->SetDirectionTolerance(0.001);
+ labelStatisticsImageFilterParenchyma->Update();
+
+ //Write them in the csv file
+ csvFile << timeFrame << ";" << labelStatisticsImageFilterParenchyma->GetSum(1) << "\n";
+ }
+
+ csvFile.close();
+
+ // this is the end my friend
+ return 0;
+}
+//-------------------------------------------------------------------=
+
+#endif /* end #define CLITKSCINTIVOLSTATS_CXX */
git://git.creatis.insa-lyon.fr / clitk.git / commitdiff