Slicer  4.11
Slicer is a multi-platform, free and open source software package for visualization and medical image computing
coordTypes.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Extract Skeleton
4  Module: $HeadURL$
5  Language: C++
6  Date: $Date$
7  Version: $Revision$
8 
9  Copyright (c) Brigham and Women's Hospital (BWH) All Rights Reserved.
10 
11  See License.txt or http://www.slicer.org/copyright/copyright.txt for details.
12 
13 ==========================================================================*/
14 
15 #ifndef COORD_TYPES_H
16 #define COORD_TYPES_H
17 
18 #include "math.h"
19 
20 #include "misc.h"
21 
22 class Coord3i
23 {
24  int v[3];
25 public:
27  {
28  v[0] = v[1] = v[2] = -1;
29  }
30  inline int & operator[](int i)
31  {
32  return v[i];
33  }
34  inline void conv(double * i)
35  {
36  i[0] = v[0]; i[1] = v[1]; i[2] = v[2];
37  };
38 };
39 
40 class Coord3f
41 {
42  float v[3];
43 public:
44  inline float & operator[](int i)
45  {
46  return v[i];
47  }
48  inline void conv(float * i)
49  {
50  i[0] = v[0]; i[1] = v[1]; i[2] = v[2];
51  };
52  inline void conv(double * i)
53  {
54  i[0] = v[0]; i[1] = v[1]; i[2] = v[2];
55  };
56 };
57 
58 class Coord3d
59 {
60  double v[3];
61 public:
62  inline double & operator[](int i)
63  {
64  return v[i];
65  };
66  inline void conv(int * i)
67  {
68  i[0] = (int) v[0]; i[1] = (int) v[1]; i[2] = (int) v[2];
69  };
70  inline void conv(float * i)
71  {
72  i[0] = static_cast<float>(v[0]); i[1] = static_cast<float>(v[1]); i[2] = static_cast<float>(v[2]);
73  };
74  inline void conv(double * i)
75  {
76  i[0] = v[0]; i[1] = v[1]; i[2] = v[2];
77  };
78 };
79 
80 // Euclidean distance between two points
81 // TODO: use image pixel spacing
82 inline double pointdistance(Coord3i &p1, Coord3i &p2)
83 {
84  return sqrt((double)sqr(p1[0] - p2[0]) + sqr(p1[1] - p2[1]) + sqr(p1[2] - p2[2]));
85 }
86 
87 inline void normcrossprod(double *v1, double *v2, double *norm)
88 // calculate normalized crossproduct
89 {
90  double absval;
91 
92  norm[0] = v1[1] * v2[2] - v1[2] * v2[1];
93  norm[1] = v1[2] * v2[0] - v1[0] * v2[2];
94  norm[2] = v1[0] * v2[1] - v1[1] * v2[0];
95 
96  absval = sqrt(norm[0] * norm[0] + norm[1] * norm[1] + norm[2] * norm[2]);
97  norm[0] /= absval;
98  norm[1] /= absval;
99  norm[2] /= absval;
100 }
101 
102 // calculate angle between two vectors (0..M_PI), you might want to adjust to 0..M_PI/2
103 // range after call via 'if (angle > M_PI/2) angle = M_PI-angle;'
104 inline double vectorangle(double *v1, double *v2)
105 {
106  double prod = 0, length1 = 0, length2 = 0;
107 
108  for( int k = 0; k < 3; k++ )
109  {
110  prod += v1[k] * v2[k];
111  length1 += v1[k] * v1[k];
112  length2 += v2[k] * v2[k];
113  }
114  return acos(prod / sqrt(length1 * length2) );
115 }
116 
117 // calculate angle between two vectors (0..M_PI), you might want to adjust to 0..M_PI/2
118 // range after call via 'if (angle > M_PI/2) angle = M_PI-angle;'
119 inline double vectorangle(Coord3d v1, Coord3d v2)
120 {
121  double prod = 0, length1 = 0, length2 = 0;
122 
123  for( int k = 0; k < 3; k++ )
124  {
125  prod += v1[k] * v2[k];
126  length1 += v1[k] * v1[k];
127  length2 += v2[k] * v2[k];
128  }
129  return acos(prod / sqrt(length1 * length2) );
130 }
131 
132 inline double vec_length(Coord3d x)
133 {
134  return sqrt(sqr(x[0]) + sqr(x[1]) + sqr(x[2]) );
135 }
136 
137 inline double vec_length(double *x)
138 {
139  return sqrt(sqr(x[0]) + sqr(x[1]) + sqr(x[2]) );
140 }
141 
142 inline double vec_length(double *x, double *y)
143 {
144  return sqrt(sqr(x[0] - y[0]) + sqr(x[1] - y[1]) + sqr(x[2] - y[2]) );
145 }
146 
147 // transform and check index, returns 0 on success and 1 if corrected location
148 inline int transWorldToImage(Coord3d loc_world, int *loc_img,
149  double *origin, int *dims, double voxelsize)
150 {
151  int adjust = 0;
152 
153  for( int i = 0; i < 3; i++ )
154  {
155  loc_img[i] = (int) ( (loc_world[i] - origin[i]) / voxelsize);
156  if( loc_img[i] < 0 )
157  {
158  adjust = 1; loc_img[i] = 0;
159  }
160  if( loc_img[i] >= dims[i] )
161  {
162  loc_img[i] = dims[i] - 1; adjust = 1;
163  }
164  }
165 
166  return adjust;
167 }
168 
169 inline int transWorldToImage(double *loc_world, int *loc_img,
170  double *origin, int *dims, double voxelsize)
171 // transform and check index, returns 0 on success and 1 if corrected location
172 {
173  int adjust = 0;
174 
175  for( int i = 0; i < 3; i++ )
176  {
177  loc_img[i] = (int) ( (loc_world[i] - origin[i]) / voxelsize);
178  if( loc_img[i] < 0 )
179  {
180  adjust = 1; loc_img[i] = 0;
181  }
182  if( loc_img[i] >= dims[i] )
183  {
184  loc_img[i] = dims[i] - 1; adjust = 1;
185  }
186  }
187 
188  return adjust;
189 }
190 
191 #endif
void conv(double *i)
Definition: coordTypes.h:34
void conv(float *i)
Definition: coordTypes.h:48
void conv(int *i)
Definition: coordTypes.h:66
double vectorangle(double *v1, double *v2)
Definition: coordTypes.h:104
double vec_length(Coord3d x)
Definition: coordTypes.h:132
void conv(double *i)
Definition: coordTypes.h:74
double & operator[](int i)
Definition: coordTypes.h:62
void conv(double *i)
Definition: coordTypes.h:52
void normcrossprod(double *v1, double *v2, double *norm)
Definition: coordTypes.h:87
int & operator[](int i)
Definition: coordTypes.h:30
int transWorldToImage(Coord3d loc_world, int *loc_img, double *origin, int *dims, double voxelsize)
Definition: coordTypes.h:148
Coord3i()
Definition: coordTypes.h:26
void conv(float *i)
Definition: coordTypes.h:70
float & operator[](int i)
Definition: coordTypes.h:44
double pointdistance(Coord3i &p1, Coord3i &p2)
Definition: coordTypes.h:82
T sqr(T x)
Definition: misc.h:64