00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 #ifndef __itkPhasedArray3DSpecialCoordinatesImage_h
00018 #define __itkPhasedArray3DSpecialCoordinatesImage_h
00019
00020 #include "itkSpecialCoordinatesImage.h"
00021 #include "itkImageRegion.h"
00022 #include "itkPoint.h"
00023 #include "itkContinuousIndex.h"
00024 #include "vnl/vnl_math.h"
00025
00026
00027 namespace itk
00028 {
00029
00085 template <class TPixel>
00086 class ITK_EXPORT PhasedArray3DSpecialCoordinatesImage :
00087 public SpecialCoordinatesImage<TPixel,3>
00088 {
00089 public:
00091 typedef PhasedArray3DSpecialCoordinatesImage Self;
00092 typedef SpecialCoordinatesImage<TPixel,3> Superclass;
00093 typedef SmartPointer<Self> Pointer;
00094 typedef SmartPointer<const Self> ConstPointer;
00095 typedef WeakPointer<const Self> ConstWeakPointer;
00096
00098 itkNewMacro(Self);
00099
00101 itkTypeMacro(PhasedArray3DSpecialCoordinatesImage, SpecialCoordinatesImage);
00102
00105 typedef TPixel PixelType;
00106
00108 typedef TPixel ValueType ;
00109
00114 typedef TPixel InternalPixelType;
00115
00118 typedef DefaultPixelAccessor< PixelType > AccessorType;
00119
00124 itkStaticConstMacro(ImageDimension, unsigned int, 3);
00125
00127 typedef ImportImageContainer<unsigned long, PixelType> PixelContainer;
00128
00130 typedef typename Superclass::IndexType IndexType;
00131
00133 typedef typename Superclass::OffsetType OffsetType;
00134
00136 typedef typename Superclass::SizeType SizeType;
00137
00139 typedef typename Superclass::RegionType RegionType;
00140
00145 typedef typename Superclass::SpacingType SpacingType;
00146
00149 typedef typename Superclass::PointType PointType;
00150
00152 typedef typename PixelContainer::Pointer PixelContainerPointer;
00153 typedef typename PixelContainer::ConstPointer PixelContainerConstPointer;
00154
00159 template<class TCoordRep>
00160 bool TransformPhysicalPointToContinuousIndex(
00161 const Point<TCoordRep, 3>& point,
00162 ContinuousIndex<TCoordRep, 3>& index ) const
00163 {
00164 RegionType region = this->GetLargestPossibleRegion();
00165 double maxAzimuth = region.GetSize(0) - 1;
00166 double maxElevation = region.GetSize(1) - 1;
00167
00168
00169 TCoordRep azimuth = atan(point[0] / point[2]);
00170 TCoordRep elevation = atan(point[1] / point[2]);
00171 TCoordRep radius = sqrt( point[0] * point[0]
00172 + point[1] * point[1]
00173 + point[2] * point[2] );
00174
00175
00176 index[0] = static_cast<TCoordRep>( (azimuth/m_AzimuthAngularSeparation)
00177 + (maxAzimuth/2.0) );
00178 index[1] = static_cast<TCoordRep>( (elevation/m_ElevationAngularSeparation)
00179 + (maxElevation/2.0) );
00180 index[2] = static_cast<TCoordRep>( ( (radius-m_FirstSampleDistance)
00181 / m_RadiusSampleSize) );
00182
00183
00184 const bool isInside = region.IsInside( index );
00185
00186 return isInside;
00187 }
00188
00193 template<class TCoordRep>
00194 bool TransformPhysicalPointToIndex(
00195 const Point<TCoordRep, 3>& point,
00196 IndexType & index ) const
00197 {
00198 typedef typename IndexType::IndexValueType IndexValueType;
00199
00200 RegionType region = this->GetLargestPossibleRegion();
00201 double maxAzimuth = region.GetSize(0) - 1;
00202 double maxElevation = region.GetSize(1) - 1;
00203
00204
00205 TCoordRep azimuth = atan(point[0] / point[2]);
00206 TCoordRep elevation = atan(point[1] / point[2]);
00207 TCoordRep radius = sqrt( point[0] * point[0]
00208 + point[1] * point[1]
00209 + point[2] * point[2] );
00210
00211
00212 index[0] = static_cast<IndexValueType>( (azimuth/m_AzimuthAngularSeparation)
00213 + (maxAzimuth/2.0) );
00214 index[1] = static_cast<IndexValueType>( (elevation/m_ElevationAngularSeparation)
00215 + (maxElevation/2.0) );
00216 index[2] = static_cast<IndexValueType>( ( (radius-m_FirstSampleDistance)
00217 / m_RadiusSampleSize ) );
00218
00219
00220 const bool isInside = region.IsInside( index );
00221
00222 return isInside;
00223 }
00224
00229 template<class TCoordRep>
00230 void TransformContinuousIndexToPhysicalPoint(
00231 const ContinuousIndex<TCoordRep, 3>& index,
00232 Point<TCoordRep, 3>& point ) const
00233 {
00234 RegionType region = this->GetLargestPossibleRegion();
00235 double maxAzimuth = region.GetSize(0) - 1;
00236 double maxElevation = region.GetSize(1) - 1;
00237
00238
00239 TCoordRep azimuth = ( index[0] - (maxAzimuth/2.0) )
00240 * m_AzimuthAngularSeparation;
00241 TCoordRep elevation = ( index[1] - (maxElevation/2.0) )
00242 * m_ElevationAngularSeparation;
00243 TCoordRep radius = (index[2]*m_RadiusSampleSize)+m_FirstSampleDistance;
00244
00245
00246 TCoordRep tanOfAzimuth = tan(azimuth);
00247 TCoordRep tanOfElevation = tan(elevation);
00248 point[2] = static_cast<TCoordRep>( radius /
00249 sqrt(1 + tanOfAzimuth*tanOfAzimuth + tanOfElevation*tanOfElevation));
00250 point[1] = static_cast<TCoordRep>( point[2] * tanOfElevation );
00251 point[0] = static_cast<TCoordRep>( point[2] * tanOfAzimuth );
00252 }
00253
00259 template<class TCoordRep>
00260 void TransformIndexToPhysicalPoint(
00261 const IndexType & index,
00262 Point<TCoordRep, 3>& point ) const
00263 {
00264 RegionType region = this->GetLargestPossibleRegion();
00265 double maxAzimuth = region.GetSize(0) - 1;
00266 double maxElevation = region.GetSize(1) - 1;
00267
00268
00269 TCoordRep azimuth = ( static_cast<double>(index[0]) - (maxAzimuth/2.0) )
00270 * m_AzimuthAngularSeparation;
00271 TCoordRep elevation = ( static_cast<double>(index[1]) - (maxElevation/2.0) )
00272 * m_ElevationAngularSeparation;
00273 TCoordRep radius = (static_cast<double>(index[2]) * m_RadiusSampleSize)
00274 + m_FirstSampleDistance;
00275
00276
00277 TCoordRep tanOfAzimuth = tan(azimuth);
00278 TCoordRep tanOfElevation = tan(elevation);
00279 point[2] = static_cast<TCoordRep>( radius / sqrt(
00280 1.0 + tanOfAzimuth*tanOfAzimuth + tanOfElevation*tanOfElevation) );
00281 point[1] = static_cast<TCoordRep>( point[2] * tanOfElevation );
00282 point[0] = static_cast<TCoordRep>( point[2] * tanOfAzimuth );
00283 }
00284
00285
00287 itkSetMacro(AzimuthAngularSeparation, double);
00288
00290 itkSetMacro(ElevationAngularSeparation, double);
00291
00293 itkSetMacro(RadiusSampleSize, double);
00294
00296 itkSetMacro(FirstSampleDistance, double);
00297
00298 protected:
00299 PhasedArray3DSpecialCoordinatesImage()
00300 {
00301 m_RadiusSampleSize = 1;
00302 m_AzimuthAngularSeparation = 1 * (2.0*vnl_math::pi/360.0);
00303 m_ElevationAngularSeparation = 1 * (2.0*vnl_math::pi/360.0);
00304 m_FirstSampleDistance = 0;
00305 }
00306 virtual ~PhasedArray3DSpecialCoordinatesImage() {};
00307 void PrintSelf(std::ostream& os, Indent indent) const;
00308
00309 private:
00310 PhasedArray3DSpecialCoordinatesImage(const Self&);
00311 void operator=(const Self&);
00312
00313 double m_AzimuthAngularSeparation;
00314 double m_ElevationAngularSeparation;
00315 double m_RadiusSampleSize;
00316 double m_FirstSampleDistance;
00317
00318 };
00319 #ifdef ITK_EXPLICIT_INSTANTIATION
00320 extern template class PhasedArray3DSpecialCoordinatesImage<float >;
00321 extern template class PhasedArray3DSpecialCoordinatesImage<double >;
00322 extern template class PhasedArray3DSpecialCoordinatesImage<unsigned char >;
00323 extern template class PhasedArray3DSpecialCoordinatesImage<unsigned short>;
00324 extern template class PhasedArray3DSpecialCoordinatesImage<unsigned int >;
00325 extern template class PhasedArray3DSpecialCoordinatesImage<signed char >;
00326 extern template class PhasedArray3DSpecialCoordinatesImage<signed short >;
00327 extern template class PhasedArray3DSpecialCoordinatesImage<signed int >;
00328 #endif
00329 }
00330 #ifndef ITK_MANUAL_INSTANTIATION
00331 #include "itkPhasedArray3DSpecialCoordinatesImage.txx"
00332 #endif
00333
00334 #endif
00335