1 : /******************************************************************************
2 : * $Id: ntfstroke.cpp 10645 2007-01-18 02:22:39Z warmerdam $
3 : *
4 : * Project: NTF Translator
5 : * Purpose: NTF Arc to polyline stroking code. This code is really generic,
6 : * and might be moved into an OGR module at some point in the
7 : * future.
8 : * Author: Frank Warmerdam, warmerdam@pobox.com
9 : *
10 : ******************************************************************************
11 : * Copyright (c) 2001, Frank Warmerdam
12 : *
13 : * Permission is hereby granted, free of charge, to any person obtaining a
14 : * copy of this software and associated documentation files (the "Software"),
15 : * to deal in the Software without restriction, including without limitation
16 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
17 : * and/or sell copies of the Software, and to permit persons to whom the
18 : * Software is furnished to do so, subject to the following conditions:
19 : *
20 : * The above copyright notice and this permission notice shall be included
21 : * in all copies or substantial portions of the Software.
22 : *
23 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29 : * DEALINGS IN THE SOFTWARE.
30 : ****************************************************************************/
31 :
32 : #include <stdarg.h>
33 : #include "ntf.h"
34 : #include "cpl_conv.h"
35 : #include "cpl_string.h"
36 :
37 : CPL_CVSID("$Id: ntfstroke.cpp 10645 2007-01-18 02:22:39Z warmerdam $");
38 :
39 : #ifndef PI
40 : #define PI 3.14159265358979323846
41 : #endif
42 :
43 : /************************************************************************/
44 : /* NTFArcCenterFromEdgePoints() */
45 : /* */
46 : /* Compute the center of an arc/circle from three edge points. */
47 : /************************************************************************/
48 :
49 0 : int NTFArcCenterFromEdgePoints( double x_c0, double y_c0,
50 : double x_c1, double y_c1,
51 : double x_c2, double y_c2,
52 : double *x_center, double *y_center )
53 :
54 : {
55 :
56 : /* -------------------------------------------------------------------- */
57 : /* Handle a degenerate case that occurs in OSNI products by */
58 : /* making some assumptions. If the first and third points are */
59 : /* the same assume they are intended to define a full circle, */
60 : /* and that the second point is on the opposite side of the */
61 : /* circle. */
62 : /* -------------------------------------------------------------------- */
63 0 : if( x_c0 == x_c2 && y_c0 == y_c2 )
64 : {
65 0 : *x_center = (x_c0 + x_c1) * 0.5;
66 0 : *y_center = (y_c0 + y_c1) * 0.5;
67 :
68 0 : return TRUE;
69 : }
70 :
71 : /* -------------------------------------------------------------------- */
72 : /* Compute the inverse of the slopes connecting the first and */
73 : /* second points. Also compute the center point of the two */
74 : /* lines ... the point our crossing line will go through. */
75 : /* -------------------------------------------------------------------- */
76 : double m1, x1, y1;
77 :
78 0 : if( (y_c1 - y_c0) != 0.0 )
79 0 : m1 = (x_c0 - x_c1) / (y_c1 - y_c0);
80 : else
81 0 : m1 = 1e+10;
82 :
83 0 : x1 = (x_c0 + x_c1) * 0.5;
84 0 : y1 = (y_c0 + y_c1) * 0.5;
85 :
86 : /* -------------------------------------------------------------------- */
87 : /* Compute the same for the second point compared to the third */
88 : /* point. */
89 : /* -------------------------------------------------------------------- */
90 : double m2, x2, y2;
91 :
92 0 : if( (y_c2 - y_c1) != 0.0 )
93 0 : m2 = (x_c1 - x_c2) / (y_c2 - y_c1);
94 : else
95 0 : m2 = 1e+10;
96 :
97 0 : x2 = (x_c1 + x_c2) * 0.5;
98 0 : y2 = (y_c1 + y_c2) * 0.5;
99 :
100 : /* -------------------------------------------------------------------- */
101 : /* Turn these into the Ax+By+C = 0 form of the lines. */
102 : /* -------------------------------------------------------------------- */
103 : double a1, a2, b1, b2, c1, c2;
104 :
105 0 : a1 = m1;
106 0 : a2 = m2;
107 :
108 0 : b1 = -1.0;
109 0 : b2 = -1.0;
110 :
111 0 : c1 = (y1 - m1*x1);
112 0 : c2 = (y2 - m2*x2);
113 :
114 : /* -------------------------------------------------------------------- */
115 : /* Compute the intersection of the two lines through the center */
116 : /* of the circle, using Kramers rule. */
117 : /* -------------------------------------------------------------------- */
118 : double det_inv;
119 :
120 0 : if( a1*b2 - a2*b1 == 0.0 )
121 0 : return FALSE;
122 :
123 0 : det_inv = 1 / (a1*b2 - a2*b1);
124 :
125 0 : *x_center = (b1*c2 - b2*c1) * det_inv;
126 0 : *y_center = (a2*c1 - a1*c2) * det_inv;
127 :
128 0 : return TRUE;
129 : }
130 :
131 : /************************************************************************/
132 : /* NTFStrokeArcToOGRGeometry_Points() */
133 : /************************************************************************/
134 :
135 : OGRGeometry *
136 0 : NTFStrokeArcToOGRGeometry_Points( double dfStartX, double dfStartY,
137 : double dfAlongX, double dfAlongY,
138 : double dfEndX, double dfEndY,
139 : int nVertexCount )
140 :
141 : {
142 : double dfStartAngle, dfEndAngle, dfAlongAngle;
143 : double dfCenterX, dfCenterY, dfRadius;
144 :
145 0 : if( !NTFArcCenterFromEdgePoints( dfStartX, dfStartY, dfAlongX, dfAlongY,
146 : dfEndX, dfEndY, &dfCenterX, &dfCenterY ) )
147 0 : return NULL;
148 :
149 0 : if( dfStartX == dfEndX && dfStartY == dfEndY )
150 : {
151 0 : dfStartAngle = 0.0;
152 0 : dfEndAngle = 360.0;
153 : }
154 : else
155 : {
156 : double dfDeltaX, dfDeltaY;
157 :
158 0 : dfDeltaX = dfStartX - dfCenterX;
159 0 : dfDeltaY = dfStartY - dfCenterY;
160 0 : dfStartAngle = atan2(dfDeltaY,dfDeltaX) * 180.0 / PI;
161 :
162 0 : dfDeltaX = dfAlongX - dfCenterX;
163 0 : dfDeltaY = dfAlongY - dfCenterY;
164 0 : dfAlongAngle = atan2(dfDeltaY,dfDeltaX) * 180.0 / PI;
165 :
166 0 : dfDeltaX = dfEndX - dfCenterX;
167 0 : dfDeltaY = dfEndY - dfCenterY;
168 0 : dfEndAngle = atan2(dfDeltaY,dfDeltaX) * 180.0 / PI;
169 :
170 : #ifdef notdef
171 : if( dfStartAngle > dfAlongAngle && dfAlongAngle > dfEndAngle )
172 : {
173 : double dfTempAngle;
174 :
175 : dfTempAngle = dfStartAngle;
176 : dfStartAngle = dfEndAngle;
177 : dfEndAngle = dfTempAngle;
178 : }
179 : #endif
180 :
181 0 : while( dfAlongAngle < dfStartAngle )
182 0 : dfAlongAngle += 360.0;
183 :
184 0 : while( dfEndAngle < dfAlongAngle )
185 0 : dfEndAngle += 360.0;
186 :
187 0 : if( dfEndAngle - dfStartAngle > 360.0 )
188 : {
189 : double dfTempAngle;
190 :
191 0 : dfTempAngle = dfStartAngle;
192 0 : dfStartAngle = dfEndAngle;
193 0 : dfEndAngle = dfTempAngle;
194 :
195 0 : while( dfEndAngle < dfStartAngle )
196 0 : dfStartAngle -= 360.0;
197 : }
198 : }
199 :
200 : dfRadius = sqrt( (dfCenterX - dfStartX) * (dfCenterX - dfStartX)
201 0 : + (dfCenterY - dfStartY) * (dfCenterY - dfStartY) );
202 :
203 : return NTFStrokeArcToOGRGeometry_Angles( dfCenterX, dfCenterY,
204 : dfRadius,
205 : dfStartAngle, dfEndAngle,
206 0 : nVertexCount );
207 : }
208 :
209 : /************************************************************************/
210 : /* NTFStrokeArcToOGRGeometry_Angles() */
211 : /************************************************************************/
212 :
213 : OGRGeometry *
214 0 : NTFStrokeArcToOGRGeometry_Angles( double dfCenterX, double dfCenterY,
215 : double dfRadius,
216 : double dfStartAngle, double dfEndAngle,
217 : int nVertexCount )
218 :
219 : {
220 0 : OGRLineString *poLine = new OGRLineString;
221 : double dfArcX, dfArcY, dfSlice;
222 : int iPoint;
223 :
224 0 : nVertexCount = MAX(2,nVertexCount);
225 0 : dfSlice = (dfEndAngle-dfStartAngle)/(nVertexCount-1);
226 :
227 0 : poLine->setNumPoints( nVertexCount );
228 :
229 0 : for( iPoint=0; iPoint < nVertexCount; iPoint++ )
230 : {
231 : double dfAngle;
232 :
233 0 : dfAngle = (dfStartAngle + iPoint * dfSlice) * PI / 180.0;
234 :
235 0 : dfArcX = dfCenterX + cos(dfAngle) * dfRadius;
236 0 : dfArcY = dfCenterY + sin(dfAngle) * dfRadius;
237 :
238 0 : poLine->setPoint( iPoint, dfArcX, dfArcY );
239 : }
240 :
241 0 : return poLine;
242 : }
243 :
244 :
|