/*=========================================================================
Program: vv
Language: C++
Author : Joel Schaerer (joel.schaerer@insa-lyon.fr)
Copyright (C) 2008
Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr
CREATIS-LRMN http://www.creatis.insa-lyon.fr
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, version 3 of the License.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program. If not, see .
=========================================================================*/
#include
#include
#include
#include
#include
#include
#include
#include
#include "vvToITK.h"
#include "vvFromITK.h"
#include "vvImageWarp.h"
#include "clitkCommon.h"
//====================================================================
vvImageWarp::vvImageWarp(vvImage::Pointer input,vvImage::Pointer vf,unsigned int ref_image,QWidget* parent):
mRefImage(ref_image),
parent_window(parent),
mInputImage(input),
mVF(vf)
{}
//====================================================================
template
void vvImageWarp::Update_WithDim() {
#define TRY_TYPE(TYPE) \
if (clitk::IsSameType(mInputImage->GetScalarTypeAsString())) { this->Update_WithDimAndPixelType(); return; }
// TRY_TYPE(signed char);
// TRY_TYPE(uchar);
TRY_TYPE(short);
// TRY_TYPE(ushort);
// TRY_TYPE(int); // no uint ...
// TRY_TYPE(float);
//TRY_TYPE(double);
#undef TRY_TYPE
std::string list = clitk::CreateListOfTypes();
std::cerr << "Error, I don't know the type '" << mInputImage->GetScalarTypeAsString() << "' for the input image. "
<< std::endl << "Known types are " << list << std::endl;
exit(0);
}
//====================================================================
template
void vvImageWarp::Update_WithDimAndPixelType() {
QProgressDialog progress(parent_window);
progress.setCancelButton(0);
progress.setLabelText("Computing warped and ventilation images...");
progress.show();
typedef itk::Image ImageType;
typedef itk::Image JacobianImageType;
typedef itk::Image,Dim> VectorImageType;
typedef std::vector ImageSeriesType;
ImageSeriesType input = vvImageToITKImageVector(mInputImage);
ImageSeriesType output;
typename itk::VTKImageToImageFilter::Pointer vf_connector=
itk::VTKImageToImageFilter::New();
//Warp, then join, then convert to vv
typedef itk::WarpImageFilter WarpFilterType;
typedef itk::DisplacementFieldJacobianDeterminantFilter JacobianFilterType;
vvImage::Pointer result=vvImage::New();
typedef itk::JoinSeriesImageFilter< ImageType,itk::Image > JoinFilterType;
typedef itk::JoinSeriesImageFilter< JacobianImageType,itk::Image > JacobianJoinFilterType;
typename JoinFilterType::Pointer join=JoinFilterType::New();
typename JoinFilterType::Pointer diff_join=JoinFilterType::New();
typename JacobianJoinFilterType::Pointer jacobian_join=JacobianJoinFilterType::New();
join->SetSpacing(1);
join->SetOrigin(0); //Set the temporal origin
diff_join->SetSpacing(1);
diff_join->SetOrigin(0);
jacobian_join->SetSpacing(1);
jacobian_join->SetOrigin(0);
typedef itk::SubtractImageFilter DiffFilter;
std::vector warped_images;
for (unsigned int num = 0; num < input.size(); num++)
{
typename WarpFilterType::Pointer warp_filter=WarpFilterType::New();
typename JacobianFilterType::Pointer jacobian_filter=JacobianFilterType::New();
jacobian_filter->SetUseImageSpacingOn();
vf_connector->SetInput(mVF->GetVTKImages()[num]);
warp_filter->SetInput(input[num]);
warp_filter->SetDeformationField(vf_connector->GetOutput());
jacobian_filter->SetInput(vf_connector->GetOutput());
warp_filter->SetOutputSpacing(input[num]->GetSpacing());
warp_filter->SetOutputOrigin(input[num]->GetOrigin());
warp_filter->SetEdgePaddingValue(-1000);
warp_filter->Update();
jacobian_filter->Update();
warped_images.push_back(warp_filter->GetOutput());
jacobian_join->PushBackInput(jacobian_filter->GetOutput());
join->PushBackInput(warp_filter->GetOutput());
progress.setValue(progress.value()+1);
}
for (typename std::vector::const_iterator i = warped_images.begin(); i!=warped_images.end();i++)
{
typename DiffFilter::Pointer diff_filter = DiffFilter::New();
diff_filter->SetInput2(*i);
diff_filter->SetInput1(*(warped_images.begin()+mRefImage));
diff_filter->Update();
diff_join->PushBackInput(diff_filter->GetOutput());
progress.setValue(progress.value()+1);
}
join->Update();
diff_join->Update();
jacobian_join->Update();
mWarpedImage = vvImageFromITK(join->GetOutput());
mDiffImage = vvImageFromITK(diff_join->GetOutput());
mJacobianImage = vvImageFromITK(jacobian_join->GetOutput());
//mJacobianImage = vvImageFromITK(temporal_filter->GetOutput());
}
//====================================================================
//
bool vvImageWarp::ComputeWarpedImage()
{
for (int i=0;iGetNumberOfDimensions();i++)
{
if (mInputImage->GetSpacing()[i] != mVF->GetSpacing()[i])
return false;
}
switch (mInputImage->GetNumberOfDimensions())
{
// case 2: this->Update_WithDim<2>(); break;;
// case 3: this->Update_WithDim<3>(); break;;
case 4:
this->Update_WithDim<3>();
break;;
default:
DD("Error: dimension not handled.");
}
return true;
}