C++ Pointer::GetOrigin方法代码示例

bool imageSort(const ImageType::Pointer &im1, const ImageType::Pointer &im2)
{
	ImageType::PointType or1 = im1->GetOrigin();
	ImageType::PointType or2 = im2->GetOrigin();

	ImageType::DirectionType direction = im1->GetDirection();
	itk::Vector<double, 3> dir, v1, v2, normed;
   	for(unsigned int i = 0; i < 3; i++)
   	{
   		dir[i] = direction(i,2);
		v1[i] = or1[i];
		v2[i] = or2[i];
   	}

	std::cout << dir << std::endl;
	double mul1 = v1[0] / dir[0];
	double mul2 = v2[0] / dir[0];

	
	



	return (mul1<mul2);

}

// ------------------------------------------------------------------------
void transformPointToPatient(const ImageType::Pointer &reference,
		const SeriesTransform &series, const std::vector<int> &roiOffset,
		const ImageType::PointType &input, ImageType::PointType &output)
{
	// rotation is easy
	MatrixType rotTemp = reference->GetDirection().GetVnlMatrix();
	MatrixType rotation;
	rotation.set_size(3,3);
	for(unsigned int i = 0; i < 3; i++)
	{
		rotation(i,0) = rotTemp(i,1);
		rotation(i,1) = rotTemp(i,0);
		rotation(i,2) = rotTemp(i,2);
	}

	VectorType refOrig = reference->GetOrigin().GetVnlVector();
	VectorType offset  = vnl_inverse(rotation) * refOrig;

	vnl_vector<double> vnl_output;

	// translate the point
	vnl_output = vnl_inverse(rotation) * input.GetVnlVector();
	vnl_output = vnl_output - offset;
	vnl_output /= reference->GetSpacing()[0];
	for(unsigned int i = 0; i < 3; i++)
	{
		vnl_output[i] -= roiOffset[i];
		vnl_output[i] -= series.translation[i];
	}

	output[0] = vnl_output[0];
	output[1] = vnl_output[1];
	output[2] = vnl_output[2];
}

void setup_snapshot_image(string basename, ImageType::Pointer model)
{
   ImageType::SpacingType spacing = model->GetSpacing();
   spacing[0] /= constants::SnapshotImageZoom;
   spacing[1] /= constants::SnapshotImageZoom;
   spacing[2] /= constants::SnapshotImageZoom;

   g_snapshot_image = ImageType::New();
   g_snapshot_image->SetSpacing( spacing );
   g_snapshot_image->SetOrigin( model->GetOrigin() );

   ImageType::RegionType const& region = model->GetLargestPossibleRegion();
   ImageType::RegionType::SizeType size = region.GetSize();

   // size is in pixels
   ImageType::SizeType doubled_size(size);
   doubled_size[0] *= constants::SnapshotImageZoom;
   doubled_size[1] *= constants::SnapshotImageZoom;
   doubled_size[2] *= constants::SnapshotImageZoom;
   g_snapshot_image->SetRegions( doubled_size );

   g_snapshot_image->Allocate();

   s_snapshot_basename = basename;
}

// ------------------------------------------------------------------------
void buildOutput(const SeriesTransform::List &series, ImageType::Pointer &outputImage,
		ImageType::Pointer &outputLabel, const unsigned int & timestep)
{
	// get the output parameters
	unsigned int slices = series.size();
	unsigned int timeSteps = series.front().images.size();
	ImageType::Pointer ref = series.front().images.front();




	
	ImageType::SpacingType spacing = ref->GetSpacing();
	spacing[2] = series.front().sliceThickness;
	ImageType::DirectionType direction = ref->GetDirection();
	ImageType::PointType origin = ref->GetOrigin();
	ImageType::RegionType region = ref->GetLargestPossibleRegion();
	region.SetSize(2,slices);


	// create the outputs
	outputImage->SetSpacing(spacing);
	outputImage->SetDirection(direction);
	outputImage->SetOrigin(origin);
	outputImage->SetRegions(region);
	outputImage->Allocate();
	

	outputLabel->SetSpacing(spacing);
	outputLabel->SetDirection(direction);
	outputLabel->SetOrigin(origin);
	outputLabel->SetRegions(region);
	outputLabel->Allocate();

	itk::ImageRegionIterator<ImageType> outLIt(outputLabel, outputLabel->GetLargestPossibleRegion());
	itk::ImageRegionIterator<ImageType> outImIt(outputImage, outputImage->GetLargestPossibleRegion());



	// loop through the slices
	for(unsigned int i = 0; i < slices; i++)
	{
		ImageType::Pointer im = series[i].images[timestep];
		ImageType::Pointer label = series[i].labelImages[timestep];

		itk::ImageRegionConstIterator<ImageType> imIt(im, im->GetLargestPossibleRegion());
		itk::ImageRegionConstIterator<ImageType> lIt(label, label->GetLargestPossibleRegion());

		while(!imIt.IsAtEnd())
		{
			outLIt.Set(lIt.Get());
			outImIt.Set(imIt.Get());

			++imIt; ++lIt;
			++outLIt; ++outImIt;
		}
	}
}

// ------------------------------------------------------------------------
void createOutput(const ImageType::Pointer &input, ImageType::Pointer &output)
{
	output->SetDirection(input->GetDirection());
	output->SetSpacing(input->GetSpacing());
	output->SetRegions(input->GetLargestPossibleRegion());
	output->SetOrigin(input->GetOrigin());
	output->Allocate();
	output->FillBuffer(0);
}

// ------------------------------------------------------------------------
void buildShortAxisVolume(const SeriesTransform::Map &transforms, 
		const unsigned int seriesNumber, ImageType::Pointer &saVolume)
{
	// get the short axis transforms
	SeriesTransform::Map::const_iterator mapIt = transforms.begin();
	std::vector<SeriesTransform> saSlices;
	while(mapIt != transforms.end())
	{
		if(mapIt->second.series == seriesNumber)
		{
			unsigned int sliceNum = mapIt->second.slice;

			if(saSlices.size() < (sliceNum+1))
				saSlices.resize(sliceNum+1);

			saSlices[sliceNum] = mapIt->second;
		}
		
		++mapIt;
	}


	// get the dimensions of the output image
	ImageType::Pointer reference = saSlices[0].images[0];
	unsigned int x,y,z;
	x = reference->GetLargestPossibleRegion().GetSize()[0];
	y = reference->GetLargestPossibleRegion().GetSize()[1];
	z = saSlices.size();

	ImageType::RegionType region;
	ImageType::SizeType size;
	ImageType::IndexType index;
	size[0] = x;
	size[1] = y;
	size[2] = z;

	index.Fill(0);

	region.SetSize(size);
	region.SetIndex(index);
	
	// get the other parts
	ImageType::SpacingType spacing = reference->GetSpacing();
	spacing[2] = saSlices[0].sliceThickness;
	ImageType::DirectionType direction = reference->GetDirection();
	ImageType::PointType origin = reference->GetOrigin();


	saVolume->SetOrigin(origin);
	saVolume->SetDirection(direction);
	saVolume->SetSpacing(spacing);
	saVolume->SetRegions(region);
	saVolume->Allocate();
	saVolume->FillBuffer(0);
	
}

// ------------------------------------------------------------------------
void DicomParser::writeSeries(const DicomSeries &series, const QString &outputFolder)
{
	// create the series name
	QString name = QString::fromStdString(series.description);
	name.replace(" ","_");
	name += "_" + QString::number(series.images.front().seriesNumber);
	name += ".nrrd";
	
	QDir path(outputFolder);
	QString fullPath = path.absoluteFilePath(name);
	

	std::vector<DicomImage> images = series.images;
	std::sort(images.begin(), images.end());

	// write and build the output images 
	typedef itk::Image<unsigned short, 3> ImageType;
	typedef itk::Image<unsigned short, 4> OutputImageType;
	typedef itk::ImageFileReader<ImageType> ReaderType;
	typedef itk::JoinSeriesImageFilter<ImageType, OutputImageType> JoinerType;

	JoinerType::Pointer joiner = JoinerType::New();
	ImageType::Pointer orig;
	for(unsigned int i = 0; i < images.size(); i++)
	{
		ReaderType::Pointer reader = ReaderType::New();
		reader->SetFileName(images[i].filename);
		std::cout << images[i].filename << std::endl;
		reader->SetImageIO(itk::GDCMImageIO::New());
		reader->Update();

		ImageType::Pointer im = reader->GetOutput();
		if(i == 0) orig = im;


		im->SetOrigin(orig->GetOrigin());
		im->SetDirection(orig->GetDirection());
		im->SetSpacing(orig->GetSpacing());

		joiner->SetInput(i, reader->GetOutput());
	}


	std::cout << joiner->GetOutput()->GetDirection() << std::endl;


	typedef itk::ImageFileWriter<OutputImageType> WriterType;
	WriterType::Pointer writer = WriterType::New();
	writer->SetInput(joiner->GetOutput());
	writer->SetFileName(fullPath.toStdString());
	writer->SetImageIO(itk::NrrdImageIO::New());
	writer->Update();

}

// ------------------------------------------------------------------------
InitialTransformExtractor::TranslationType InitialTransformExtractor::ComputeTranslation(
    const ImageType::Pointer &image,
    const OptionsData &options)
{
    // set the offsets
    TranslationType translation;
    translation.Fill(0);
    translation[0] -= options.roiOffset[0];
    translation[1] -= options.roiOffset[1];
    translation[2] -= options.roiOffset[2];


    ImageType::PointType origin = image->GetOrigin();
    translation[0] -= origin[0];
    translation[1] -= origin[1];
    translation[2] -= origin[2];

    return translation;
}

// ------------------------------------------------------------------------
void InitialTransformExtractor::Compute()
{
    GetFilenames(m_Filenames);
    BuildMapping(m_Mapping);
    std::string referenceFilename = GetReferenceImageFilename();


    // load the reference file
    typedef itk::ImageFileReader<ImageType> ReaderType;
    ReaderType::Pointer reader = ReaderType::New();
    reader->SetFileName(referenceFilename);
    reader->SetImageIO(itk::GDCMImageIO::New());
    reader->Update();
    ImageType::Pointer refImage = reader->GetOutput();


    // flip the x and y axis
    typedef itk::PermuteAxesImageFilter<ImageType> FlipperType;
    FlipperType::Pointer flipper = FlipperType::New();
    itk::FixedArray<unsigned int, 3> order;
    order[0] = 1;
    order[1] = 0;
    order[2] = 2;

    flipper->SetOrder(order);
    flipper->SetInput(refImage);
    flipper->Update();

    //refImage = flipper->GetOutput();

    // transformation is the negative origin of the reference
    ImageType::PointType origin = refImage->GetOrigin();
    m_Translation[0] = -origin[0];
    m_Translation[1] = -origin[1];
    m_Translation[2] = -origin[2];

    m_Rotation = refImage->GetDirection().GetInverse();
    m_Reference = refImage;
}

  void NrrdQBallImageReader
    ::GenerateData()
  {
    if ( m_FileName == "")
    {
      throw itk::ImageFileReaderException(__FILE__, __LINE__, "Sorry, the filename of the vessel tree to be read is empty!");
    }
    else
    {
      try
      {
        const std::string& locale = "C";
        const std::string& currLocale = setlocale( LC_ALL, NULL );

        if ( locale.compare(currLocale)!=0 )
        {
          try
          {
            setlocale(LC_ALL, locale.c_str());
          }
          catch(...)
          {
            MITK_INFO << "Could not set locale " << locale;
          }
        }

        typedef itk::VectorImage<float,3> ImageType;
        itk::NrrdImageIO::Pointer io = itk::NrrdImageIO::New();
        typedef itk::ImageFileReader<ImageType> FileReaderType;
        FileReaderType::Pointer reader = FileReaderType::New();
        reader->SetImageIO(io);
        reader->SetFileName(this->m_FileName);
        reader->Update();
        ImageType::Pointer img = reader->GetOutput();

        typedef itk::Image<itk::Vector<float,QBALL_ODFSIZE>,3> VecImgType;
        VecImgType::Pointer vecImg = VecImgType::New();
        vecImg->SetSpacing( img->GetSpacing() );   // Set the image spacing
        vecImg->SetOrigin( img->GetOrigin() );     // Set the image origin
        vecImg->SetDirection( img->GetDirection() );  // Set the image direction
        vecImg->SetLargestPossibleRegion( img->GetLargestPossibleRegion());
        vecImg->SetBufferedRegion( img->GetLargestPossibleRegion() );
        vecImg->Allocate();

        itk::ImageRegionIterator<VecImgType> ot (vecImg, vecImg->GetLargestPossibleRegion() );
        ot = ot.Begin();

        itk::ImageRegionIterator<ImageType> it (img, img->GetLargestPossibleRegion() );

        typedef ImageType::PixelType  VarPixType;
        typedef VecImgType::PixelType FixPixType;

        for (it = it.Begin(); !it.IsAtEnd(); ++it)
        {
          VarPixType vec = it.Get();
          FixPixType fixVec(vec.GetDataPointer());
          ot.Set(fixVec);
          ++ot;
        }

        this->GetOutput()->InitializeByItk(vecImg.GetPointer());
        this->GetOutput()->SetVolume(vecImg->GetBufferPointer());

        try
        {
          setlocale(LC_ALL, currLocale.c_str());
        }
        catch(...)
        {
          MITK_INFO << "Could not reset locale " << currLocale;
        }
      }
      catch(std::exception& e)
      {
        throw itk::ImageFileReaderException(__FILE__, __LINE__, e.what());
      }
      catch(...)
      {
        throw itk::ImageFileReaderException(__FILE__, __LINE__, "Sorry, an error occurred while reading the requested vessel tree file!");
      }
    }
  }

  std::vector<itk::SmartPointer<BaseData> > NrrdTensorImageReader::Read()
  {
    std::vector<itk::SmartPointer<mitk::BaseData> > result;
    std::string location = GetInputLocation();

    if ( location == "")
    {
      throw itk::ImageFileReaderException(__FILE__, __LINE__, "Sorry, the filename is empty!");
    }
    else
    {
      try
      {
        mitk::LocaleSwitch localeSwitch("C");

        try
        {
          std::string fname3 = mitk::IOUtil::GetTempPath()+"/temp_dti.nii.gz";

          int c = 0;
          while( itksys::SystemTools::FileExists(fname3) )
          {
            fname3 = mitk::IOUtil::GetTempPath()+"/temp_dti_" + boost::lexical_cast<std::string>(c) + ".nii.gz";
            ++c;
          }

          itksys::SystemTools::CopyAFile(location.c_str(), fname3.c_str());

          typedef itk::VectorImage<float,3> ImageType;
          itk::NiftiImageIO::Pointer io = itk::NiftiImageIO::New();
          typedef itk::ImageFileReader<ImageType> FileReaderType;
          FileReaderType::Pointer reader = FileReaderType::New();
          reader->SetImageIO(io);
          reader->SetFileName(fname3);
          reader->Update();
          ImageType::Pointer img = reader->GetOutput();

          TensorImage::ItkTensorImageType::Pointer vecImg = TensorImage::ItkTensorImageType::New();
          vecImg->SetSpacing( img->GetSpacing() );   // Set the image spacing
          vecImg->SetOrigin( img->GetOrigin() );     // Set the image origin
          vecImg->SetDirection( img->GetDirection() );  // Set the image direction
          vecImg->SetRegions( img->GetLargestPossibleRegion());
          vecImg->Allocate();

          itk::ImageRegionIterator<TensorImage::ItkTensorImageType> ot (vecImg, vecImg->GetLargestPossibleRegion() );
          ot.GoToBegin();

          itk::ImageRegionIterator<ImageType> it (img, img->GetLargestPossibleRegion() );
          it.GoToBegin();

          typedef ImageType::PixelType  VarPixType;
          typedef TensorImage::PixelType FixPixType;
          int numComponents = img->GetNumberOfComponentsPerPixel();

          if (numComponents==6)
          {
            MITK_INFO << "Trying to load dti as 6-comp nifti ...";
            while (!it.IsAtEnd())
            {
              VarPixType vec = it.Get();
              FixPixType fixVec(vec.GetDataPointer());

              TensorImage::PixelType tensor;
              tensor.SetElement(0, vec.GetElement(0));
              tensor.SetElement(1, vec.GetElement(1));
              tensor.SetElement(2, vec.GetElement(2));
              tensor.SetElement(3, vec.GetElement(3));
              tensor.SetElement(4, vec.GetElement(4));
              tensor.SetElement(5, vec.GetElement(5));

              fixVec = tensor;

              ot.Set(fixVec);
              ++ot;
              ++it;
            }
          }
          else if(numComponents==9)
          {
            MITK_INFO << "Trying to load dti as 9-comp nifti ...";
            while (!it.IsAtEnd())
            {
              VarPixType vec = it.Get();
              TensorImage::PixelType tensor;
              tensor.SetElement(0, vec.GetElement(0));
              tensor.SetElement(1, vec.GetElement(1));
              tensor.SetElement(2, vec.GetElement(2));
              tensor.SetElement(3, vec.GetElement(4));
              tensor.SetElement(4, vec.GetElement(5));
              tensor.SetElement(5, vec.GetElement(8));

              FixPixType fixVec(tensor);
              ot.Set(fixVec);
              ++ot;
              ++it;
            }
          }
          else if (numComponents==1)
          {
            MITK_INFO << "Trying to load dti as 4D nifti ...";
//.........这里部分代码省略.........

        virtual ReadResult readImage(const std::string& file, const osgDB::ReaderWriter::Options* options) const
        {
            std::string ext = osgDB::getLowerCaseFileExtension(file);
            if (!acceptsExtension(ext)) return ReadResult::FILE_NOT_HANDLED;

            std::string fileName = osgDB::findDataFile( file, options );
            if (fileName.empty()) return ReadResult::FILE_NOT_FOUND;

            notice()<<"Reading DICOM file "<<fileName<<std::endl;

            typedef unsigned short PixelType;
            const unsigned int Dimension = 3;
            typedef itk::Image< PixelType, Dimension > ImageType;


            typedef itk::ImageFileReader< ImageType > ReaderType;
            ReaderType::Pointer reader = ReaderType::New();
            reader->SetFileName( fileName.c_str() );

            typedef itk::GDCMImageIO           ImageIOType;
            ImageIOType::Pointer gdcmImageIO = ImageIOType::New();
            reader->SetImageIO( gdcmImageIO );

            try
            {
                reader->Update();
            }
            catch (itk::ExceptionObject & e)
            {
                std::cerr << "exception in file reader " << std::endl;
                std::cerr << e.GetDescription() << std::endl;
                std::cerr << e.GetLocation() << std::endl;
                return ReadResult::ERROR_IN_READING_FILE;
            }
            
            ImageType::Pointer inputImage = reader->GetOutput();
            
            ImageType::RegionType region = inputImage->GetBufferedRegion();
            ImageType::SizeType size = region.GetSize();
            ImageType::IndexType start = region.GetIndex();
            
            //inputImage->GetSpacing();
            //inputImage->GetOrigin();
            
            unsigned int width = size[0];
            unsigned int height = size[1];
            unsigned int depth = size[2];
            
            osg::RefMatrix* matrix = new osg::RefMatrix;
            
            notice()<<"width = "<<width<<" height = "<<height<<" depth = "<<depth<<std::endl;
            for(unsigned int i=0; i<Dimension; ++i)
            {
                (*matrix)(i,i) = inputImage->GetSpacing()[i];
                (*matrix)(3,i) = inputImage->GetOrigin()[i];
            }
            
            osg::Image* image = new osg::Image;
            image->allocateImage(width, height, depth, GL_LUMINANCE, GL_UNSIGNED_BYTE, 1);

            unsigned char* data = image->data();            
            typedef itk::ImageRegionConstIterator< ImageType > IteratorType;
            IteratorType it(inputImage, region);
            
            it.GoToBegin();
            while (!it.IsAtEnd())
            {
                *data = it.Get();
                ++data;
                ++it;
            }
            
            image->setUserData(matrix);
            
            matrix->preMult(osg::Matrix::scale(double(image->s()), double(image->t()), double(image->r())));

            return image;
        }

// ------------------------------------------------------------------------
void NormaliseImage(const ImageType::Pointer &input, 
		const ImageType::Pointer &reference,	
		ImageType::Pointer &output,
		SeriesTransform &trans)
{
	// flip the x and y axis 
	typedef itk::PermuteAxesImageFilter<ImageType> FlipperType;
	FlipperType::Pointer flipper = FlipperType::New();
	itk::FixedArray<unsigned int, 3> order;
	order[0] = 1;
	order[1] = 0;
	order[2] = 2;

	flipper->SetOrder(order);
	flipper->SetInput(input);
	flipper->Update();

	output = flipper->GetOutput();

	// get the reference offset
	vnl_vector<double> refOrigin(3);
	refOrigin[0] = reference->GetOrigin()[0];// + trans.translation[0];
	refOrigin[1] = reference->GetOrigin()[1];// + trans.translation[1];
	refOrigin[2] = reference->GetOrigin()[2];// + trans.translation[2];

	vnl_matrix<double> refRot = reference->GetDirection().GetVnlMatrix();
	vnl_matrix<double> refRotInv = vnl_inverse(refRot);
	vnl_vector<double> refOffset = refRotInv * refOrigin;

	// get the image origin
	vnl_vector<double> origin(3);
	origin[0] = output->GetOrigin()[0];
	origin[1] = output->GetOrigin()[1];
	origin[2] = output->GetOrigin()[2];
	
	// apply the rotation to the origin
	origin = refRotInv * origin;

	// apply the offset to the origin
	origin = origin - refOffset;

	// apply the scaling 
	origin /= reference->GetSpacing()[0];

	// put the values into the output image
	ImageType::PointType itkOrigin;
	itkOrigin[0] = origin[0];
	itkOrigin[1] = origin[1];
	itkOrigin[2] = origin[2];

	ImageType::SpacingType spacing;
	spacing.Fill(output->GetSpacing()[0] / reference->GetSpacing()[0]);

	// get the new direction
	ImageType::DirectionType dirOut = output->GetDirection();
	dirOut = reference->GetDirection().GetInverse() * dirOut.GetVnlMatrix();
	
	output->SetSpacing(spacing);
	output->SetDirection(dirOut);
	output->SetOrigin(itkOrigin);
}

//.........这里部分代码省略.........
      nameAddition << "_Threshold";
      std::cout << "Thresholding successful." << std::endl;
      break;
    }

  case INVERSION:
    {
      InversionFilterType::Pointer invFilter = InversionFilterType::New();
      mitk::ScalarType min = newImage->GetScalarValueMin();
      mitk::ScalarType max = newImage->GetScalarValueMax();
      invFilter->SetMaximum( max + min );
      invFilter->SetInput(itkImage);
      invFilter->UpdateLargestPossibleRegion();
      newImage = mitk::ImportItkImage(invFilter->GetOutput())->Clone();
      nameAddition << "_Inverted";
      std::cout << "Image inversion successful." << std::endl;
      break;
    }

  case DOWNSAMPLING:
    {
      ResampleImageFilterType::Pointer downsampler = ResampleImageFilterType::New();
      downsampler->SetInput( itkImage );

      NearestInterpolatorType::Pointer interpolator = NearestInterpolatorType::New();
      downsampler->SetInterpolator( interpolator );

      downsampler->SetDefaultPixelValue( 0 );

      ResampleImageFilterType::SpacingType spacing = itkImage->GetSpacing();
      spacing *= (double) param1;
      downsampler->SetOutputSpacing( spacing );

      downsampler->SetOutputOrigin( itkImage->GetOrigin() );
      downsampler->SetOutputDirection( itkImage->GetDirection() );

      ResampleImageFilterType::SizeType size = itkImage->GetLargestPossibleRegion().GetSize();
      for ( int i = 0; i < 3; ++i )
      {
        size[i] /= param1;
      }
      downsampler->SetSize( size );
      downsampler->UpdateLargestPossibleRegion();

      newImage = mitk::ImportItkImage(downsampler->GetOutput())->Clone();
      nameAddition << "_Downsampled_by_" << param1;
      std::cout << "Downsampling successful." << std::endl;
      break;
    }

  case FLIPPING:
    {
      FlipImageFilterType::Pointer flipper = FlipImageFilterType::New();
      flipper->SetInput( itkImage );
      itk::FixedArray<bool, 3> flipAxes;
      for(int i=0; i<3; ++i)
      {
        if(i == param1)
        {
          flipAxes[i] = true;
        }
        else
        {
          flipAxes[i] = false;
        }
      }

void NrrdTensorImageReader
::GenerateData()
{
    if ( m_FileName == "")
    {
        throw itk::ImageFileReaderException(__FILE__, __LINE__, "Sorry, the filename is empty!");
    }
    else
    {
        try
        {
            const std::string& locale = "C";
            const std::string& currLocale = setlocale( LC_ALL, NULL );

            if ( locale.compare(currLocale)!=0 )
            {
                try
                {
                    setlocale(LC_ALL, locale.c_str());
                }
                catch(...)
                {
                    MITK_INFO << "Could not set locale " << locale;
                }
            }

            try
            {
                MITK_INFO << "Trying to load dti as nifti ...";
                std::string fname3 = "temp_dti.nii";
                itksys::SystemTools::CopyAFile(m_FileName.c_str(), fname3.c_str());

                typedef itk::Image<float,4> ImageType;
                itk::NiftiImageIO::Pointer io = itk::NiftiImageIO::New();
                typedef itk::ImageFileReader<ImageType> FileReaderType;
                FileReaderType::Pointer reader = FileReaderType::New();
                reader->SetImageIO(io);
                reader->SetFileName(fname3);
                reader->Update();
                itksys::SystemTools::RemoveFile(fname3.c_str());

                ImageType::Pointer img = reader->GetOutput();
                itk::Size<4> size = img->GetLargestPossibleRegion().GetSize();

                itk::ImageRegion< 3 > region;
                region.SetSize(0,size[0]);
                region.SetSize(1,size[1]);
                region.SetSize(2,size[2]);
                itk::Vector<double,3> spacing;
                spacing[0] = img->GetSpacing()[0];
                spacing[1] = img->GetSpacing()[1];
                spacing[2] = img->GetSpacing()[2];
                itk::Point<double,3> origin;
                origin[0] = img->GetOrigin()[0];
                origin[1] = img->GetOrigin()[1];
                origin[2] = img->GetOrigin()[2];
                itk::Matrix<double,3,3> direction;
                direction[0][0] = img->GetDirection()[0][0];
                direction[1][0] = img->GetDirection()[1][0];
                direction[2][0] = img->GetDirection()[2][0];
                direction[0][1] = img->GetDirection()[0][1];
                direction[1][1] = img->GetDirection()[1][1];
                direction[2][1] = img->GetDirection()[2][1];
                direction[0][2] = img->GetDirection()[0][2];
                direction[1][2] = img->GetDirection()[1][2];
                direction[2][2] = img->GetDirection()[2][2];

                typedef itk::Image<itk::DiffusionTensor3D<float>,3> VecImgType;
                typedef VecImgType::PixelType FixPixType;
                VecImgType::Pointer vecImg = VecImgType::New();
                vecImg->SetSpacing( spacing );
                vecImg->SetOrigin( origin );
                vecImg->SetDirection( direction );
                vecImg->SetRegions( region );
                vecImg->Allocate();

                itk::ImageRegionIterator<VecImgType> ot (vecImg, vecImg->GetLargestPossibleRegion() );
                ot = ot.Begin();

                while (!ot.IsAtEnd())
                {
                    itk::DiffusionTensor3D<float> tensor;
                    ImageType::IndexType idx;
                    idx[0] = ot.GetIndex()[0]; idx[1] = ot.GetIndex()[1]; idx[2] = ot.GetIndex()[2];

                    if (size[3]==6)
                    {
                        for (int te=0; te<size[3]; te++)
                        {
                            idx[3] = te;
                            tensor.SetElement(te, img->GetPixel(idx));
                        }
                    }
                    else if (size[3]==9)
                    {
                        idx[3] = 0;
                        tensor.SetElement(0, img->GetPixel(idx));
                        idx[3] = 1;
                        tensor.SetElement(1, img->GetPixel(idx));
                        idx[3] = 2;
//.........这里部分代码省略.........