selxWBIRDemoTest.cxx 29.5 KB
Newer Older
Floris Berendsen's avatar
Floris Berendsen committed
1
2
/*=========================================================================
 *
3
 *  Copyright Leiden University Medical Center, Erasmus University Medical
Floris Berendsen's avatar
Floris Berendsen committed
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
 *  Center and contributors
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *        http://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 *
 *=========================================================================*/

20
#include "selxSuperElastixFilter.h"
Floris Berendsen's avatar
Floris Berendsen committed
21
22
23
24
25
26

#include "selxItkSmoothingRecursiveGaussianImageFilterComponent.h"
#include "selxDisplacementFieldItkImageFilterSink.h"
#include "selxItkImageSource.h"

#include "selxElastixComponent.h"
27
28
29
#include "selxMonolithicElastix.h"
#include "selxMonolithicTransformix.h"

Floris Berendsen's avatar
Floris Berendsen committed
30
31
32
33
34
#include "selxItkImageSink.h"

#include "selxItkImageRegistrationMethodv4Component.h"
#include "selxItkANTSNeighborhoodCorrelationImageToImageMetricv4.h"
#include "selxItkMeanSquaresImageToImageMetricv4.h"
35
36
37
38
#include "selxItkGradientDescentOptimizerv4.h"
#include "selxItkGaussianExponentialDiffeomorphicTransform.h"
#include "selxItkTransformDisplacementFilter.h"
#include "selxItkResampleFilter.h"
39
#include "selxRegistrationController.h"
Floris Berendsen's avatar
Floris Berendsen committed
40
41
42
#include "selxItkImageSourceFixed.h"
#include "selxItkImageSourceMoving.h"

43
#include "selxDefaultComponents.h"
44

45
46
47
48
49
#include "itkTestingComparisonImageFilter.h"
#include "itkSubtractImageFilter.h"
#include "itkVectorMagnitudeImageFilter.h"
#include "itkStatisticsImageFilter.h"

50
#include "selxDataManager.h"
Floris Berendsen's avatar
Floris Berendsen committed
51
52
53
54
55
#include "gtest/gtest.h"

#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"

56
57
namespace selx
{
58
/** A demo for our WBIR paper written as a Unit Test in the Google Test Framework */
59
60
class WBIRDemoTest : public ::testing::Test
{
Floris Berendsen's avatar
Floris Berendsen committed
61
public:
62

63
  /** Fill SUPERelastix' component data base by registering various components */
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
  typedef TypeList<
    DisplacementFieldItkImageFilterSinkComponent< 2, float >,
    ItkImageSinkComponent< 2, float >,
    ItkImageSourceFixedComponent< 2, float >,
    ItkImageSourceMovingComponent< 2, float >,
    ElastixComponent< 2, float >,
    MonolithicElastixComponent< 2, float >,
    MonolithicTransformixComponent< 2, float >,
    ItkImageRegistrationMethodv4Component< 2, float >,
    ItkANTSNeighborhoodCorrelationImageToImageMetricv4Component< 2, float >,
    ItkMeanSquaresImageToImageMetricv4Component< 2, float >,
    ItkGradientDescentOptimizerv4Component< double >,
    ItkGaussianExponentialDiffeomorphicTransformComponent< double, 2 >,
    ItkTransformDisplacementFilterComponent< 2, float, double >,
    ItkResampleFilterComponent< 2, float, double >,
79
    RegistrationControllerComponent< > > RegisterComponents;
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98

  typedef SuperElastixFilter< RegisterComponents > SuperElastixFilterType;

  typedef Blueprint::Pointer            BlueprintPointerType;
  typedef Blueprint::ConstPointer       BlueprintConstPointerType;
  typedef Blueprint::ParameterMapType   ParameterMapType;
  typedef Blueprint::ParameterValueType ParameterValueType;
  typedef DataManager                   DataManagerType;

  typedef itk::Image< float, 2 >              Image2DType;
  typedef itk::ImageFileReader< Image2DType > ImageReader2DType;
  typedef itk::ImageFileWriter< Image2DType > ImageWriter2DType;

  typedef itk::Image< itk::Vector< float, 2 >, 2 >  VectorImage2DType;
  typedef itk::ImageFileReader< VectorImage2DType > VectorImageReader2DType;
  typedef itk::ImageFileWriter< VectorImage2DType > VectorImageWriter2DType;

  typedef itk::Testing::ComparisonImageFilter< Image2DType, Image2DType > ComparisonImageFilterType;
  // Unfortunately comparing Vector Images cannot by done by this filter, therefore we define our own comparison pipeline.
99
100
  // typedef itk::Testing::ComparisonImageFilter<VectorImage2DType, VectorImage2DType> ComparisonVectorImageFilterType;
  //typedef itk::VectorImageToImageAdaptor<float, 2> VectorImageToImageAdaptorType;
101
102
103
  typedef itk::SubtractImageFilter< VectorImage2DType, VectorImage2DType, VectorImage2DType > SubtractVectorImageFilterType;
  typedef itk::VectorMagnitudeImageFilter< VectorImage2DType, Image2DType >                   VectorMagnitudeImageFilterType;
  typedef itk::StatisticsImageFilter< Image2DType >                                           StatisticsImageFilterType;
104

105
106
  virtual void SetUp()
  {
107
    baselineImageReader = ImageReader2DType::New();
108
109
110
    compareImageFilter  = ComparisonImageFilterType::New();
    compareImageFilter->SetValidInput( baselineImageReader->GetOutput() );
    compareImageFilter->SetDifferenceThreshold( 0.0 );
111
112
113
114
115
116

    baselineVectorImageReader = VectorImageReader2DType::New();
    subtractVectorImageFilter = SubtractVectorImageFilterType::New();
    VectorMagnitudeImageFilterType::Pointer vectorMagnitudeImageFilter = VectorMagnitudeImageFilterType::New();
    statisticsImageFilter = StatisticsImageFilterType::New();

117
118
119
120
    subtractVectorImageFilter->SetInput1( baselineVectorImageReader->GetOutput() );

    vectorMagnitudeImageFilter->SetInput( subtractVectorImageFilter->GetOutput() );
    statisticsImageFilter->SetInput( vectorMagnitudeImageFilter->GetOutput() );
Floris Berendsen's avatar
Floris Berendsen committed
121
122
  }

123
124
125

  virtual void TearDown()
  {
Floris Berendsen's avatar
Floris Berendsen committed
126
127
128
    itk::ObjectFactoryBase::UnRegisterAllFactories();
  }

129
130

  BlueprintPointerType            blueprint;
131
  SuperElastixFilterType::Pointer superElastixFilter;
132

133
134
135
  ImageReader2DType::Pointer             baselineImageReader;
  ComparisonImageFilterType::Pointer     compareImageFilter;
  VectorImageReader2DType::Pointer       baselineVectorImageReader;
136
  SubtractVectorImageFilterType::Pointer subtractVectorImageFilter;
137
  StatisticsImageFilterType::Pointer     statisticsImageFilter;
Floris Berendsen's avatar
Floris Berendsen committed
138
139
};

140
/** Experiment 2a: ITKv4 framework, stationary velocity field transform, ANTs neighborhood correlation metric */
141
TEST_F( WBIRDemoTest, itkv4_SVF_ANTSCC )
Floris Berendsen's avatar
Floris Berendsen committed
142
143
144
145
{
  /** make example blueprint configuration */
  blueprint = Blueprint::New();

146
147
148
149
150
151
  blueprint->AddComponent( "RegistrationMethod", { { "NameOfClass", { "ItkImageRegistrationMethodv4Component" } } } );
  blueprint->AddComponent( "Metric", { { "NameOfClass", { "ItkANTSNeighborhoodCorrelationImageToImageMetricv4Component" } } } );
  blueprint->AddComponent( "Optimizer", { { "NameOfClass", { "ItkGradientDescentOptimizerv4Component" } },
                                          { "NumberOfIterations", { "100" } },
                                          { "LearningRate", { "100" } } } );
  blueprint->AddComponent( "Transform", { { "NameOfClass", { "ItkGaussianExponentialDiffeomorphicTransformComponent" } } } );
Floris Berendsen's avatar
Floris Berendsen committed
152

153
154
  blueprint->AddComponent( "ResampleFilter", { { "NameOfClass", { "ItkResampleFilterComponent" } } } );
  blueprint->AddComponent( "TransformDisplacementFilter", { { "NameOfClass", { "ItkTransformDisplacementFilterComponent" } } } );
Floris Berendsen's avatar
Floris Berendsen committed
155

156
157
158
159
160
  blueprint->AddComponent( "FixedImageSource", { { "NameOfClass", { "ItkImageSourceFixedComponent" } } } );
  blueprint->AddComponent( "MovingImageSource", { { "NameOfClass", { "ItkImageSourceMovingComponent" } } } );
  blueprint->AddComponent( "ResultImageSink", { { "NameOfClass", { "ItkImageSinkComponent" } } } );
  blueprint->AddComponent( "ResultDisplacementFieldSink", { { "NameOfClass", { "DisplacementFieldItkImageFilterSinkComponent" } } } );
  blueprint->AddComponent( "Controller", { { "NameOfClass", { "RegistrationControllerComponent" } } } );
Floris Berendsen's avatar
Floris Berendsen committed
161

162
  //optionally, tie properties to connection to avoid ambiguities
163
  //blueprint->AddConnection("FixedImageSource", "RegistrationMethod", { { "NameOfInterface", { "itkImageFixedInterface" } } });
164
  blueprint->AddConnection( "FixedImageSource", "RegistrationMethod", { {} } );
Floris Berendsen's avatar
Floris Berendsen committed
165

166
  //optionally, tie properties to connection to avoid ambiguities
167
  //blueprint->AddConnection("MovingImageSource", "RegistrationMethod", { { "NameOfInterface", { "itkImageMovingInterface" } } });
168
  blueprint->AddConnection( "MovingImageSource", "RegistrationMethod", { {} } );
Floris Berendsen's avatar
Floris Berendsen committed
169

170
  //optionally, tie properties to connection to avoid ambiguities
171
  //blueprint->AddConnection("RegistrationMethod", "ResultImageSink", { { "NameOfInterface", { "itkImageSourceInterface" } } });
172
  blueprint->AddConnection( "ResampleFilter", "ResultImageSink", { {} } );
Floris Berendsen's avatar
Floris Berendsen committed
173

174
  //optionally, tie properties to connection to avoid ambiguities
175
  //blueprint->AddConnection("RegistrationMethod", "ResultDisplacementFieldSink", { { "NameOfInterface", { "DisplacementFieldItkImageSourceInterface" } } });
176
  blueprint->AddConnection( "TransformDisplacementFilter", "ResultDisplacementFieldSink", { {} } );
Floris Berendsen's avatar
Floris Berendsen committed
177

178
  //optionally, tie properties to connection to avoid ambiguities
179
  //blueprint->AddConnection("Metric", "RegistrationMethod", { { "NameOfInterface", { "itkMetricv4Interface" } } });
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
  blueprint->AddConnection( "Metric", "RegistrationMethod", { {} } );

  blueprint->AddConnection( "FixedImageSource", "Transform", { {} } );
  blueprint->AddConnection( "Transform", "RegistrationMethod", { {} } );
  blueprint->AddConnection( "Optimizer", "RegistrationMethod", { {} } );
  blueprint->AddConnection( "RegistrationMethod", "TransformDisplacementFilter", { {} } );
  blueprint->AddConnection( "FixedImageSource", "TransformDisplacementFilter", { {} } );
  blueprint->AddConnection( "RegistrationMethod", "ResampleFilter", { {} } );
  blueprint->AddConnection( "FixedImageSource", "ResampleFilter", { {} } );
  blueprint->AddConnection( "MovingImageSource", "ResampleFilter", { {} } );

  blueprint->AddConnection( "RegistrationMethod", "Controller", { {} } );          //RunRegistrationInterface
  blueprint->AddConnection( "ResampleFilter", "Controller", { {} } );              //ReconnectTransformInterface
  blueprint->AddConnection( "TransformDisplacementFilter", "Controller", { {} } ); //ReconnectTransformInterface
  blueprint->AddConnection( "ResultImageSink", "Controller", { {} } );             //AfterRegistrationInterface
  blueprint->AddConnection( "ResultDisplacementFieldSink", "Controller", { {} } ); //AfterRegistrationInterface
196
197
  // Data manager provides the paths to the input and output data for unit tests
  DataManagerType::Pointer dataManager = DataManagerType::New();
Floris Berendsen's avatar
Floris Berendsen committed
198

199
  blueprint->WriteBlueprint( dataManager->GetOutputFile( "itkv4_SVF_ANTSCC.dot" ) );
200

201
  // Instantiate SuperElastix
202
  EXPECT_NO_THROW( superElastixFilter = SuperElastixFilterType::New() );
Floris Berendsen's avatar
Floris Berendsen committed
203

204
  // Set up the readers and writers
205
  ImageReader2DType::Pointer fixedImageReader = ImageReader2DType::New();
206
  fixedImageReader->SetFileName( dataManager->GetInputFile( "coneA2d64.mhd" ) );
Floris Berendsen's avatar
Floris Berendsen committed
207

208
  ImageReader2DType::Pointer movingImageReader = ImageReader2DType::New();
209
210
  movingImageReader->SetFileName( dataManager->GetInputFile( "coneB2d64.mhd" ) );

211
  ImageWriter2DType::Pointer resultImageWriter = ImageWriter2DType::New();
212
  resultImageWriter->SetFileName( dataManager->GetOutputFile( "itkv4_SVF_ANTSCC_Image.mhd" ) );
Floris Berendsen's avatar
Floris Berendsen committed
213

214
  VectorImageWriter2DType::Pointer vectorImageWriter = VectorImageWriter2DType::New();
215
216
  vectorImageWriter->SetFileName( dataManager->GetOutputFile( "itkv4_SVF_ANTSCC_Displacement.mhd" ) );

217
  // Connect SuperElastix in an itk pipeline
218
219
  superElastixFilter->SetInput( "FixedImageSource", fixedImageReader->GetOutput() );
  superElastixFilter->SetInput( "MovingImageSource", movingImageReader->GetOutput() );
220

221
222
  resultImageWriter->SetInput( superElastixFilter->GetOutput< Image2DType >( "ResultImageSink" ) );
  vectorImageWriter->SetInput( superElastixFilter->GetOutput< VectorImage2DType >( "ResultDisplacementFieldSink" ) );
223

224
225
  baselineImageReader->SetFileName( dataManager->GetBaselineFile( "itkv4_SVF_ANTSCC_Image.mhd" ) );
  baselineVectorImageReader->SetFileName( dataManager->GetBaselineFile( "itkv4_SVF_ANTSCC_Displacement.mhd" ) );
226

227
228
  compareImageFilter->SetTestInput( superElastixFilter->GetOutput< Image2DType >( "ResultImageSink" ) );
  subtractVectorImageFilter->SetInput2( superElastixFilter->GetOutput< VectorImage2DType >( "ResultDisplacementFieldSink" ) );
229

230
  EXPECT_NO_THROW( superElastixFilter->SetBlueprint( blueprint ) );
231
232
233
234
235

  //Optional Update call
  //superElastixFilter->Update();

  // Update call on the writers triggers SuperElastix to configure and execute
236
237
  EXPECT_NO_THROW( resultImageWriter->Update() );
  EXPECT_NO_THROW( vectorImageWriter->Update() );
238

239
240
  EXPECT_NO_THROW( compareImageFilter->Update() );
  EXPECT_NO_THROW( statisticsImageFilter->Update() );
241

242
243
  EXPECT_EQ( 0, compareImageFilter->GetNumberOfPixelsWithDifferences() );
  EXPECT_LT( statisticsImageFilter->GetSumOutput()->Get(), 1e-16 );
244
245
246
}

/** Experiment 2b: ITKv4 framework, stationary velocity field transform, mean squared differences metric */
247
TEST_F( WBIRDemoTest, itkv4_SVF_MSD )
248
249
250
251
{
  /** make example blueprint configuration */
  blueprint = Blueprint::New();

252
253
254
255
256
257
  blueprint->AddComponent( "RegistrationMethod", { { "NameOfClass", { "ItkImageRegistrationMethodv4Component" } } } );
  blueprint->AddComponent( "Metric", { { "NameOfClass", { "ItkMeanSquaresImageToImageMetricv4Component" } } } );
  blueprint->AddComponent( "Optimizer", { { "NameOfClass", { "ItkGradientDescentOptimizerv4Component" } },
                                          { "NumberOfIterations", { "100" } },
                                          { "LearningRate", { "0.001" } } } );
  blueprint->AddComponent( "Transform", { { "NameOfClass", { "ItkGaussianExponentialDiffeomorphicTransformComponent" } } } );
258

259
260
  blueprint->AddComponent( "ResampleFilter", { { "NameOfClass", { "ItkResampleFilterComponent" } } } );
  blueprint->AddComponent( "TransformDisplacementFilter", { { "NameOfClass", { "ItkTransformDisplacementFilterComponent" } } } );
261

262
263
264
265
266
  blueprint->AddComponent( "FixedImageSource", { { "NameOfClass", { "ItkImageSourceFixedComponent" } } } );
  blueprint->AddComponent( "MovingImageSource", { { "NameOfClass", { "ItkImageSourceMovingComponent" } } } );
  blueprint->AddComponent( "ResultImageSink", { { "NameOfClass", { "ItkImageSinkComponent" } } } );
  blueprint->AddComponent( "ResultDisplacementFieldSink", { { "NameOfClass", { "DisplacementFieldItkImageFilterSinkComponent" } } } );
  blueprint->AddComponent( "Controller", { { "NameOfClass", { "RegistrationControllerComponent" } } } );
267
268

  //optionally, tie properties to connection to avoid ambiguities
269
  //blueprint->AddConnection("FixedImageSource", "RegistrationMethod", { { "NameOfInterface", { "itkImageFixedInterface" } } });
270
  blueprint->AddConnection( "FixedImageSource", "RegistrationMethod", { {} } );
271
272

  //optionally, tie properties to connection to avoid ambiguities
273
  //blueprint->AddConnection("MovingImageSource", "RegistrationMethod", { { "NameOfInterface", { "itkImageMovingInterface" } } });
274
  blueprint->AddConnection( "MovingImageSource", "RegistrationMethod", { {} } );
275
276

  //optionally, tie properties to connection to avoid ambiguities
277
  //blueprint->AddConnection("RegistrationMethod", "ResultImageSink", { { "NameOfInterface", { "itkImageSourceInterface" } } });
278
  blueprint->AddConnection( "ResampleFilter", "ResultImageSink", { {} } );
279
280

  //optionally, tie properties to connection to avoid ambiguities
281
  //blueprint->AddConnection("RegistrationMethod", "ResultDisplacementFieldSink", { { "NameOfInterface", { "DisplacementFieldItkImageSourceInterface" } } });
282
  blueprint->AddConnection( "TransformDisplacementFilter", "ResultDisplacementFieldSink", { {} } );
283
284

  //optionally, tie properties to connection to avoid ambiguities
285
  //blueprint->AddConnection("Metric", "RegistrationMethod", { { "NameOfInterface", { "itkMetricv4Interface" } } });
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
  blueprint->AddConnection( "Metric", "RegistrationMethod", { {} } );

  blueprint->AddConnection( "FixedImageSource", "Transform", { {} } );
  blueprint->AddConnection( "Transform", "RegistrationMethod", { {} } );
  blueprint->AddConnection( "Optimizer", "RegistrationMethod", { {} } );
  blueprint->AddConnection( "RegistrationMethod", "TransformDisplacementFilter", { {} } );
  blueprint->AddConnection( "FixedImageSource", "TransformDisplacementFilter", { {} } );
  blueprint->AddConnection( "RegistrationMethod", "ResampleFilter", { {} } );
  blueprint->AddConnection( "FixedImageSource", "ResampleFilter", { {} } );
  blueprint->AddConnection( "MovingImageSource", "ResampleFilter", { {} } );

  blueprint->AddConnection( "RegistrationMethod", "Controller", { {} } );          //RunRegistrationInterface
  blueprint->AddConnection( "ResampleFilter", "Controller", { {} } );              //ReconnectTransformInterface
  blueprint->AddConnection( "TransformDisplacementFilter", "Controller", { {} } ); //ReconnectTransformInterface
  blueprint->AddConnection( "ResultImageSink", "Controller", { {} } );             //AfterRegistrationInterface
  blueprint->AddConnection( "ResultDisplacementFieldSink", "Controller", { {} } ); //AfterRegistrationInterface
302

303
304
305
  // Data manager provides the paths to the input and output data for unit tests
  DataManagerType::Pointer dataManager = DataManagerType::New();

306
  blueprint->WriteBlueprint( dataManager->GetOutputFile( "itkv4_SVF_MSD.dot" ) );
307
308

  // Instantiate SuperElastix
309
  EXPECT_NO_THROW( superElastixFilter = SuperElastixFilterType::New() );
310
311
312

  // Set up the readers and writers
  ImageReader2DType::Pointer fixedImageReader = ImageReader2DType::New();
313
  fixedImageReader->SetFileName( dataManager->GetInputFile( "coneA2d64.mhd" ) );
314
315

  ImageReader2DType::Pointer movingImageReader = ImageReader2DType::New();
316
  movingImageReader->SetFileName( dataManager->GetInputFile( "coneB2d64.mhd" ) );
317

318
  ImageWriter2DType::Pointer resultImageWriter = ImageWriter2DType::New();
319
  resultImageWriter->SetFileName( dataManager->GetOutputFile( "itkv4_SVF_MSD_Image.mhd" ) );
320

321
  VectorImageWriter2DType::Pointer vectorImageWriter = VectorImageWriter2DType::New();
322
  vectorImageWriter->SetFileName( dataManager->GetOutputFile( "itkv4_SVF_MSD_Displacement.mhd" ) );
323
324

  // Connect SuperElastix in an itk pipeline
325
326
  superElastixFilter->SetInput( "FixedImageSource", fixedImageReader->GetOutput() );
  superElastixFilter->SetInput( "MovingImageSource", movingImageReader->GetOutput() );
Floris Berendsen's avatar
Floris Berendsen committed
327

328
329
  resultImageWriter->SetInput( superElastixFilter->GetOutput< Image2DType >( "ResultImageSink" ) );
  vectorImageWriter->SetInput( superElastixFilter->GetOutput< VectorImage2DType >( "ResultDisplacementFieldSink" ) );
Floris Berendsen's avatar
Floris Berendsen committed
330

331
332
  baselineImageReader->SetFileName( dataManager->GetBaselineFile( "itkv4_SVF_MSD_Image.mhd" ) );
  baselineVectorImageReader->SetFileName( dataManager->GetBaselineFile( "itkv4_SVF_MSD_Displacement.mhd" ) );
Floris Berendsen's avatar
Floris Berendsen committed
333

334
335
336
337
  compareImageFilter->SetTestInput( superElastixFilter->GetOutput< Image2DType >( "ResultImageSink" ) );
  subtractVectorImageFilter->SetInput2( superElastixFilter->GetOutput< VectorImage2DType >( "ResultDisplacementFieldSink" ) );

  EXPECT_NO_THROW( superElastixFilter->SetBlueprint( blueprint ) );
338

339
  //Optional Update call
340
341
  //superElastixFilter->Update();

342
  // Update call on the writers triggers SuperElastix to configure and execute
343
344
  EXPECT_NO_THROW( resultImageWriter->Update() );
  EXPECT_NO_THROW( vectorImageWriter->Update() );
345

346
347
  EXPECT_NO_THROW( compareImageFilter->Update() );
  EXPECT_NO_THROW( statisticsImageFilter->Update() );
348

349
350
351
  EXPECT_EQ( 0, compareImageFilter->GetNumberOfPixelsWithDifferences() );
  EXPECT_LT( statisticsImageFilter->GetSumOutput()->Get(), 1e-16 );
}
Floris Berendsen's avatar
Floris Berendsen committed
352

353
/** Experiment 1a: elastix framework, B-spline transform, normalized correlation metric */
354
TEST_F( WBIRDemoTest, elastix_BS_NCC )
355
{
356
  /** make blueprint configuration */
357
358
  blueprint = Blueprint::New();

359
360
  blueprint->AddComponent( "RegistrationMethod", { { "NameOfClass", { "MonolithicElastixComponent" } },
                                                   { "Transform", { "BSplineTransform" } }, { "Metric", { "AdvancedNormalizedCorrelation" } } } );
361

362
  blueprint->AddComponent( "TransformDisplacementField", { { "NameOfClass", { "MonolithicTransformixComponent" } } } );
363

364
  blueprint->AddComponent( "FixedImageSource", { { "NameOfClass", { "ItkImageSourceFixedComponent" } } } );
365

366
  blueprint->AddComponent( "MovingImageSource", { { "NameOfClass", { "ItkImageSourceMovingComponent" } } } );
367

368
  blueprint->AddComponent( "ResultImageSink", { { "NameOfClass", { "ItkImageSinkComponent" } } } );
369

370
  blueprint->AddComponent( "Controller", { { "NameOfClass", { "RegistrationControllerComponent" } } } );
371
372

  //optionally, tie properties to connection to avoid ambiguities
373
  blueprint->AddConnection( "FixedImageSource", "RegistrationMethod", { {} } ); // {{"NameOfInterface", { "itkImageFixedInterface" }}};
374
375

  //optionally, tie properties to connection to avoid ambiguities
376
377
378
  blueprint->AddConnection( "MovingImageSource", "RegistrationMethod", { {} } ); // {{"NameOfInterface", { "itkImageMovingInterface" }}};

  blueprint->AddConnection( "RegistrationMethod", "TransformDisplacementField", { {} } ); // { { "NameOfInterface", { "elastixTransformParameterObjectInterface" } } } ;
379

380
  blueprint->AddConnection( "FixedImageSource", "TransformDisplacementField", { {} } ); // { { "NameOfInterface", { "itkImageDomainFixedInterface" } } } ;
381

382
  blueprint->AddConnection( "MovingImageSource", "TransformDisplacementField", { {} } ); // { { "NameOfInterface", { "itkImageMovingInterface" } } };
383

384
  blueprint->AddConnection( "TransformDisplacementField", "ResultImageSink", { {} } ); // { { "NameOfInterface", { "itkImageInterface" } } } ;
385

386
  blueprint->AddConnection( "RegistrationMethod", "Controller", { {} } ); // { { "NameOfInterface", { "RunRegistrationInterface" } } } ;
387

388
  blueprint->AddConnection( "TransformDisplacementField", "Controller", { {} } ); // { { "NameOfInterface", { "ReconnectTransformInterface" } } } ;
389

390
  blueprint->AddConnection( "ResultImageSink", "Controller", { {} } ); // { { "NameOfInterface", { "AfterRegistrationInterface" } } } ;
391

392
393
394
395
  // Data manager provides the paths to the input and output data for unit tests
  DataManagerType::Pointer dataManager = DataManagerType::New();

  blueprint->WriteBlueprint(dataManager->GetOutputFile("elastix_BS_NCC.dot"));
396
397

  // Instantiate SuperElastix
398
  EXPECT_NO_THROW( superElastixFilter = SuperElastixFilterType::New() );
399
400
401

  // Set up the readers and writers
  ImageReader2DType::Pointer fixedImageReader = ImageReader2DType::New();
402
  fixedImageReader->SetFileName( dataManager->GetInputFile( "coneA2d64.mhd" ) );
403
404

  ImageReader2DType::Pointer movingImageReader = ImageReader2DType::New();
405
  movingImageReader->SetFileName( dataManager->GetInputFile( "coneB2d64.mhd" ) );
406

407
  ImageWriter2DType::Pointer resultImageWriter = ImageWriter2DType::New();
408
  resultImageWriter->SetFileName( dataManager->GetOutputFile( "elastix_BS_NCC_Image.mhd" ) );
409

410
  VectorImageWriter2DType::Pointer vectorImageWriter = VectorImageWriter2DType::New();
411
  vectorImageWriter->SetFileName( dataManager->GetOutputFile( "elastix_BS_NCC_Displacement.mhd" ) );
412
413

  // Connect SuperElastix in an itk pipeline
414
415
  superElastixFilter->SetInput( "FixedImageSource", fixedImageReader->GetOutput() );
  superElastixFilter->SetInput( "MovingImageSource", movingImageReader->GetOutput() );
416

417
  resultImageWriter->SetInput( superElastixFilter->GetOutput< Image2DType >( "ResultImageSink" ) );
418
  //elastix component does not have an deformation output, but writes deformationfield to disk.
419
420
  //vectorImageWriter->SetInput(superElastixFilter->GetOutput<VectorImage2DType>("ResultDisplacementFieldSink"));

421
422
423
424
  baselineImageReader->SetFileName( dataManager->GetBaselineFile( "elastix_BS_NCC_Image.mhd" ) );
  baselineVectorImageReader->SetFileName( dataManager->GetBaselineFile( "elastix_BS_NCC_Displacement.mhd" ) );

  compareImageFilter->SetTestInput( superElastixFilter->GetOutput< Image2DType >( "ResultImageSink" ) );
425
426
427
  //elastix component does not have an deformation output, but writes deformationfield to disk.
  //subtractVectorImageFilter->SetInput2(superElastixFilter->GetOutput<VectorImage2DType>("ResultDisplacementFieldSink"));
  VectorImageReader2DType::Pointer elastixDeformationFieldReader = VectorImageReader2DType::New();
428
429
430
431
  elastixDeformationFieldReader->SetFileName( "deformationField.nii" );
  subtractVectorImageFilter->SetInput2( elastixDeformationFieldReader->GetOutput() );

  EXPECT_NO_THROW( superElastixFilter->SetBlueprint( blueprint ) );
432
433
434
435
436

  //Optional Update call
  //superElastixFilter->Update();

  // Update call on the writers triggers SuperElastix to configure and execute
437
  EXPECT_NO_THROW( resultImageWriter->Update() );
438
439
  //EXPECT_NO_THROW(vectorImageWriter->Update());

440
441
  EXPECT_NO_THROW( compareImageFilter->Update() );
  EXPECT_NO_THROW( statisticsImageFilter->Update() );
442

443
444
  EXPECT_EQ( 0, compareImageFilter->GetNumberOfPixelsWithDifferences() );
  EXPECT_LT( statisticsImageFilter->GetSumOutput()->Get(), 1e-16 );
445
446
447
}

/** Experiment 1b: elastix framework, B-spline transform, mean squared differences metric */
448
TEST_F( WBIRDemoTest, elastix_BS_MSD )
449
{
450
  /** make blueprint configuration */
451
452
  blueprint = Blueprint::New();

453
454
  blueprint->AddComponent( "RegistrationMethod", { { "NameOfClass", { "MonolithicElastixComponent" } },
                                                   { "Transform", { "BSplineTransform" } }, { "Metric", { "AdvancedMeanSquares" } } } );
455

456
  blueprint->AddComponent( "TransformDisplacementField", { { "NameOfClass", { "MonolithicTransformixComponent" } } } );
457

458
  blueprint->AddComponent( "FixedImageSource", { { "NameOfClass", { "ItkImageSourceFixedComponent" } } } );
459

460
  blueprint->AddComponent( "MovingImageSource", { { "NameOfClass", { "ItkImageSourceMovingComponent" } } } );
461

462
  blueprint->AddComponent( "ResultImageSink", { { "NameOfClass", { "ItkImageSinkComponent" } } } );
463

464
  blueprint->AddComponent( "Controller", { { "NameOfClass", { "RegistrationControllerComponent" } } } );
465
466

  //optionally, tie properties to connection to avoid ambiguities
467
  blueprint->AddConnection( "FixedImageSource", "RegistrationMethod", { {} } ); // {{"NameOfInterface", { "itkImageFixedInterface" }}};
468
469

  //optionally, tie properties to connection to avoid ambiguities
470
  blueprint->AddConnection( "MovingImageSource", "RegistrationMethod", { {} } ); // {{"NameOfInterface", { "itkImageMovingInterface" }}};
471

472
  blueprint->AddConnection( "RegistrationMethod", "TransformDisplacementField", { {} } ); // { { "NameOfInterface", { "elastixTransformParameterObjectInterface" } } } ;
473

474
  blueprint->AddConnection( "FixedImageSource", "TransformDisplacementField", { {} } ); // { { "NameOfInterface", { "itkImageDomainFixedInterface" } } } ;
475

476
  blueprint->AddConnection( "MovingImageSource", "TransformDisplacementField", { {} } ); // { { "NameOfInterface", { "itkImageMovingInterface" } } };
477

478
  blueprint->AddConnection( "TransformDisplacementField", "ResultImageSink", { {} } ); // { { "NameOfInterface", { "itkImageInterface" } } } ;
479

480
  blueprint->AddConnection( "RegistrationMethod", "Controller", { {} } ); // { { "NameOfInterface", { "RunRegistrationInterface" } } } ;
481

482
  blueprint->AddConnection( "TransformDisplacementField", "Controller", { {} } ); // { { "NameOfInterface", { "ReconnectTransformInterface" } } } ;
483

484
  blueprint->AddConnection( "ResultImageSink", "Controller", { {} } ); // { { "NameOfInterface", { "AfterRegistrationInterface" } } } ;
485

486
487
488
489
  // Data manager provides the paths to the input and output data for unit tests
  DataManagerType::Pointer dataManager = DataManagerType::New();

  blueprint->WriteBlueprint(dataManager->GetOutputFile("elastix_BS_MSD.dot"));
490
491

  // Instantiate SuperElastix
492
  EXPECT_NO_THROW( superElastixFilter = SuperElastixFilterType::New() );
493
494
495

  // Set up the readers and writers
  ImageReader2DType::Pointer fixedImageReader = ImageReader2DType::New();
496
  fixedImageReader->SetFileName( dataManager->GetInputFile( "coneA2d64.mhd" ) );
497
498

  ImageReader2DType::Pointer movingImageReader = ImageReader2DType::New();
499
  movingImageReader->SetFileName( dataManager->GetInputFile( "coneB2d64.mhd" ) );
500

501
  ImageWriter2DType::Pointer resultImageWriter = ImageWriter2DType::New();
502
  resultImageWriter->SetFileName( dataManager->GetOutputFile( "elastix_BS_MSD_Image.mhd" ) );
503

504
  VectorImageWriter2DType::Pointer vectorImageWriter = VectorImageWriter2DType::New();
505
  vectorImageWriter->SetFileName( dataManager->GetOutputFile( "elastix_BS_MSD_Displacement.mhd" ) );
506
507

  // Connect SuperElastix in an itk pipeline
508
509
  superElastixFilter->SetInput( "FixedImageSource", fixedImageReader->GetOutput() );
  superElastixFilter->SetInput( "MovingImageSource", movingImageReader->GetOutput() );
510

511
  resultImageWriter->SetInput( superElastixFilter->GetOutput< Image2DType >( "ResultImageSink" ) );
512
513
  //elastix component does not have an deformation output, but writes deformationfield to disk.
  //vectorImageWriter->SetInput(superElastixFilter->GetOutput<VectorImage2DType>("ResultDisplacementFieldSink"));
514

515
516
  baselineImageReader->SetFileName( dataManager->GetBaselineFile( "elastix_BS_MSD_Image.mhd" ) );
  baselineVectorImageReader->SetFileName( dataManager->GetBaselineFile( "elastix_BS_MSD_Displacement.mhd" ) );
517

518
  compareImageFilter->SetTestInput( superElastixFilter->GetOutput< Image2DType >( "ResultImageSink" ) );
519
520
521
  //elastix component does not have an deformation output, but writes deformationfield to disk.
  //subtractVectorImageFilter->SetInput2(superElastixFilter->GetOutput<VectorImage2DType>("ResultDisplacementFieldSink"));
  VectorImageReader2DType::Pointer elastixDeformationFieldReader = VectorImageReader2DType::New();
522
523
  elastixDeformationFieldReader->SetFileName( "deformationField.nii" );
  subtractVectorImageFilter->SetInput2( elastixDeformationFieldReader->GetOutput() );
524

525
  EXPECT_NO_THROW( superElastixFilter->SetBlueprint( blueprint ) );
526
527
528
529
530

  //Optional Update call
  //superElastixFilter->Update();

  // Update call on the writers triggers SuperElastix to configure and execute
531
  EXPECT_NO_THROW( resultImageWriter->Update() );
532
  //EXPECT_NO_THROW(vectorImageWriter->Update());
533

534
535
  EXPECT_NO_THROW( compareImageFilter->Update() );
  EXPECT_NO_THROW( statisticsImageFilter->Update() );
536

537
538
  EXPECT_EQ( 0, compareImageFilter->GetNumberOfPixelsWithDifferences() );
  EXPECT_LT( statisticsImageFilter->GetSumOutput()->Get(), 1e-16 );
539
}
540
} // namespace selx