[Insight-users] Image Intensities Changed During DICOM Read-Write-Read in ITK 3.20 vs. ITK 4.4.1

Constantine Z mnkz at leeds.ac.uk
Mon Aug 12 14:45:27 EDT 2013


Thanks. Yes, I can see it all looks completely legit.

Perhaps I can send you my data? MR from a Philips scanner. Or maybe on 
Win it might work differently?

Regards,
Constartine

On 12/08/2013 7:38 p.m., Bill Lorensen wrote:
>
> I just tried your code on two different dicom datasets, a ct and an mr.
>
> Both worked for me.
> ./DICOMTest ~/Data/itkTestcases/Back\ to\ Front/ .
> ITK Version: 4.5
> Stats for file: /home/lorensen/Data/itkTestcases/Back to
> Front//i24389.MRDC.34
> Min: 0, Max: 293, Mean: 28.7654, Sigma: 37.2513, Variance: 1387.66, Sum:
> 1.88517e+06
> Stats for file: ./file00.dcm
> Min: 0, Max: 293, Mean: 28.7654, Sigma: 37.2513, Variance: 1387.66, Sum:
> 1.88517e+06
> [build] ./DICOMTest ~/Data/NAV/CT/ ./
> ITK Version: 4.5
> Stats for file: /home/lorensen/Data/NAV/CT//ct.1
> Min: -1024, Max: 31744, Mean: 6277.24, Sigma: 13249.9, Variance:
> 1.75559e+08, Sum: 1.64554e+09
> Stats for file: .//file00.dcm
> Min: -1024, Max: 31744, Mean: 6277.24, Sigma: 13249.9, Variance:
> 1.75559e+08, Sum: 1.64554e+09
>
>
>
>
> On Mon, Aug 12, 2013 at 11:46 AM, Constantine Zakkaroff
> <mnkz at leeds.ac.uk <mailto:mnkz at leeds.ac.uk>> wrote:
>
>     Hi Brad,
>
>     Yes, here's the minimal code:
>
>     ******************************__******************************__************
>     #include <cstdlib>
>     #include <itkImage.h>
>     #include <itkMetaImageIO.h>
>     #include <itkGDCMImageIO.h>
>     #include <itkGDCMSeriesFileNames.h>
>     #include <itkImageFileReader.h>
>     #include <itkImageFileWriter.h>
>     #include <itkStatisticsImageFilter.h>
>
>     int main(int argc, char *argv[])
>     {
>              std::clog << "ITK Version: " << ITK_VERSION_STRING <<
>     std::endl;
>
>              typedef itk::GDCMSeriesFileNames NamesGeneratorType;
>              typedef itk::Image<short, 3> ImageType;
>
>              NamesGeneratorType::Pointer namesGenerator =
>     NamesGeneratorType::New();
>              namesGenerator->__SetUseSeriesDetails(true);
>              namesGenerator->__SetInputDirectory(argv[1]);
>
>              typedef itk::StatisticsImageFilter<__ImageType>
>     StatsFilterType;
>              typedef itk::ImageFileReader<__ImageType> ImageReaderType;
>              typedef itk::ImageFileWriter<__ImageType> ImageWriterType;
>
>              ImageReaderType::Pointer imageReader0 = ImageReaderType::New();
>              imageReader0->SetImageIO(itk::__GDCMImageIO::New());
>              StatsFilterType::Pointer statsFilter0 = StatsFilterType::New();
>
>     imageReader0->SetFileName(__namesGenerator->__GetInputFileNames()[0]);
>              statsFilter0->SetInput(__imageReader0->GetOutput());
>              statsFilter0->Update();
>
>              std::clog << "Stats for file: " <<
>     namesGenerator->__GetInputFileNames()[0] << std::endl;
>              std::clog <<
>                      "Min: " << statsFilter0->GetMinimum() << ", "
>                      "Max: " << statsFilter0->GetMaximum() << ", "
>                      "Mean: " << statsFilter0->GetMean() << ", "
>                      "Sigma: " << statsFilter0->GetSigma() << ", "
>                      "Variance: " << statsFilter0->GetVariance() << ", "
>                      "Sum: " << statsFilter0->GetSum()
>                      << std::endl;
>
>              ImageWriterType::Pointer imageWriter0 = ImageWriterType::New();
>              imageWriter0->SetInput(__imageReader0->GetOutput());
>              imageWriter0->SetImageIO(itk::__GDCMImageIO::New());
>              std::string fileName0 = std::string(argv[2]) + "/file00.dcm";
>              imageWriter0->SetFileName(__fileName0);
>              imageWriter0->Update();
>
>              ImageReaderType::Pointer imageReader1 = ImageReaderType::New();
>              imageReader1->SetImageIO(itk::__GDCMImageIO::New());
>              StatsFilterType::Pointer statsFilter1 = StatsFilterType::New();
>              imageReader1->SetFileName(__fileName0);
>              statsFilter1->SetInput(__imageReader1->GetOutput());
>              statsFilter1->Update();
>
>              std::clog << "Stats for file: " << fileName0 << std::endl;
>              std::clog <<
>                      "Min: " << statsFilter1->GetMinimum() << ", "
>                      "Max: " << statsFilter1->GetMaximum() << ", "
>                      "Mean: " << statsFilter1->GetMean() << ", "
>                      "Sigma: " << statsFilter1->GetSigma() << ", "
>                      "Variance: " << statsFilter1->GetVariance() << ", "
>                      "Sum: " << statsFilter1->GetSum()
>                      << std::endl;
>
>              system("pause");
>
>              return EXIT_SUCCESS;
>     }
>     ******************************__******************************__************
>
>     Many thanks,
>     Constantine
>
>
>     On 12/08/2013 14:48, Bradley Lowekamp wrote:
>
>
>         Do you have the code to share so that it can be easily reproduced?
>
>
>         Thanks,
>         Brad
>
>         On Aug 12, 2013, at 9:33 AM, Constantine Zakkaroff
>         <mnkz at leeds.ac.uk <mailto:mnkz at leeds.ac.uk>> wrote:
>
>             Hello ALL,
>
>             I'm wondering if anyone can explain what's going on with the
>             changing image intensities in a read-write-read round-trip.
>             I think it's a bug, but the issues.itk.org
>             <http://issues.itk.org> is unresponsive at the moment.
>             Another bug? ;)
>
>             Anyhow, the strange thing is that the results for the round
>             trip read-write-read test are very different for ITK 3.20
>             and ITK 4.4.1.
>
>             I use to StatisticsImageFilter to examine the intensity
>             values. Here's what happens, when I simply read a DICOM
>             image, print stats, save the image and read it again,
>             followed by printing the stats:
>
>             ITK Version: 3.20
>             Stats for original file:
>             Min: 0, Max: 2474, Mean: 169.236, Sigma: 256.719, Variance:
>             65904.7, Sum: 1.40371e+007
>             Stats for saved file:
>             Min: 0, Max: 2471, Mean: 167.132, Sigma: 256.125, Variance:
>             65600.1, Sum: 1.38626e+007
>
>             ITK Version: 4.4
>             Stats for original file:
>             Min: 0, Max: 2474, Mean: 169.236, Sigma: 256.719, Variance:
>             65904.7, Sum: 1.40371e+007
>             Stats for saved file:
>             Min: 0, Max: 7446, Mean: 509.114, Sigma: 772.59, Variance:
>             596895, Sum: 4.22279e+007
>
>             In both cases the code is identical. Just compiled against
>             different versions of ITK. The images are read and written
>             as itk::Image<unsigned short, 3>, with no processing, no
>             changes to any aspect of the images.
>
>             So, ITK 3.20 introduces a small change, which I'm pretty
>             sure shouldn't be there. But in ITK 4.4.1 the scale of
>             intensities changes three-fold.
>
>             I'd not be surprised by a change of this magnitude if I was
>             photocopying those images, but not when using digital
>             software and hardware. :^) Is there something I'm not
>             understanding here?
>
>             None of this happens with
>             read(DICOM)-write(MHD)-read(__MHD), where the original and
>             round-trip values are identical:
>             Stats for saved file: file0.mhd
>             Min: 0, Max: 2474, Mean: 169.236, Sigma: 256.719, Variance:
>             65904.7, Sum: 1.40371e+007
>
>             Many thanks,
>             Constantine
>
>             _______________________________________
>             Powered by www.kitware.com <http://www.kitware.com>
>
>             Visit other Kitware open-source projects at
>             http://www.kitware.com/__opensource/opensource.html
>             <http://www.kitware.com/opensource/opensource.html>
>
>             Kitware offers ITK Training Courses, for more information visit:
>             http://www.kitware.com/__products/protraining.php
>             <http://www.kitware.com/products/protraining.php>
>
>             Please keep messages on-topic and check the ITK FAQ at:
>             http://www.itk.org/Wiki/ITK___FAQ
>             <http://www.itk.org/Wiki/ITK_FAQ>
>
>             Follow this link to subscribe/unsubscribe:
>             http://www.itk.org/mailman/__listinfo/insight-users
>             <http://www.itk.org/mailman/listinfo/insight-users>
>
>
>     _______________________________________
>     Powered by www.kitware.com <http://www.kitware.com>
>
>     Visit other Kitware open-source projects at
>     http://www.kitware.com/__opensource/opensource.html
>     <http://www.kitware.com/opensource/opensource.html>
>
>     Kitware offers ITK Training Courses, for more information visit:
>     http://www.kitware.com/__products/protraining.php
>     <http://www.kitware.com/products/protraining.php>
>
>     Please keep messages on-topic and check the ITK FAQ at:
>     http://www.itk.org/Wiki/ITK___FAQ <http://www.itk.org/Wiki/ITK_FAQ>
>
>     Follow this link to subscribe/unsubscribe:
>     http://www.itk.org/mailman/__listinfo/insight-users
>     <http://www.itk.org/mailman/listinfo/insight-users>
>
>
>
>
> --
> Unpaid intern in BillsBasement at noware dot com



More information about the Insight-users mailing list