← Back to branch summary

~ubuntu-branches/ubuntu/oneiric/qwt/oneiric-proposed

« back to all changes in this revision

Viewing changes to src/qwt_spline.cpp

  • Committer: Bazaar Package Importer
  • Author(s): Fathi Boudra
  • Date: 2011-06-10 11:22:47 UTC
  • mfrom: (1.1.6 upstream) (2.1.4 sid)
  • Revision ID: james.westby@ubuntu.com-20110610112247-0i1019vvmzaq6p86
Tags: 6.0.0-1
http://bugs.debian.org/624107
http://bugs.debian.org/604379
http://bugs.debian.org/626868
http://bugs.debian.org/626567
* New upstream release (Closes: #624107):
  - drop Qt3 support. (Closes: #604379, #626868)
* Register documentation with doc-base. (Closes: #626567)
* Drop patches:
  - 01_makefiles.diff
  - 02_add_missing_warnings.diff
  - 03_qwt_branch_pull_r544.diff
* Add qwt_install_paths.patch to fix the hardcoded installation paths.
* Update debian/control:
  - drop libqt3-mt-dev build dependency.
  - bump Standards-Version to 3.9.2 (no changes).
  - drop Qt3 related packages.
  - due to bump soname (and as we dropper Qt3 support):
    - libqwt5-qt4-dev -> libqwt-dev
    - libqwt5-qt4 -> libqwt6
    - libqwt5-doc -> libqwt-doc
* Update debian/copyright file.
* Update debian/rules: drop Qt3 packages support.

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/qwt_spline.cpp
 
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_spline.h"
 
11
#include "qwt_math.h"
 
12
 
 
13
class QwtSpline::PrivateData
 
14
{
 
15
public:
 
16
    PrivateData():
 
17
        splineType( QwtSpline::Natural )
 
18
    {
 
19
    }
 
20
 
 
21
    QwtSpline::SplineType splineType;
 
22
 
 
23
    // coefficient vectors
 
24
    QVector<double> a;
 
25
    QVector<double> b;
 
26
    QVector<double> c;
 
27
 
 
28
    // control points
 
29
    QPolygonF points;
 
30
};
 
31
 
 
32
static int lookup( double x, const QPolygonF &values )
 
33
{
 
34
#if 0
 
35
//qLowerBound/qHigherBound ???
 
36
#endif
 
37
    int i1;
 
38
    const int size = ( int )values.size();
 
39
 
 
40
    if ( x <= values[0].x() )
 
41
        i1 = 0;
 
42
    else if ( x >= values[size - 2].x() )
 
43
        i1 = size - 2;
 
44
    else
 
45
    {
 
46
        i1 = 0;
 
47
        int i2 = size - 2;
 
48
        int i3 = 0;
 
49
 
 
50
        while ( i2 - i1 > 1 )
 
51
        {
 
52
            i3 = i1 + ( ( i2 - i1 ) >> 1 );
 
53
 
 
54
            if ( values[i3].x() > x )
 
55
                i2 = i3;
 
56
            else
 
57
                i1 = i3;
 
58
        }
 
59
    }
 
60
    return i1;
 
61
}
 
62
 
 
63
//! Constructor
 
64
QwtSpline::QwtSpline()
 
65
{
 
66
    d_data = new PrivateData;
 
67
}
 
68
 
 
69
/*!
 
70
   Copy constructor
 
71
   \param other Spline used for initilization
 
72
*/
 
73
QwtSpline::QwtSpline( const QwtSpline& other )
 
74
{
 
75
    d_data = new PrivateData( *other.d_data );
 
76
}
 
77
 
 
78
/*!
 
79
   Assignment operator
 
80
   \param other Spline used for initilization
 
81
*/
 
82
QwtSpline &QwtSpline::operator=( const QwtSpline & other )
 
83
{
 
84
    *d_data = *other.d_data;
 
85
    return *this;
 
86
}
 
87
 
 
88
//! Destructor
 
89
QwtSpline::~QwtSpline()
 
90
{
 
91
    delete d_data;
 
92
}
 
93
 
 
94
/*!
 
95
   Select the algorithm used for calculating the spline
 
96
 
 
97
   \param splineType Spline type
 
98
   \sa splineType()
 
99
*/
 
100
void QwtSpline::setSplineType( SplineType splineType )
 
101
{
 
102
    d_data->splineType = splineType;
 
103
}
 
104
 
 
105
/*!
 
106
   \return the spline type
 
107
   \sa setSplineType()
 
108
*/
 
109
QwtSpline::SplineType QwtSpline::splineType() const
 
110
{
 
111
    return d_data->splineType;
 
112
}
 
113
 
 
114
/*!
 
115
  \brief Calculate the spline coefficients
 
116
 
 
117
  Depending on the value of \a periodic, this function
 
118
  will determine the coefficients for a natural or a periodic
 
119
  spline and store them internally.
 
120
 
 
121
  \param points Points
 
122
  \return true if successful
 
123
  \warning The sequence of x (but not y) values has to be strictly monotone
 
124
           increasing, which means <code>points[i].x() < points[i+1].x()</code>.
 
125
       If this is not the case, the function will return false
 
126
*/
 
127
bool QwtSpline::setPoints( const QPolygonF& points )
 
128
{
 
129
    const int size = points.size();
 
130
    if ( size <= 2 )
 
131
    {
 
132
        reset();
 
133
        return false;
 
134
    }
 
135
 
 
136
    d_data->points = points;
 
137
 
 
138
    d_data->a.resize( size - 1 );
 
139
    d_data->b.resize( size - 1 );
 
140
    d_data->c.resize( size - 1 );
 
141
 
 
142
    bool ok;
 
143
    if ( d_data->splineType == Periodic )
 
144
        ok = buildPeriodicSpline( points );
 
145
    else
 
146
        ok = buildNaturalSpline( points );
 
147
 
 
148
    if ( !ok )
 
149
        reset();
 
150
 
 
151
    return ok;
 
152
}
 
153
 
 
154
/*!
 
155
   Return points passed by setPoints
 
156
*/
 
157
QPolygonF QwtSpline::points() const
 
158
{
 
159
    return d_data->points;
 
160
}
 
161
 
 
162
//! \return A coefficients
 
163
const QVector<double> &QwtSpline::coefficientsA() const
 
164
{
 
165
    return d_data->a;
 
166
}
 
167
 
 
168
//! \return B coefficients
 
169
const QVector<double> &QwtSpline::coefficientsB() const
 
170
{
 
171
    return d_data->b;
 
172
}
 
173
 
 
174
//! \return C coefficients
 
175
const QVector<double> &QwtSpline::coefficientsC() const
 
176
{
 
177
    return d_data->c;
 
178
}
 
179
 
 
180
 
 
181
//! Free allocated memory and set size to 0
 
182
void QwtSpline::reset()
 
183
{
 
184
    d_data->a.resize( 0 );
 
185
    d_data->b.resize( 0 );
 
186
    d_data->c.resize( 0 );
 
187
    d_data->points.resize( 0 );
 
188
}
 
189
 
 
190
//! True if valid
 
191
bool QwtSpline::isValid() const
 
192
{
 
193
    return d_data->a.size() > 0;
 
194
}
 
195
 
 
196
/*!
 
197
  Calculate the interpolated function value corresponding
 
198
  to a given argument x.
 
199
*/
 
200
double QwtSpline::value( double x ) const
 
201
{
 
202
    if ( d_data->a.size() == 0 )
 
203
        return 0.0;
 
204
 
 
205
    const int i = lookup( x, d_data->points );
 
206
 
 
207
    const double delta = x - d_data->points[i].x();
 
208
    return( ( ( ( d_data->a[i] * delta ) + d_data->b[i] )
 
209
        * delta + d_data->c[i] ) * delta + d_data->points[i].y() );
 
210
}
 
211
 
 
212
/*!
 
213
  \brief Determines the coefficients for a natural spline
 
214
  \return true if successful
 
215
*/
 
216
bool QwtSpline::buildNaturalSpline( const QPolygonF &points )
 
217
{
 
218
    int i;
 
219
 
 
220
    const QPointF *p = points.data();
 
221
    const int size = points.size();
 
222
 
 
223
    double *a = d_data->a.data();
 
224
    double *b = d_data->b.data();
 
225
    double *c = d_data->c.data();
 
226
 
 
227
    //  set up tridiagonal equation system; use coefficient
 
228
    //  vectors as temporary buffers
 
229
    QVector<double> h( size - 1 );
 
230
    for ( i = 0; i < size - 1; i++ )
 
231
    {
 
232
        h[i] = p[i+1].x() - p[i].x();
 
233
        if ( h[i] <= 0 )
 
234
            return false;
 
235
    }
 
236
 
 
237
    QVector<double> d( size - 1 );
 
238
    double dy1 = ( p[1].y() - p[0].y() ) / h[0];
 
239
    for ( i = 1; i < size - 1; i++ )
 
240
    {
 
241
        b[i] = c[i] = h[i];
 
242
        a[i] = 2.0 * ( h[i-1] + h[i] );
 
243
 
 
244
        const double dy2 = ( p[i+1].y() - p[i].y() ) / h[i];
 
245
        d[i] = 6.0 * ( dy1 - dy2 );
 
246
        dy1 = dy2;
 
247
    }
 
248
 
 
249
    //
 
250
    // solve it
 
251
    //
 
252
 
 
253
    // L-U Factorization
 
254
    for ( i = 1; i < size - 2; i++ )
 
255
    {
 
256
        c[i] /= a[i];
 
257
        a[i+1] -= b[i] * c[i];
 
258
    }
 
259
 
 
260
    // forward elimination
 
261
    QVector<double> s( size );
 
262
    s[1] = d[1];
 
263
    for ( i = 2; i < size - 1; i++ )
 
264
        s[i] = d[i] - c[i-1] * s[i-1];
 
265
 
 
266
    // backward elimination
 
267
    s[size - 2] = - s[size - 2] / a[size - 2];
 
268
    for ( i = size - 3; i > 0; i-- )
 
269
        s[i] = - ( s[i] + b[i] * s[i+1] ) / a[i];
 
270
    s[size - 1] = s[0] = 0.0;
 
271
 
 
272
    //
 
273
    // Finally, determine the spline coefficients
 
274
    //
 
275
    for ( i = 0; i < size - 1; i++ )
 
276
    {
 
277
        a[i] = ( s[i+1] - s[i] ) / ( 6.0 * h[i] );
 
278
        b[i] = 0.5 * s[i];
 
279
        c[i] = ( p[i+1].y() - p[i].y() ) / h[i]
 
280
            - ( s[i+1] + 2.0 * s[i] ) * h[i] / 6.0;
 
281
    }
 
282
 
 
283
    return true;
 
284
}
 
285
 
 
286
/*!
 
287
  \brief Determines the coefficients for a periodic spline
 
288
  \return true if successful
 
289
*/
 
290
bool QwtSpline::buildPeriodicSpline( const QPolygonF &points )
 
291
{
 
292
    int i;
 
293
 
 
294
    const QPointF *p = points.data();
 
295
    const int size = points.size();
 
296
 
 
297
    double *a = d_data->a.data();
 
298
    double *b = d_data->b.data();
 
299
    double *c = d_data->c.data();
 
300
 
 
301
    QVector<double> d( size - 1 );
 
302
    QVector<double> h( size - 1 );
 
303
    QVector<double> s( size );
 
304
 
 
305
    //
 
306
    //  setup equation system; use coefficient
 
307
    //  vectors as temporary buffers
 
308
    //
 
309
    for ( i = 0; i < size - 1; i++ )
 
310
    {
 
311
        h[i] = p[i+1].x() - p[i].x();
 
312
        if ( h[i] <= 0.0 )
 
313
            return false;
 
314
    }
 
315
 
 
316
    const int imax = size - 2;
 
317
    double htmp = h[imax];
 
318
    double dy1 = ( p[0].y() - p[imax].y() ) / htmp;
 
319
    for ( i = 0; i <= imax; i++ )
 
320
    {
 
321
        b[i] = c[i] = h[i];
 
322
        a[i] = 2.0 * ( htmp + h[i] );
 
323
        const double dy2 = ( p[i+1].y() - p[i].y() ) / h[i];
 
324
        d[i] = 6.0 * ( dy1 - dy2 );
 
325
        dy1 = dy2;
 
326
        htmp = h[i];
 
327
    }
 
328
 
 
329
    //
 
330
    // solve it
 
331
    //
 
332
 
 
333
    // L-U Factorization
 
334
    a[0] = qSqrt( a[0] );
 
335
    c[0] = h[imax] / a[0];
 
336
    double sum = 0;
 
337
 
 
338
    for ( i = 0; i < imax - 1; i++ )
 
339
    {
 
340
        b[i] /= a[i];
 
341
        if ( i > 0 )
 
342
            c[i] = - c[i-1] * b[i-1] / a[i];
 
343
        a[i+1] = qSqrt( a[i+1] - qwtSqr( b[i] ) );
 
344
        sum += qwtSqr( c[i] );
 
345
    }
 
346
    b[imax-1] = ( b[imax-1] - c[imax-2] * b[imax-2] ) / a[imax-1];
 
347
    a[imax] = qSqrt( a[imax] - qwtSqr( b[imax-1] ) - sum );
 
348
 
 
349
 
 
350
    // forward elimination
 
351
    s[0] = d[0] / a[0];
 
352
    sum = 0;
 
353
    for ( i = 1; i < imax; i++ )
 
354
    {
 
355
        s[i] = ( d[i] - b[i-1] * s[i-1] ) / a[i];
 
356
        sum += c[i-1] * s[i-1];
 
357
    }
 
358
    s[imax] = ( d[imax] - b[imax-1] * s[imax-1] - sum ) / a[imax];
 
359
 
 
360
 
 
361
    // backward elimination
 
362
    s[imax] = - s[imax] / a[imax];
 
363
    s[imax-1] = -( s[imax-1] + b[imax-1] * s[imax] ) / a[imax-1];
 
364
    for ( i = imax - 2; i >= 0; i-- )
 
365
        s[i] = - ( s[i] + b[i] * s[i+1] + c[i] * s[imax] ) / a[i];
 
366
 
 
367
    //
 
368
    // Finally, determine the spline coefficients
 
369
    //
 
370
    s[size-1] = s[0];
 
371
    for ( i = 0; i < size - 1; i++ )
 
372
    {
 
373
        a[i] = ( s[i+1] - s[i] ) / ( 6.0 * h[i] );
 
374
        b[i] = 0.5 * s[i];
 
375
        c[i] = ( p[i+1].y() - p[i].y() )
 
376
            / h[i] - ( s[i+1] + 2.0 * s[i] ) * h[i] / 6.0;
 
377
    }
 
378
 
 
379
    return true;
 
380
}