source: ntrip/trunk/BNC/qwt/qwt_curve_fitter.cpp@ 5209

Last change on this file since 5209 was 4271, checked in by mervart, 12 years ago
File size: 9.0 KB
Line 
1/* -*- mode: C++ ; c-file-style: "stroustrup" -*- *****************************
2 * Qwt Widget Library
3 * Copyright (C) 1997 Josef Wilgen
4 * Copyright (C) 2002 Uwe Rathmann
5 *
6 * This library is free software; you can redistribute it and/or
7 * modify it under the terms of the Qwt License, Version 1.0
8 *****************************************************************************/
9
10#include "qwt_curve_fitter.h"
11#include "qwt_math.h"
12#include "qwt_spline.h"
13#include <qstack.h>
14#include <qvector.h>
15
16#if QT_VERSION < 0x040601
17#define qFabs(x) ::fabs(x)
18#endif
19
20//! Constructor
21QwtCurveFitter::QwtCurveFitter()
22{
23}
24
25//! Destructor
26QwtCurveFitter::~QwtCurveFitter()
27{
28}
29
30class QwtSplineCurveFitter::PrivateData
31{
32public:
33 PrivateData():
34 fitMode( QwtSplineCurveFitter::Auto ),
35 splineSize( 250 )
36 {
37 }
38
39 QwtSpline spline;
40 QwtSplineCurveFitter::FitMode fitMode;
41 int splineSize;
42};
43
44//! Constructor
45QwtSplineCurveFitter::QwtSplineCurveFitter()
46{
47 d_data = new PrivateData;
48}
49
50//! Destructor
51QwtSplineCurveFitter::~QwtSplineCurveFitter()
52{
53 delete d_data;
54}
55
56/*!
57 Select the algorithm used for building the spline
58
59 \param mode Mode representing a spline algorithm
60 \sa fitMode()
61*/
62void QwtSplineCurveFitter::setFitMode( FitMode mode )
63{
64 d_data->fitMode = mode;
65}
66
67/*!
68 \return Mode representing a spline algorithm
69 \sa setFitMode()
70*/
71QwtSplineCurveFitter::FitMode QwtSplineCurveFitter::fitMode() const
72{
73 return d_data->fitMode;
74}
75
76/*!
77 Assign a spline
78
79 \param spline Spline
80 \sa spline()
81*/
82void QwtSplineCurveFitter::setSpline( const QwtSpline &spline )
83{
84 d_data->spline = spline;
85 d_data->spline.reset();
86}
87
88/*!
89 \return Spline
90 \sa setSpline()
91*/
92const QwtSpline &QwtSplineCurveFitter::spline() const
93{
94 return d_data->spline;
95}
96
97/*!
98 \return Spline
99 \sa setSpline()
100*/
101QwtSpline &QwtSplineCurveFitter::spline()
102{
103 return d_data->spline;
104}
105
106/*!
107 Assign a spline size ( has to be at least 10 points )
108
109 \param splineSize Spline size
110 \sa splineSize()
111*/
112void QwtSplineCurveFitter::setSplineSize( int splineSize )
113{
114 d_data->splineSize = qMax( splineSize, 10 );
115}
116
117/*!
118 \return Spline size
119 \sa setSplineSize()
120*/
121int QwtSplineCurveFitter::splineSize() const
122{
123 return d_data->splineSize;
124}
125
126/*!
127 Find a curve which has the best fit to a series of data points
128
129 \param points Series of data points
130 \return Curve points
131*/
132QPolygonF QwtSplineCurveFitter::fitCurve( const QPolygonF &points ) const
133{
134 const int size = points.size();
135 if ( size <= 2 )
136 return points;
137
138 FitMode fitMode = d_data->fitMode;
139 if ( fitMode == Auto )
140 {
141 fitMode = Spline;
142
143 const QPointF *p = points.data();
144 for ( int i = 1; i < size; i++ )
145 {
146 if ( p[i].x() <= p[i-1].x() )
147 {
148 fitMode = ParametricSpline;
149 break;
150 }
151 };
152 }
153
154 if ( fitMode == ParametricSpline )
155 return fitParametric( points );
156 else
157 return fitSpline( points );
158}
159
160QPolygonF QwtSplineCurveFitter::fitSpline( const QPolygonF &points ) const
161{
162 d_data->spline.setPoints( points );
163 if ( !d_data->spline.isValid() )
164 return points;
165
166 QPolygonF fittedPoints( d_data->splineSize );
167
168 const double x1 = points[0].x();
169 const double x2 = points[int( points.size() - 1 )].x();
170 const double dx = x2 - x1;
171 const double delta = dx / ( d_data->splineSize - 1 );
172
173 for ( int i = 0; i < d_data->splineSize; i++ )
174 {
175 QPointF &p = fittedPoints[i];
176
177 const double v = x1 + i * delta;
178 const double sv = d_data->spline.value( v );
179
180 p.setX( v );
181 p.setY( sv );
182 }
183 d_data->spline.reset();
184
185 return fittedPoints;
186}
187
188QPolygonF QwtSplineCurveFitter::fitParametric( const QPolygonF &points ) const
189{
190 int i;
191 const int size = points.size();
192
193 QPolygonF fittedPoints( d_data->splineSize );
194 QPolygonF splinePointsX( size );
195 QPolygonF splinePointsY( size );
196
197 const QPointF *p = points.data();
198 QPointF *spX = splinePointsX.data();
199 QPointF *spY = splinePointsY.data();
200
201 double param = 0.0;
202 for ( i = 0; i < size; i++ )
203 {
204 const double x = p[i].x();
205 const double y = p[i].y();
206 if ( i > 0 )
207 {
208 const double delta = qSqrt( qwtSqr( x - spX[i-1].y() )
209 + qwtSqr( y - spY[i-1].y() ) );
210 param += qMax( delta, 1.0 );
211 }
212 spX[i].setX( param );
213 spX[i].setY( x );
214 spY[i].setX( param );
215 spY[i].setY( y );
216 }
217
218 d_data->spline.setPoints( splinePointsX );
219 if ( !d_data->spline.isValid() )
220 return points;
221
222 const double deltaX =
223 splinePointsX[size - 1].x() / ( d_data->splineSize - 1 );
224 for ( i = 0; i < d_data->splineSize; i++ )
225 {
226 const double dtmp = i * deltaX;
227 fittedPoints[i].setX( d_data->spline.value( dtmp ) );
228 }
229
230 d_data->spline.setPoints( splinePointsY );
231 if ( !d_data->spline.isValid() )
232 return points;
233
234 const double deltaY =
235 splinePointsY[size - 1].x() / ( d_data->splineSize - 1 );
236 for ( i = 0; i < d_data->splineSize; i++ )
237 {
238 const double dtmp = i * deltaY;
239 fittedPoints[i].setY( d_data->spline.value( dtmp ) );
240 }
241
242 return fittedPoints;
243}
244
245class QwtWeedingCurveFitter::PrivateData
246{
247public:
248 PrivateData():
249 tolerance( 1.0 )
250 {
251 }
252
253 double tolerance;
254};
255
256class QwtWeedingCurveFitter::Line
257{
258public:
259 Line( int i1 = 0, int i2 = 0 ):
260 from( i1 ),
261 to( i2 )
262 {
263 }
264
265 int from;
266 int to;
267};
268
269/*!
270 Constructor
271
272 \param tolerance Tolerance
273 \sa setTolerance(), tolerance()
274*/
275QwtWeedingCurveFitter::QwtWeedingCurveFitter( double tolerance )
276{
277 d_data = new PrivateData;
278 setTolerance( tolerance );
279}
280
281//! Destructor
282QwtWeedingCurveFitter::~QwtWeedingCurveFitter()
283{
284 delete d_data;
285}
286
287/*!
288 Assign the tolerance
289
290 The tolerance is the maximum distance, that is accaptable
291 between the original curve and the smoothed curve.
292
293 Increasing the tolerance will reduce the number of the
294 resulting points.
295
296 \param tolerance Tolerance
297
298 \sa tolerance()
299*/
300void QwtWeedingCurveFitter::setTolerance( double tolerance )
301{
302 d_data->tolerance = qMax( tolerance, 0.0 );
303}
304
305/*!
306 \return Tolerance
307 \sa setTolerance()
308*/
309double QwtWeedingCurveFitter::tolerance() const
310{
311 return d_data->tolerance;
312}
313
314/*!
315 \param points Series of data points
316 \return Curve points
317*/
318QPolygonF QwtWeedingCurveFitter::fitCurve( const QPolygonF &points ) const
319{
320 QStack<Line> stack;
321 stack.reserve( 500 );
322
323 const QPointF *p = points.data();
324 const int nPoints = points.size();
325
326 QVector<bool> usePoint( nPoints, false );
327
328 double distToSegment;
329
330 stack.push( Line( 0, nPoints - 1 ) );
331
332 while ( !stack.isEmpty() )
333 {
334 const Line r = stack.pop();
335
336 // initialize line segment
337 const double vecX = p[r.to].x() - p[r.from].x();
338 const double vecY = p[r.to].y() - p[r.from].y();
339
340 const double vecLength = qSqrt( vecX * vecX + vecY * vecY );
341
342 const double unitVecX = ( vecLength != 0.0 ) ? vecX / vecLength : 0.0;
343 const double unitVecY = ( vecLength != 0.0 ) ? vecY / vecLength : 0.0;
344
345 double maxDist = 0.0;
346 int nVertexIndexMaxDistance = r.from + 1;
347 for ( int i = r.from + 1; i < r.to; i++ )
348 {
349 //compare to anchor
350 const double fromVecX = p[i].x() - p[r.from].x();
351 const double fromVecY = p[i].y() - p[r.from].y();
352 const double fromVecLength =
353 qSqrt( fromVecX * fromVecX + fromVecY * fromVecY );
354
355 if ( fromVecX * unitVecX + fromVecY * unitVecY < 0.0 )
356 {
357 distToSegment = fromVecLength;
358 }
359 if ( fromVecX * unitVecX + fromVecY * unitVecY < 0.0 )
360 {
361 distToSegment = fromVecLength;
362 }
363 else
364 {
365 const double toVecX = p[i].x() - p[r.to].x();
366 const double toVecY = p[i].y() - p[r.to].y();
367 const double toVecLength = qSqrt( toVecX * toVecX + toVecY * toVecY );
368 const double s = toVecX * ( -unitVecX ) + toVecY * ( -unitVecY );
369 if ( s < 0.0 )
370 distToSegment = toVecLength;
371 else
372 {
373 distToSegment = qSqrt( qFabs( toVecLength * toVecLength - s * s ) );
374 }
375 }
376
377 if ( maxDist < distToSegment )
378 {
379 maxDist = distToSegment;
380 nVertexIndexMaxDistance = i;
381 }
382 }
383 if ( maxDist <= d_data->tolerance )
384 {
385 usePoint[r.from] = true;
386 usePoint[r.to] = true;
387 }
388 else
389 {
390 stack.push( Line( r.from, nVertexIndexMaxDistance ) );
391 stack.push( Line( nVertexIndexMaxDistance, r.to ) );
392 }
393 }
394
395 int cnt = 0;
396
397 QPolygonF stripped( nPoints );
398 for ( int i = 0; i < nPoints; i++ )
399 {
400 if ( usePoint[i] )
401 stripped[cnt++] = p[i];
402 }
403 stripped.resize( cnt );
404 return stripped;
405}
Note: See TracBrowser for help on using the repository browser.