- //===========================================
- // Create image
- gdcm::File * header = clitk::readDicomHeader(input_files[sliceIndex.front()]);
- typedef itk::Image<float, 3> ImageType;
- typedef itk::Image<signed short, 3> ImageTypeOutput;
- ImageType::SpacingType spacing;
- ImageType::SizeType size;
- ImageType::PointType origin;
- spacing[0] = header->GetXSpacing();
- spacing[1] = header->GetYSpacing();
- spacing[2] = dist;
- size[0] = header->GetXSize();
- size[1] = header->GetYSize();
- size[2] = sliceIndex.size();
- ///Optional use of special origin scheme used at CLB
- if (args_info.focal_origin_flag)
- {
- origin[0] = -spacing[0]*size[0]/2;
- origin[1] = -spacing[1]*size[1]/2;
- }
- else
- {
- origin[0] = header->GetXOrigin();
- origin[1] = header->GetYOrigin();
- }
- origin[2] = header->GetZOrigin();
- ImageTypeOutput::Pointer output = ImageTypeOutput::New();
- itk::ImageRegion<3> region;
- region.SetSize(size);
- output->SetRegions(region);
- output->SetOrigin(origin);
- output->SetSpacing(spacing);
- output->Allocate();
-
- //===========================================
- // Fill image
- ImageType::Pointer slice;
- typedef itk::ImageRegionIterator<ImageTypeOutput> IteratorType;
- typedef itk::ImageRegionConstIterator<ImageType> ConstIteratorType;
- IteratorType po(output, output->GetLargestPossibleRegion());
- po.Begin();
- for(unsigned int i=0; i<sliceIndex.size(); i++) {
- slice = clitk::readImage<ImageType>(input_files[sliceIndex[i]]);
- ConstIteratorType pi(slice, slice->GetLargestPossibleRegion());
- while ( !pi.IsAtEnd() ) {
- po.Set((signed short)pi.Get());
- ++po;
- ++pi;
+ //===========================================
+ // Analyze slices locations ...
+ double currentDist;
+ double dist=0;
+ double tolerance = args_info.tolerance_arg;
+ double previous = sliceLocations[sliceIndex[0]];
+ for(unsigned int i=1; i<sliceIndex.size(); i++) {
+ currentDist = sliceLocations[sliceIndex[i]]-previous;
+ if (i!=1) {
+ if (fabs(dist-currentDist) > tolerance) {
+ std::cout << "ERROR : " << std::endl
+ << "Current slice pos is = " << sliceLocations[sliceIndex[i]] << std::endl
+ << "Previous slice pos is = " << previous << std::endl
+ << "Current file is = " << input_files[sliceIndex[i]] << std::endl
+ << "Current index is = " << i << std::endl
+ << "Current sortindex is = " << sliceIndex[i] << std::endl
+ << "Current slice diff = " << dist << std::endl
+ << "Current error = " << fabs(dist-currentDist) << std::endl;
+ exit(1);