1 /*=========================================================================
2 *
3 * Copyright Insight Software Consortium
4 *
5 * Licensed under the Apache License, Version 2.0 (the "License");
6 * you may not use this file except in compliance with the License.
7 * You may obtain a copy of the License at
8 *
9 * http://www.apache.org/licenses/LICENSE-2.0.txt
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 *
17 *=========================================================================*/
18 #include "itkMeanSquaresImageToImageMetricv4.h"
19 #include "itkTranslationTransform.h"
20 #include "itkTestingMacros.h"
21
22 /*
23 * Simple test to run using unix 'time' function for speed test.
24 */
25
itkMeanSquaresImageToImageMetricv4SpeedTest(int argc,char * argv[])26 int itkMeanSquaresImageToImageMetricv4SpeedTest(int argc, char *argv[] )
27 {
28 if( argc < 3 )
29 {
30 std::cerr << "usage: " << itkNameOfTestExecutableMacro(argv) << ": image-dimension number-of-reps" << std::endl;
31 return EXIT_FAILURE;
32 }
33 int imageSize = std::stoi( argv[1] );
34 int numberOfReps = std::stoi( argv[2] );
35
36 std::cout << "image dim: " << imageSize << ", reps: " << numberOfReps << std::endl;
37
38 constexpr unsigned int imageDimensionality = 3;
39 using ImageType = itk::Image< double, imageDimensionality >;
40
41 ImageType::SizeType size;
42 size.Fill( imageSize );
43 ImageType::IndexType index;
44 index.Fill( 0 );
45 ImageType::RegionType region;
46 region.SetSize( size );
47 region.SetIndex( index );
48 ImageType::SpacingType spacing;
49 spacing.Fill(1.0);
50 ImageType::PointType origin;
51 origin.Fill(0);
52 ImageType::DirectionType direction;
53 direction.SetIdentity();
54
55 /* Create simple test images. */
56 ImageType::Pointer fixedImage = ImageType::New();
57 fixedImage->SetRegions( region );
58 fixedImage->SetSpacing( spacing );
59 fixedImage->SetOrigin( origin );
60 fixedImage->SetDirection( direction );
61 fixedImage->Allocate();
62
63 ImageType::Pointer movingImage = ImageType::New();
64 movingImage->SetRegions( region );
65 movingImage->SetSpacing( spacing );
66 movingImage->SetOrigin( origin );
67 movingImage->SetDirection( direction );
68 movingImage->Allocate();
69
70 /* Fill images */
71 itk::ImageRegionIterator<ImageType> itFixed( fixedImage, region );
72 itFixed.GoToBegin();
73 unsigned int count = 1;
74 while( !itFixed.IsAtEnd() )
75 {
76 itFixed.Set( count );
77 count++;
78 ++itFixed;
79 }
80
81 itk::ImageRegionIteratorWithIndex<ImageType> itMoving( movingImage, region );
82
83 itMoving.GoToBegin();
84 count = 1;
85
86 while( !itMoving.IsAtEnd() )
87 {
88 itMoving.Set( count*count );
89 count++;
90 ++itMoving;
91 }
92
93 /* Transforms */
94 using FixedTransformType = itk::TranslationTransform<double,imageDimensionality>;
95 using MovingTransformType = itk::TranslationTransform<double,imageDimensionality>;
96
97 FixedTransformType::Pointer fixedTransform = FixedTransformType::New();
98 MovingTransformType::Pointer movingTransform = MovingTransformType::New();
99
100 fixedTransform->SetIdentity();
101 movingTransform->SetIdentity();
102
103 /* The metric */
104 using MetricType = itk::MeanSquaresImageToImageMetricv4< ImageType, ImageType, ImageType >;
105
106 MetricType::Pointer metric = MetricType::New();
107
108 /* Assign images and transforms.
109 * By not setting a virtual domain image or virtual domain settings,
110 * the metric will use the fixed image for the virtual domain. */
111 metric->SetFixedImage( fixedImage );
112 metric->SetMovingImage( movingImage );
113 metric->SetFixedTransform( fixedTransform );
114 metric->SetMovingTransform( movingTransform );
115
116 /* Initialize. */
117 std::cout << "Calling Initialize..." << std::endl;
118 metric->Initialize();
119
120 // Evaluate with GetValueAndDerivative
121 MetricType::MeasureType valueReturn1;
122 MetricType::DerivativeType derivativeReturn;
123
124 MetricType::MeasureType sum = itk::NumericTraits<MetricType::MeasureType>::ZeroValue();
125 for( int r=0; r < numberOfReps; r++ )
126 {
127 metric->GetValueAndDerivative( valueReturn1, derivativeReturn );
128 //Sum results to prevent optimizations
129 sum += valueReturn1 + derivativeReturn[0];
130 }
131
132 std::cout << "sum: " << sum << std::endl;
133
134 return EXIT_SUCCESS;
135 }
136