This file is indexed.

/usr/include/ITK-4.5/itkImageBase.h is in libinsighttoolkit4-dev 4.5.0-3.

This file is owned by root:root, with mode 0o644.

The actual contents of the file can be viewed below.

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
/*=========================================================================
 *
 *  Copyright Insight Software Consortium
 *
 *  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.
 *
 *=========================================================================*/
/*=========================================================================
 *
 *  Portions of this file are subject to the VTK Toolkit Version 3 copyright.
 *
 *  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
 *
 *  For complete copyright, license and disclaimer of warranty information
 *  please refer to the NOTICE file at the top of the ITK source tree.
 *
 *=========================================================================*/
#ifndef __itkImageBase_h
#define __itkImageBase_h

#include "itkDataObject.h"

#include "itkImageRegion.h"
#include "itkMatrix.h"
#include "itkObjectFactory.h"
#include "itkOffset.h"
#include "itkFixedArray.h"
#include "itkImageHelper.h"
#include "itkFloatTypes.h"

//HACK:  vnl/vnl_matrix_fixed.txx is needed here?
//      to avoid undefined symbol vnl_matrix_fixed<double, 8u, 8u>::set_identity()", referenced from
#include "vnl/vnl_matrix_fixed.txx"

#include "itkImageTransformHelper.h"

namespace itk
{

/* Forward declaration (ImageTransformHelper include's ImageBase) */
template< unsigned int NImageDimension, unsigned int R, unsigned int C, typename TPointValue, typename TMatrixValue >
class ImageTransformHelper;

/** \class ImageBase
 * \brief Base class for templated image classes.
 *
 * ImageBase is the base class for the templated Image
 * classes. ImageBase is templated over the dimension of the image. It
 * provides the API and ivars that depend solely on the dimension of
 * the image.  ImageBase does not store any of the image (pixel) data.
 * Storage for the pixels and the pixel access methods are defined in
 * subclasses of ImageBase, namely Image and ImageAdaptor.
 *
 * ImageBase manages the geometry of an image. The geometry of an
 * image is defined by its position, orientation, spacing, and extent.
 *
 * The position and orientation of an image is defined by its "Origin"
 * and its "Directions".  The "Origin" is the physical position of the
 * pixel whose "Index" is all zeros. The "Direction" of an image is a
 * matrix whose columns indicate the direction in physical space that
 * each dimension of the image traverses. The first column defines the
 * direction that the fastest moving index in the image traverses in
 * physical space while the last column defines the direction that the
 * slowest moving index in the image traverses in physical space.
 *
 * The extent of an image is defined by the pixel spacing and a set of
 * regions. The "Spacing" is the size of a pixel in physical space
 * along each dimension.  Regions describe a portion of an image grid
 * via a starting index for the image array and a size (or number of
 * pixels) in each dimension. The ivar LargestPossibleRegion defines
 * the size and starting index of the image dataset. The entire image
 * dataset, however, may not be resident in memory. The region of the
 * image that is resident in memory is defined by the
 * "BufferedRegion". The Buffer is a contiguous block of memory.  The
 * third set of meta-data defines a region of interest, called the
 * "RequestedRegion". The RequestedRegion is used by the pipeline
 * execution model to define what a filter is requested to produce.
 *
 * [RegionIndex, RegionSize] C [BufferIndex, BufferSize]
 *                           C [ImageIndex, ImageSize]
 *
 * ImageBase provides all the methods for converting between the
 * physical space and index coordinate
 * frames. TransformIndexToPhysicalPoint() converts an Index in the
 * pixel array into its coordinates in physical space.
 * TransformPhysicalPointToIndex() converts a position in physical
 * space into an Index into the pixel array (using
 * rounding). Subpixel locations are supported by methods that
 * convert to and from ContinuousIndex types.
 *
 * ImageBase also provides helper routines for the ImageIterators
 * which convert an Index to an offset in memory from the first pixel
 * address as well as covert an offset in memory from the first pixel
 * address to an Index.
 *
 * \ingroup ImageObjects
 * \ingroup ITKSystemObjects
 *
 * \ingroup ITKCommon
 */
template< unsigned int VImageDimension = 2 >
class ImageBase:public DataObject
{
public:
  /** Standard typedefs. */
  typedef ImageBase                  Self;
  typedef DataObject                 Superclass;
  typedef SmartPointer< Self >       Pointer;
  typedef SmartPointer< const Self > ConstPointer;

  /** Method for creation through the object factory. */
  itkNewMacro(Self);

  /** Run-time type information (and related methods). */
  itkTypeMacro(ImageBase, DataObject);

  /** Dimension of the image.  This constant is used by functions that are
   * templated over image type (as opposed to being templated over pixel
   * type and dimension) when they need compile time access to the dimension
   * of the image. */
  itkStaticConstMacro(ImageDimension, unsigned int, VImageDimension);

  /** Index typedef support. An index is used to access pixel values. */
  typedef Index< VImageDimension >           IndexType;
  typedef typename IndexType::IndexValueType IndexValueType;

  /** Offset typedef support. An offset represent relative position
   * between indices. */
  typedef Offset< VImageDimension >            OffsetType;
  typedef typename OffsetType::OffsetValueType OffsetValueType;

  /** Size typedef support. A size is used to define region bounds. */
  typedef Size< VImageDimension >          SizeType;
  typedef typename SizeType::SizeValueType SizeValueType;

  /** Region typedef support. A region is used to specify a subset of an image. */
  typedef ImageRegion< VImageDimension > RegionType;

  /** Spacing typedef support.  Spacing holds the size of a pixel.
   * The spacing is the geometric distance between image samples along
   * each dimension. ITK only supports positive spacing value:
   * negative values may cause undesirable results.  */
  typedef SpacePrecisionType                          SpacingValueType;
  typedef Vector< SpacingValueType, VImageDimension > SpacingType;

  /** Origin typedef support.  The origin is the geometric coordinates
   * of the index (0,0). */
  typedef SpacePrecisionType                       PointValueType;
  typedef Point< PointValueType, VImageDimension > PointType;

  /** Direction typedef support.  The Direction is a matix of
   * direction cosines that specify the direction in physical space
   * between samples along each dimension. */
  typedef Matrix< SpacePrecisionType, VImageDimension, VImageDimension > DirectionType;

  /** Restore object to initialized state. */
  void Initialize();

  /** Image dimension. The dimension of an image is fixed at construction. */
  static unsigned int GetImageDimension()
  { return VImageDimension; }

  /** Set the origin of the image. The origin is the geometric
   * coordinates of the image origin (pixel [0,0]).  It is stored internally
   * as SpacePrecisionType but may be set from float or double.
   * \sa GetOrigin() */
  itkSetMacro(Origin, PointType);
  virtual void SetOrigin(const double origin[VImageDimension]);
  virtual void SetOrigin(const float origin[VImageDimension]);

  /** Set the direction cosines of the image. The direction cosines
   * are vectors that point from one pixel to the next.
   *
   * Each column of the matrix indicates the direction cosines of the unit vector
   * that is parallel to the lines of the image grid corresponding to that
   * dimension. For example, an image with Direction matrix
   *
   *    0.866   0.500
   *   -0.500   0.866
   *
   * has an image grid were the fastest changing index (dimension[0]) walks
   * over a line that in physical space is oriented parallel to the vector
   * (0.866, -0.5). The second fastest changing index (dimension[1]) walks along
   * a line that in Physical space is oriented parallel to the vector
   * (0.5, 0.866)
   *
   * The columns of the Direction matrix are expected to form an
   * orthogonal right handed coordinate syste.  But this is not
   * checked nor enforced in itk::ImageBase.
   *
   * For details, please see:
   *
   * http://www.itk.org/Wiki/Proposals:Orientation#Some_notes_on_the_DICOM_convention_and_current_ITK_usage
   *
   * \sa GetDirection() */
  virtual void SetDirection(const DirectionType & direction);

  /** Get the direction cosines of the image. The direction cosines
   * are vectors that point from one pixel to the next.
   * For ImageBase and Image, the default direction is identity. */
  itkGetConstReferenceMacro(Direction, DirectionType);

  /** Get the inverse direction cosines of the image.
   * These are calculated automatically in SetDirection, thus there
   * is no Set accessor. */
  itkGetConstReferenceMacro(InverseDirection, DirectionType);

  /** Get the spacing (size of a pixel) `of the image. The
   * spacing is the geometric distance between image samples along
   * each dimension. The value returned is a Vector<double, VImageDimension>.
   * For ImageBase and Image, the default data spacing is unity. */
  itkGetConstReferenceMacro(Spacing, SpacingType);

  /** Get the origin of the image. The origin is the geometric
   * coordinates of the index (0,0).  The value returned is a
   * Point<double, VImageDimension>. For ImageBase and Image, the
   * default origin is 0. */
  itkGetConstReferenceMacro(Origin, PointType);

  /** Allocate the image memory. The size of the image must
   * already be set, e.g. by calling SetRegions() or SetBufferedRegion().
   *
   * This method should be pure virtual, if backwards compatibility
   * was not required.
   */
  virtual void Allocate() {}

  /** Set the region object that defines the size and starting index
   * for the largest possible region this image could represent.  This
   * is used in determining how much memory would be needed to load an
   * entire dataset.  It is also used to determine boundary true
   * conditions.
   * \sa ImageRegion, SetBufferedRegion(), SetRequestedRegion() */
  virtual void SetLargestPossibleRegion(const RegionType & region);

  /** Get the region object that defines the size and starting index
   * for the largest possible region this image could represent.  This
   * is used in determining how much memory would be needed to load an
   * entire dataset.  It is also used to determine boundary true
   * conditions.
   * \sa ImageRegion, GetBufferedRegion(), GetRequestedRegion() */
  virtual const RegionType & GetLargestPossibleRegion() const
  { return m_LargestPossibleRegion; }

  /** Set the region object that defines the size and starting index
   * of the region of the image currently loaded in memory.
   * \sa ImageRegion, SetLargestPossibleRegion(), SetRequestedRegion() */
  virtual void SetBufferedRegion(const RegionType & region);

  /** Get the region object that defines the size and starting index
   * of the region of the image currently loaded in memory.
   * \sa ImageRegion, SetLargestPossibleRegion(), SetRequestedRegion() */
  virtual const RegionType & GetBufferedRegion() const
  { return m_BufferedRegion; }

  /** Set the region object that defines the size and starting index
   * for the region of the image requested (i.e., the region of the
   * image to be operated on by a filter). Setting the RequestedRegion
   * does not cause the object to be modified. This method is called
   * internally by the pipeline and therefore bypasses the modified
   * time calculation.
   * \sa ImageRegion, SetLargestPossibleRegion(), SetBufferedRegion() */
  virtual void SetRequestedRegion(const RegionType & region);

  /** Set the requested region from this data object to match the requested
   * region of the data object passed in as a parameter.  This method
   * implements the API from DataObject. The data object parameter must be
   * castable to an ImageBase. Setting the RequestedRegion does not cause
   * the object to be modified. This method is called internally by
   * the pipeline and therefore bypasses the modified time
   * calculation. */
  virtual void SetRequestedRegion( const DataObject *data );

  /** Get the region object that defines the size and starting index
   * for the region of the image requested (i.e., the region of the
   * image to be operated on by a filter).
   * \sa ImageRegion, SetLargestPossibleRegion(), SetBufferedRegion() */
  virtual const RegionType & GetRequestedRegion() const
  { return m_RequestedRegion; }

  /** Convenience methods to set the LargestPossibleRegion,
   *  BufferedRegion and RequestedRegion. Allocate must still be called.
   */
  virtual void SetRegions(const RegionType& region)
  {
    this->SetLargestPossibleRegion(region);
    this->SetBufferedRegion(region);
    this->SetRequestedRegion(region);
  }

  virtual void SetRegions(const SizeType& size)
  {
    RegionType region; region.SetSize(size);

    this->SetLargestPossibleRegion(region);
    this->SetBufferedRegion(region);
    this->SetRequestedRegion(region);
  }

  /** Get the offset table.  The offset table gives increments for
   * moving from one pixel to next in the current row, column, slice,
   * etc..  This table if of size [VImageDimension+1], because its
   * values are computed progressively as: {1, N1, N1*N2,
   * N1*N2*N3,...,(N1*...*Nn)} Where the values {N1,...,Nn} are the
   * elements of the BufferedRegion::Size array.  The last element of
   * the OffsetTable is equivalent to the BufferSize.  Having a
   * [VImageDimension+1] size array, simplifies the implementation of
   * some data accessing algorithms. The entries in the offset table
   * are only valid after the BufferedRegion is set. */
  const OffsetValueType * GetOffsetTable() const { return m_OffsetTable; }

  /** Compute an offset from the beginning of the buffer for a pixel
   * at the specified index. The index is not checked as to whether it
   * is inside the current buffer, so the computed offset could
   * conceivably be outside the buffer. If bounds checking is needed,
   * one can call ImageRegion::IsInside(ind) on the BufferedRegion
   * prior to calling ComputeOffset. */
  inline OffsetValueType ComputeOffset(const IndexType & ind) const
  {
    OffsetValueType offset = 0;

    ImageHelper< VImageDimension, VImageDimension >::ComputeOffset(this->GetBufferedRegion().GetIndex(),
                                                                   ind,
                                                                   m_OffsetTable,
                                                                   offset);
    return offset;
    /* NON TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING data version
     * Leaving here for documentation purposes
     * OffsetValueType ComputeOffset(const IndexType & ind) const
     * {
     *   // need to add bounds checking for the region/buffer?
     *   OffsetValueType   offset = 0;
     *   const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
     *   // data is arranged as [][][][slice][row][col]
     *   // with Index[0] = col, Index[1] = row, Index[2] = slice
     *   for ( int i = VImageDimension - 1; i > 0; i-- )
     *     {
     *     offset += ( ind[i] - bufferedRegionIndex[i] ) * m_OffsetTable[i];
     *     }
     *   offset += ( ind[0] - bufferedRegionIndex[0] );
     *   return offset;
     * }
     */
  }

  /** Compute the index of the pixel at a specified offset from the
   * beginning of the buffered region. Bounds checking is not
   * performed. Thus, the computed index could be outside the
   * BufferedRegion. To ensure a valid index, the parameter "offset"
   * should be between 0 and the number of pixels in the
   * BufferedRegion (the latter can be found using
   * ImageRegion::GetNumberOfPixels()). */
  inline IndexType ComputeIndex(OffsetValueType offset) const
  {
    IndexType         index;
    const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();

    ImageHelper< VImageDimension, VImageDimension >::ComputeIndex(bufferedRegionIndex,
                                                                  offset,
                                                                  m_OffsetTable,
                                                                  index);
    return index;
    /* NON TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING data version
     * Leaving here for documentation purposes
     * IndexType ComputeIndex(OffsetValueType offset) const
     * {
     *   IndexType         index;
     *   const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
     *   for ( int i = VImageDimension - 1; i > 0; i-- )
     *     {
     *     index[i] = static_cast< IndexValueType >( offset / m_OffsetTable[i] );
     *     offset -= ( index[i] * m_OffsetTable[i] );
     *     index[i] += bufferedRegionIndex[i];
     *     }
     *   index[0] = bufferedRegionIndex[0] + static_cast< IndexValueType >( offset );
     *   return index;
     * }
    */

  }

  /** Set the spacing (size of a pixel) of the image. The spacing is
   * the geometric distance between image samples along each
   * dimension. It is stored internally as double, but may be set from
   * float. These methods also pre-compute the Index to Physical point
   * transforms of the image.
   * \sa GetSpacing() */
  virtual void SetSpacing(const SpacingType & spacing);
  virtual void SetSpacing(const double spacing[VImageDimension]);
  virtual void SetSpacing(const float spacing[VImageDimension]);

  /** Get the index (discrete) of a voxel from a physical point.
   * Floating point index results are rounded to integers
   * Returns true if the resulting index is within the image, false otherwise
   * \sa Transform */
  template< typename TCoordRep >
  bool TransformPhysicalPointToIndex(
    const Point< TCoordRep, VImageDimension > & point,
    IndexType & index) const
  {
    ImageTransformHelper< VImageDimension,VImageDimension - 1, VImageDimension - 1, TCoordRep, SpacePrecisionType >
      ::TransformPhysicalPointToIndex(this->m_PhysicalPointToIndex, this->m_Origin, point, index);

    // Now, check to see if the index is within allowed bounds
    const bool isInside = this->GetLargestPossibleRegion().IsInside(index);
    return isInside;
    /* NON TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING data version
     * Leaving here for documentation purposes
     * template< typename TCoordRep >
     * bool TransformPhysicalPointToIndex(
     *   const Point< TCoordRep, VImageDimension > & point,
     *   IndexType & index) const
     * {
     *   for ( unsigned int i = 0; i < VImageDimension; i++ )
     *     {
     *     TCoordRep sum = NumericTraits< TCoordRep >::Zero;
     *     for ( unsigned int j = 0; j < VImageDimension; j++ )
     *       {
     *       sum += this->m_PhysicalPointToIndex[i][j] * ( point[j] - this->m_Origin[j] );
     *       }
     *     index[i] = Math::RoundHalfIntegerUp< IndexValueType >(sum);
     *     }
     *   // Now, check to see if the index is within allowed bounds
     *   const bool isInside = this->GetLargestPossibleRegion().IsInside(index);
     *   return isInside;
     * }
     */
  }

  /** \brief Get the continuous index from a physical point
   *
   * Returns true if the resulting index is within the image, false otherwise.
   * \sa Transform */
  template< typename TCoordRep, typename TIndexRep >
  bool TransformPhysicalPointToContinuousIndex(
    const Point< TCoordRep, VImageDimension > & point,
    ContinuousIndex< TIndexRep, VImageDimension > & index) const
  {
    Vector< SpacePrecisionType, VImageDimension > cvector;

    for ( unsigned int k = 0; k < VImageDimension; k++ )
      {
      cvector[k] = point[k] - this->m_Origin[k];
      }
    cvector = m_PhysicalPointToIndex * cvector;
    for ( unsigned int i = 0; i < VImageDimension; i++ )
      {
      index[i] = static_cast< TIndexRep >( cvector[i] );
      }

    // Now, check to see if the index is within allowed bounds
    const bool isInside = this->GetLargestPossibleRegion().IsInside(index);

    return isInside;
  }

  /** Get a physical point (in the space which
   * the origin and spacing information comes from)
   * from a continuous index (in the index space)
   * \sa Transform */
  template< typename TCoordRep, typename TIndexRep >
  void TransformContinuousIndexToPhysicalPoint(
    const ContinuousIndex< TIndexRep, VImageDimension > & index,
    Point< TCoordRep, VImageDimension > & point) const
  {
    for ( unsigned int r = 0; r < VImageDimension; r++ )
      {
      TCoordRep sum = NumericTraits< TCoordRep >::Zero;
      for ( unsigned int c = 0; c < VImageDimension; c++ )
        {
        sum += this->m_IndexToPhysicalPoint(r, c) * index[c];
        }
      point[r] = sum + this->m_Origin[r];
      }
  }

  /** Get a physical point (in the space which
   * the origin and spacing information comes from)
   * from a discrete index (in the index space)
   *
   * \sa Transform */
  template< typename TCoordRep >
  void TransformIndexToPhysicalPoint(
    const IndexType & index,
    Point< TCoordRep, VImageDimension > & point) const
  {
    ImageTransformHelper< VImageDimension, VImageDimension - 1, VImageDimension - 1,TCoordRep, SpacePrecisionType >::
      TransformIndexToPhysicalPoint(this->m_IndexToPhysicalPoint, this->m_Origin, index, point);
    /* NON TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING data version
     * Leaving here for documentation purposes
     * template< typename TCoordRep >
     * void TransformIndexToPhysicalPoint(
     *   const IndexType & index,
     *   Point< TCoordRep, VImageDimension > & point) const
     * {
     *   for ( unsigned int i = 0; i < VImageDimension; i++ )
     *     {
     *     point[i] = this->m_Origin[i];
     *     for ( unsigned int j = 0; j < VImageDimension; j++ )
     *       {
     *       point[i] += m_IndexToPhysicalPoint[i][j] * index[j];
     *       }
     *     }
     * }
     */
  }

  /** Take a vector or covariant vector that has been computed in the
   * coordinate system parallel to the image grid and rotate it by the
   * direction cosines in order to get it in terms of the coordinate system of
   * the image acquisition device.  This implementation in the Image
   * multiply the array (vector or covariant vector) by the matrix of Direction
   * Cosines. The arguments of the method are of type FixedArray to make
   * possible to use this method with both Vector and CovariantVector.
   * The Method is implemented differently in the itk::Image.
   *
   * \sa Image
   */
  template< typename TCoordRep >
  void TransformLocalVectorToPhysicalVector(
    const FixedArray< TCoordRep, VImageDimension > & inputGradient,
    FixedArray< TCoordRep, VImageDimension > & outputGradient) const
  {
    //
    //TODO: This temporary implementation should be replaced with Template
    // MetaProgramming.
    //
    const DirectionType & direction = this->GetDirection();

    for ( unsigned int i = 0; i < VImageDimension; i++ )
      {
      typedef typename NumericTraits< TCoordRep >::AccumulateType CoordSumType;
      CoordSumType sum = NumericTraits< CoordSumType >::Zero;
      for ( unsigned int j = 0; j < VImageDimension; j++ )
        {
        sum += direction[i][j] * inputGradient[j];
        }
      outputGradient[i] = static_cast< TCoordRep >( sum );
      }
  }

  /** Take a vector or covariant vector that has been computed in terms of the
   * coordinate system of the image acquisition device, and rotate it by the
   * inverse direction cosines in order to get it in the coordinate system
   * parallel to the image grid. This implementation in the Image
   * multiply the array (vector or covariant vector) by the inverse of Direction
   * Cosines. The arguments of the method are of type FixedArray to make
   * possible to use this method with both Vector and CovariantVector.
   */
  template< typename TCoordRep >
  void TransformPhysicalVectorToLocalVector(
    const FixedArray< TCoordRep, VImageDimension > & inputGradient,
    FixedArray< TCoordRep, VImageDimension > & outputGradient) const
  {
    //
    //TODO: This temporary implementation should be replaced with Template
    // MetaProgramming.
    //
    const DirectionType & inverseDirection = this->GetInverseDirection();

    for ( unsigned int i = 0; i < VImageDimension; i++ )
      {
      typedef typename NumericTraits< TCoordRep >::AccumulateType CoordSumType;
      CoordSumType sum = NumericTraits< CoordSumType >::Zero;
      for ( unsigned int j = 0; j < VImageDimension; j++ )
        {
        sum += inverseDirection[i][j] * inputGradient[j];
        }
      outputGradient[i] = static_cast< TCoordRep >( sum );
      }
  }

  /** Copy information from the specified data set.  This method is
   * part of the pipeline execution model. By default, a ProcessObject
   * will copy meta-data from the first input to all of its
   * outputs. See ProcessObject::GenerateOutputInformation().  Each
   * subclass of DataObject is responsible for being able to copy
   * whatever meta-data it needs from from another DataObject.
   * ImageBase has more meta-data than its DataObject.  Thus, it must
   * provide its own version of CopyInformation() in order to copy the
   * LargestPossibleRegion from the input parameter. */
  virtual void CopyInformation(const DataObject *data);

  /** Graft the data and information from one image to another. This
   * is a convenience method to setup a second image with all the meta
   * information of another image and use the same pixel
   * container. Note that this method is different than just using two
   * SmartPointers to the same image since separate DataObjects are
   * still maintained. This method is similar to
   * ImageSource::GraftOutput(). The implementation in ImageBase
   * simply calls CopyInformation() and copies the region ivars.
   * Subclasses of ImageBase are responsible for copying the pixel
   * container. */
  virtual void Graft(const DataObject *data);

  /** Update the information for this DataObject so that it can be used
   * as an output of a ProcessObject.  This method is used the pipeline
   * mechanism to propagate information and initialize the meta data
   * associated with a DataObject. This method calls its source's
   * ProcessObject::UpdateOutputInformation() which determines modified
   * times, LargestPossibleRegions, and any extra meta data like spacing,
   * origin, etc. */
  virtual void UpdateOutputInformation();

  /** UpdateOutputData() is part of the pipeline infrastructure to
   * communicate between ProcessObjects and DataObjects. The method of
   * the superclass is overriden to check if the requested image
   * region has zero pixels. This is needed so that filters can set an
   * input's requested region to zero, to indicate that it does not
   * need to be updated or executed.
   */
  virtual void UpdateOutputData();

  /** Set the RequestedRegion to the LargestPossibleRegion.  This
   * forces a filter to produce all of the output in one execution
   * (i.e. not streaming) on the next call to Update(). */
  virtual void SetRequestedRegionToLargestPossibleRegion();

  /** Determine whether the RequestedRegion is outside of the
   * BufferedRegion. This method returns true if the RequestedRegion
   * is outside the BufferedRegion (true if at least one pixel is
   * outside). This is used by the pipeline mechanism to determine
   * whether a filter needs to re-execute in order to satisfy the
   * current request.  If the current RequestedRegion is already
   * inside the BufferedRegion from the previous execution (and the
   * current filter is up to date), then a given filter does not need
   * to re-execute */
  virtual bool RequestedRegionIsOutsideOfTheBufferedRegion();

  /** Verify that the RequestedRegion is within the
   * LargestPossibleRegion.  If the RequestedRegion is not within the
   * LargestPossibleRegion, then the filter cannot possible satisfy
   * the request. This method returns true if the request can be
   * satisfied and returns fails if the request cannot. This method is
   * used by PropagateRequestedRegion().  PropagateRequestedRegion()
   * throws a InvalidRequestedRegionError exception is the requested
   * region is not within the LargestPossibleRegion. */
  virtual bool VerifyRequestedRegion();

  /** INTERNAL This method is used internally by filters to copy meta-data from
   * the output to the input. Users should not have a need to use this method.
   *
   * Filters that override the ProcessObject's GenerateOutputInformation()
   * should generally have the following line if they want to propagate meta-
   * data for both Image and VectorImage
   * \code
   * outputImage->SetNumberOfComponentsPerPixel(
   *    inputImage->GetNumberOfComponentsPerPixel() )
   * \endcode
   *
   * \sa ImageBase, VectorImage
   *
   * Returns/Sets the number of components in the image. Note that in the
   * ImageBase implementation, this always returns 1. Image returns the
   * # returned from NumericTraits for the pixel type, and VectorImage
   * returns the vector length set by the user.
   */
  virtual unsigned int GetNumberOfComponentsPerPixel() const;
  virtual void SetNumberOfComponentsPerPixel(unsigned int);

protected:
  ImageBase();
  ~ImageBase();
  virtual void PrintSelf(std::ostream & os, Indent indent) const;

  /** Calculate the offsets needed to move from one pixel to the next
   * along a row, column, slice, volume, etc. These offsets are based
   * on the size of the BufferedRegion. This should be called after
   * the BufferedRegion is set. */
  void ComputeOffsetTable();

  /** Compute helper matrices used to transform Index coordinates to
   * PhysicalPoint coordinates and back. This method is virtual and will be
   * overloaded in derived classes in order to provide backward compatibility
   * behavior in classes that did not used to take image orientation into
   * account.  */
  virtual void ComputeIndexToPhysicalPointMatrices();

protected:
  /** Origin, spacing, and direction in physical coordinates. This variables are
   * protected for efficiency.  They are referenced frequently by
   * inner loop calculations. */
  SpacingType   m_Spacing;
  PointType     m_Origin;
  DirectionType m_Direction;
  DirectionType m_InverseDirection;

  /** Matrices intended to help with the conversion of Index coordinates
   *  to PhysicalPoint coordinates */
  DirectionType m_IndexToPhysicalPoint;
  DirectionType m_PhysicalPointToIndex;

  /** Restores the buffered region to it's default state
   *  This method does not call Modify because Initialization is
   *  called by ReleaseData and can not modify the MTime
   * \sa  ReleaseData, Initialize, SetBufferedRegion */
  virtual void InitializeBufferedRegion(void);

private:
  ImageBase(const Self &);      //purposely not implemented
  void operator=(const Self &); //purposely not implemented

  void InternalSetSpacing(const SpacingValueType spacing[VImageDimension])
    {
      SpacingType s(spacing);
      this->SetSpacing(s);
    }

  template <typename TSpacingValue>
  void InternalSetSpacing(const TSpacingValue spacing[VImageDimension])
    {
      Vector<TSpacingValue,VImageDimension> sf(spacing);
      SpacingType                           s;
      s.CastFrom(sf);
      this->SetSpacing(s);
    }

  OffsetValueType m_OffsetTable[VImageDimension + 1];

  RegionType m_LargestPossibleRegion;
  RegionType m_RequestedRegion;
  RegionType m_BufferedRegion;
};
} // end namespace itk

#ifndef ITK_MANUAL_INSTANTIATION
#include "itkImageBase.hxx"
#endif

#endif