1 /*-------------------------------------------------------------------------
3 Copyright (c) CREATIS (Centre de Recherche et d'Applications en Traitement de
4 l'Image). All rights reserved. See Doc/License.txt or
5 http://www.creatis.insa-lyon.fr/Public/Gdcm/License.html for details.
7 This software is distributed WITHOUT ANY WARRANTY; without even
8 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
9 PURPOSE. See the above copyright notices for more information.
11 -------------------------------------------------------------------------*/
13 #ifndef CLITKIMAGEARITHMGENERICFILTER_TXX
14 #define CLITKIMAGEARITHMGENERICFILTER_TXX
16 /*-------------------------------------------------
17 * @file clitkImageArithmGenericFilter.txx
18 * @author David Sarrut <David.Sarrut@creatis.insa-lyon.fr>
21 -------------------------------------------------*/
23 //--------------------------------------------------------------------
24 template<unsigned int Dim>
25 void ImageArithmGenericFilter::Update_WithDim() {
26 #define TRY_TYPE(TYPE) \
27 if (IsSameType<TYPE>(mPixelTypeName)) { Update_WithDimAndPixelType<Dim, TYPE>(); return; }
33 TRY_TYPE(unsigned int);
38 std::string list = CreateListOfTypes<char, uchar, short, ushort, int, uint, float, double>();
39 std::cerr << "Error, I don't know the type '" << mPixelTypeName << "' for the input image '"
40 << mInputFilenames[0] << "'." << std::endl << "Known types are " << list << std::endl;
43 //--------------------------------------------------------------------
45 //--------------------------------------------------------------------
46 template<unsigned int Dim, class PixelType>
47 void ImageArithmGenericFilter::Update_WithDimAndPixelType() {
49 typedef itk::Image<PixelType,Dim> ImageType;
50 typename ImageType::Pointer input1 = clitk::readImage<ImageType>(mInputFilenames[0], mIOVerbose);
51 typename ImageType::Pointer outputImage;
53 // Read input2 (float is ok altough it could take too memory ...)
54 if (mIsOperationUseASecondImage) {
55 typedef itk::Image<float,Dim> ImageType2;
56 typename ImageType2::Pointer input2 = clitk::readImage<ImageType2>(mInputFilenames[1], mIOVerbose);
57 outputImage = ComputeImage<ImageType, ImageType2>(input1, input2);
60 outputImage = ComputeImage<ImageType>(input1);
64 this->SetNextOutput<ImageType>(outputImage);
65 //clitk::writeImage<ImageType>(outputImage, mOutputFilename, mIOVerbose);
67 //--------------------------------------------------------------------
69 //--------------------------------------------------------------------
70 template<class ImageType>
71 typename ImageType::Pointer
72 ImageArithmGenericFilter::ComputeImage(typename ImageType::Pointer inputImage) {
74 typedef typename ImageType::PixelType PixelType;
75 typedef itk::ImageRegionIterator<ImageType> IteratorType;
76 IteratorType it(inputImage, inputImage->GetLargestPossibleRegion());
79 switch (mTypeOfOperation) {
81 while (!it.IsAtEnd()) {
82 it.Set(PixelTypeDownCast<double, PixelType>((double)it.Get() + mScalar) );
87 while (!it.IsAtEnd()) {
88 it.Set(PixelTypeDownCast<double, PixelType>((double)it.Get() * mScalar) );
93 while (!it.IsAtEnd()) {
95 it.Set(PixelTypeDownCast<double, PixelType>(mScalar / (double)it.Get()));
96 else it.Set(mDefaultPixelValue);
101 while (!it.IsAtEnd()) {
102 if (it.Get() < mScalar) it.Set(PixelTypeDownCast<double, PixelType>(mScalar));
107 while (!it.IsAtEnd()) {
108 if (it.Get() > mScalar) it.Set(PixelTypeDownCast<double, PixelType>(mScalar));
112 case 5: // Absolute value
113 while (!it.IsAtEnd()) {
114 if (it.Get() <= 0) it.Set(PixelTypeDownCast<double, PixelType>(-it.Get()));
115 // <= zero to avoid warning for unsigned types
119 case 6: // Squared value
120 while (!it.IsAtEnd()) {
121 it.Set(PixelTypeDownCast<double, PixelType>((double)it.Get()*(double)it.Get()));
126 while (!it.IsAtEnd()) {
128 it.Set(PixelTypeDownCast<double, PixelType>(log((double)it.Get())));
129 else it.Set(mDefaultPixelValue);
134 while (!it.IsAtEnd()) {
135 it.Set(PixelTypeDownCast<double, PixelType>(exp((double)it.Get())));
140 while (!it.IsAtEnd()) {
142 it.Set(PixelTypeDownCast<double, PixelType>(sqrt((double)it.Get())));
144 if (it.Get() ==0) it.Set(0);
145 else it.Set(mDefaultPixelValue);
151 std::cerr << "ERROR : the operation number (" << mTypeOfOperation << ") is not known." << std::endl;
157 //--------------------------------------------------------------------
159 //--------------------------------------------------------------------
160 template<class ImageType1, class ImageType2>
161 typename ImageType1::Pointer
162 ImageArithmGenericFilter::ComputeImage(typename ImageType1::Pointer inputImage1,
163 typename ImageType2::Pointer inputImage2) {
165 typedef typename ImageType1::PixelType PixelType;
166 typedef itk::ImageRegionIterator<ImageType1> IteratorType;
167 IteratorType it1(inputImage1, inputImage1->GetLargestPossibleRegion());
170 typedef itk::ImageRegionConstIterator<ImageType2> ConstIteratorType;
171 ConstIteratorType it2(inputImage2, inputImage2->GetLargestPossibleRegion());
174 switch (mTypeOfOperation) {
176 while (!it1.IsAtEnd()) {
177 it1.Set(PixelTypeDownCast<double, PixelType>((double)it1.Get() + (double)it2.Get()) );
182 while (!it1.IsAtEnd()) {
183 it1.Set(PixelTypeDownCast<double, PixelType>((double)it1.Get() * (double)it2.Get()) );
188 while (!it1.IsAtEnd()) {
190 it1.Set(PixelTypeDownCast<double, PixelType>((double)it1.Get() / (double)it2.Get()));
191 else it1.Set(mDefaultPixelValue);
196 while (!it1.IsAtEnd()) {
197 if (it1.Get() < it2.Get()) it1.Set(PixelTypeDownCast<double, PixelType>(it2.Get()));
202 while (!it1.IsAtEnd()) {
203 if (it1.Get() > it2.Get()) it1.Set(PixelTypeDownCast<double, PixelType>(it2.Get()));
207 case 5: // Absolute difference
208 while (!it1.IsAtEnd()) {
209 it1.Set(PixelTypeDownCast<double, PixelType>(fabs(it2.Get()-it1.Get())));
213 case 6: // Squared differences
214 while (!it1.IsAtEnd()) {
215 it1.Set(PixelTypeDownCast<double, PixelType>(pow((double)it1.Get()-(double)it2.Get(),2)));
219 case 7: // Difference
220 while (!it1.IsAtEnd()) {
221 it1.Set(PixelTypeDownCast<double, PixelType>((double)it1.Get()-(double)it2.Get()));
225 case 8: // Relative Difference
226 while (!it1.IsAtEnd()) {
227 if (it1.Get() != 0) it1.Set(PixelTypeDownCast<double, PixelType>(((double)it1.Get()-(double)it2.Get()))/(double)it1.Get());
233 std::cerr << "ERROR : the operation number (" << mTypeOfOperation << ") is not known." << std::endl;
239 //--------------------------------------------------------------------
241 #endif //#define CLITKIMAGEARITHMGENERICFILTER_TXX