Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1__all__ = ['splrep', 'splprep', 'splev', 'splint', 'sproot', 'spalde', 

2 'bisplrep', 'bisplev', 'insert', 'splder', 'splantider'] 

3 

4import warnings 

5 

6import numpy as np 

7 

8# These are in the API for fitpack even if not used in fitpack.py itself. 

9from ._fitpack_impl import bisplrep, bisplev, dblint 

10from . import _fitpack_impl as _impl 

11from ._bsplines import BSpline 

12 

13 

14def splprep(x, w=None, u=None, ub=None, ue=None, k=3, task=0, s=None, t=None, 

15 full_output=0, nest=None, per=0, quiet=1): 

16 """ 

17 Find the B-spline representation of an N-D curve. 

18 

19 Given a list of N rank-1 arrays, `x`, which represent a curve in 

20 N-D space parametrized by `u`, find a smooth approximating 

21 spline curve g(`u`). Uses the FORTRAN routine parcur from FITPACK. 

22 

23 Parameters 

24 ---------- 

25 x : array_like 

26 A list of sample vector arrays representing the curve. 

27 w : array_like, optional 

28 Strictly positive rank-1 array of weights the same length as `x[0]`. 

29 The weights are used in computing the weighted least-squares spline 

30 fit. If the errors in the `x` values have standard-deviation given by 

31 the vector d, then `w` should be 1/d. Default is ``ones(len(x[0]))``. 

32 u : array_like, optional 

33 An array of parameter values. If not given, these values are 

34 calculated automatically as ``M = len(x[0])``, where 

35 

36 v[0] = 0 

37 

38 v[i] = v[i-1] + distance(`x[i]`, `x[i-1]`) 

39 

40 u[i] = v[i] / v[M-1] 

41 

42 ub, ue : int, optional 

43 The end-points of the parameters interval. Defaults to 

44 u[0] and u[-1]. 

45 k : int, optional 

46 Degree of the spline. Cubic splines are recommended. 

47 Even values of `k` should be avoided especially with a small s-value. 

48 ``1 <= k <= 5``, default is 3. 

49 task : int, optional 

50 If task==0 (default), find t and c for a given smoothing factor, s. 

51 If task==1, find t and c for another value of the smoothing factor, s. 

52 There must have been a previous call with task=0 or task=1 

53 for the same set of data. 

54 If task=-1 find the weighted least square spline for a given set of 

55 knots, t. 

56 s : float, optional 

57 A smoothing condition. The amount of smoothness is determined by 

58 satisfying the conditions: ``sum((w * (y - g))**2,axis=0) <= s``, 

59 where g(x) is the smoothed interpolation of (x,y). The user can 

60 use `s` to control the trade-off between closeness and smoothness 

61 of fit. Larger `s` means more smoothing while smaller values of `s` 

62 indicate less smoothing. Recommended values of `s` depend on the 

63 weights, w. If the weights represent the inverse of the 

64 standard-deviation of y, then a good `s` value should be found in 

65 the range ``(m-sqrt(2*m),m+sqrt(2*m))``, where m is the number of 

66 data points in x, y, and w. 

67 t : int, optional 

68 The knots needed for task=-1. 

69 full_output : int, optional 

70 If non-zero, then return optional outputs. 

71 nest : int, optional 

72 An over-estimate of the total number of knots of the spline to 

73 help in determining the storage space. By default nest=m/2. 

74 Always large enough is nest=m+k+1. 

75 per : int, optional 

76 If non-zero, data points are considered periodic with period 

77 ``x[m-1] - x[0]`` and a smooth periodic spline approximation is 

78 returned. Values of ``y[m-1]`` and ``w[m-1]`` are not used. 

79 quiet : int, optional 

80 Non-zero to suppress messages. 

81 This parameter is deprecated; use standard Python warning filters 

82 instead. 

83 

84 Returns 

85 ------- 

86 tck : tuple 

87 (t,c,k) a tuple containing the vector of knots, the B-spline 

88 coefficients, and the degree of the spline. 

89 u : array 

90 An array of the values of the parameter. 

91 fp : float 

92 The weighted sum of squared residuals of the spline approximation. 

93 ier : int 

94 An integer flag about splrep success. Success is indicated 

95 if ier<=0. If ier in [1,2,3] an error occurred but was not raised. 

96 Otherwise an error is raised. 

97 msg : str 

98 A message corresponding to the integer flag, ier. 

99 

100 See Also 

101 -------- 

102 splrep, splev, sproot, spalde, splint, 

103 bisplrep, bisplev 

104 UnivariateSpline, BivariateSpline 

105 BSpline 

106 make_interp_spline 

107 

108 Notes 

109 ----- 

110 See `splev` for evaluation of the spline and its derivatives. 

111 The number of dimensions N must be smaller than 11. 

112 

113 The number of coefficients in the `c` array is ``k+1`` less then the number 

114 of knots, ``len(t)``. This is in contrast with `splrep`, which zero-pads 

115 the array of coefficients to have the same length as the array of knots. 

116 These additional coefficients are ignored by evaluation routines, `splev` 

117 and `BSpline`. 

118 

119 References 

120 ---------- 

121 .. [1] P. Dierckx, "Algorithms for smoothing data with periodic and 

122 parametric splines, Computer Graphics and Image Processing", 

123 20 (1982) 171-184. 

124 .. [2] P. Dierckx, "Algorithms for smoothing data with periodic and 

125 parametric splines", report tw55, Dept. Computer Science, 

126 K.U.Leuven, 1981. 

127 .. [3] P. Dierckx, "Curve and surface fitting with splines", Monographs on 

128 Numerical Analysis, Oxford University Press, 1993. 

129 

130 Examples 

131 -------- 

132 Generate a discretization of a limacon curve in the polar coordinates: 

133 

134 >>> phi = np.linspace(0, 2.*np.pi, 40) 

135 >>> r = 0.5 + np.cos(phi) # polar coords 

136 >>> x, y = r * np.cos(phi), r * np.sin(phi) # convert to cartesian 

137 

138 And interpolate: 

139 

140 >>> from scipy.interpolate import splprep, splev 

141 >>> tck, u = splprep([x, y], s=0) 

142 >>> new_points = splev(u, tck) 

143 

144 Notice that (i) we force interpolation by using `s=0`, 

145 (ii) the parameterization, ``u``, is generated automatically. 

146 Now plot the result: 

147 

148 >>> import matplotlib.pyplot as plt 

149 >>> fig, ax = plt.subplots() 

150 >>> ax.plot(x, y, 'ro') 

151 >>> ax.plot(new_points[0], new_points[1], 'r-') 

152 >>> plt.show() 

153 

154 """ 

155 res = _impl.splprep(x, w, u, ub, ue, k, task, s, t, full_output, nest, per, 

156 quiet) 

157 return res 

158 

159 

160def splrep(x, y, w=None, xb=None, xe=None, k=3, task=0, s=None, t=None, 

161 full_output=0, per=0, quiet=1): 

162 """ 

163 Find the B-spline representation of a 1-D curve. 

164 

165 Given the set of data points ``(x[i], y[i])`` determine a smooth spline 

166 approximation of degree k on the interval ``xb <= x <= xe``. 

167 

168 Parameters 

169 ---------- 

170 x, y : array_like 

171 The data points defining a curve y = f(x). 

172 w : array_like, optional 

173 Strictly positive rank-1 array of weights the same length as x and y. 

174 The weights are used in computing the weighted least-squares spline 

175 fit. If the errors in the y values have standard-deviation given by the 

176 vector d, then w should be 1/d. Default is ones(len(x)). 

177 xb, xe : float, optional 

178 The interval to fit. If None, these default to x[0] and x[-1] 

179 respectively. 

180 k : int, optional 

181 The degree of the spline fit. It is recommended to use cubic splines. 

182 Even values of k should be avoided especially with small s values. 

183 1 <= k <= 5 

184 task : {1, 0, -1}, optional 

185 If task==0 find t and c for a given smoothing factor, s. 

186 

187 If task==1 find t and c for another value of the smoothing factor, s. 

188 There must have been a previous call with task=0 or task=1 for the same 

189 set of data (t will be stored an used internally) 

190 

191 If task=-1 find the weighted least square spline for a given set of 

192 knots, t. These should be interior knots as knots on the ends will be 

193 added automatically. 

194 s : float, optional 

195 A smoothing condition. The amount of smoothness is determined by 

196 satisfying the conditions: sum((w * (y - g))**2,axis=0) <= s where g(x) 

197 is the smoothed interpolation of (x,y). The user can use s to control 

198 the tradeoff between closeness and smoothness of fit. Larger s means 

199 more smoothing while smaller values of s indicate less smoothing. 

200 Recommended values of s depend on the weights, w. If the weights 

201 represent the inverse of the standard-deviation of y, then a good s 

202 value should be found in the range (m-sqrt(2*m),m+sqrt(2*m)) where m is 

203 the number of datapoints in x, y, and w. default : s=m-sqrt(2*m) if 

204 weights are supplied. s = 0.0 (interpolating) if no weights are 

205 supplied. 

206 t : array_like, optional 

207 The knots needed for task=-1. If given then task is automatically set 

208 to -1. 

209 full_output : bool, optional 

210 If non-zero, then return optional outputs. 

211 per : bool, optional 

212 If non-zero, data points are considered periodic with period x[m-1] - 

213 x[0] and a smooth periodic spline approximation is returned. Values of 

214 y[m-1] and w[m-1] are not used. 

215 quiet : bool, optional 

216 Non-zero to suppress messages. 

217 This parameter is deprecated; use standard Python warning filters 

218 instead. 

219 

220 Returns 

221 ------- 

222 tck : tuple 

223 A tuple (t,c,k) containing the vector of knots, the B-spline 

224 coefficients, and the degree of the spline. 

225 fp : array, optional 

226 The weighted sum of squared residuals of the spline approximation. 

227 ier : int, optional 

228 An integer flag about splrep success. Success is indicated if ier<=0. 

229 If ier in [1,2,3] an error occurred but was not raised. Otherwise an 

230 error is raised. 

231 msg : str, optional 

232 A message corresponding to the integer flag, ier. 

233 

234 See Also 

235 -------- 

236 UnivariateSpline, BivariateSpline 

237 splprep, splev, sproot, spalde, splint 

238 bisplrep, bisplev 

239 BSpline 

240 make_interp_spline 

241 

242 Notes 

243 ----- 

244 See `splev` for evaluation of the spline and its derivatives. Uses the 

245 FORTRAN routine ``curfit`` from FITPACK. 

246 

247 The user is responsible for assuring that the values of `x` are unique. 

248 Otherwise, `splrep` will not return sensible results. 

249 

250 If provided, knots `t` must satisfy the Schoenberg-Whitney conditions, 

251 i.e., there must be a subset of data points ``x[j]`` such that 

252 ``t[j] < x[j] < t[j+k+1]``, for ``j=0, 1,...,n-k-2``. 

253 

254 This routine zero-pads the coefficients array ``c`` to have the same length 

255 as the array of knots ``t`` (the trailing ``k + 1`` coefficients are ignored 

256 by the evaluation routines, `splev` and `BSpline`.) This is in contrast with 

257 `splprep`, which does not zero-pad the coefficients. 

258 

259 References 

260 ---------- 

261 Based on algorithms described in [1]_, [2]_, [3]_, and [4]_: 

262 

263 .. [1] P. Dierckx, "An algorithm for smoothing, differentiation and 

264 integration of experimental data using spline functions", 

265 J.Comp.Appl.Maths 1 (1975) 165-184. 

266 .. [2] P. Dierckx, "A fast algorithm for smoothing data on a rectangular 

267 grid while using spline functions", SIAM J.Numer.Anal. 19 (1982) 

268 1286-1304. 

269 .. [3] P. Dierckx, "An improved algorithm for curve fitting with spline 

270 functions", report tw54, Dept. Computer Science,K.U. Leuven, 1981. 

271 .. [4] P. Dierckx, "Curve and surface fitting with splines", Monographs on 

272 Numerical Analysis, Oxford University Press, 1993. 

273 

274 Examples 

275 -------- 

276 

277 >>> import matplotlib.pyplot as plt 

278 >>> from scipy.interpolate import splev, splrep 

279 >>> x = np.linspace(0, 10, 10) 

280 >>> y = np.sin(x) 

281 >>> spl = splrep(x, y) 

282 >>> x2 = np.linspace(0, 10, 200) 

283 >>> y2 = splev(x2, spl) 

284 >>> plt.plot(x, y, 'o', x2, y2) 

285 >>> plt.show() 

286 

287 """ 

288 res = _impl.splrep(x, y, w, xb, xe, k, task, s, t, full_output, per, quiet) 

289 return res 

290 

291 

292def splev(x, tck, der=0, ext=0): 

293 """ 

294 Evaluate a B-spline or its derivatives. 

295 

296 Given the knots and coefficients of a B-spline representation, evaluate 

297 the value of the smoothing polynomial and its derivatives. This is a 

298 wrapper around the FORTRAN routines splev and splder of FITPACK. 

299 

300 Parameters 

301 ---------- 

302 x : array_like 

303 An array of points at which to return the value of the smoothed 

304 spline or its derivatives. If `tck` was returned from `splprep`, 

305 then the parameter values, u should be given. 

306 tck : 3-tuple or a BSpline object 

307 If a tuple, then it should be a sequence of length 3 returned by 

308 `splrep` or `splprep` containing the knots, coefficients, and degree 

309 of the spline. (Also see Notes.) 

310 der : int, optional 

311 The order of derivative of the spline to compute (must be less than 

312 or equal to k, the degree of the spline). 

313 ext : int, optional 

314 Controls the value returned for elements of ``x`` not in the 

315 interval defined by the knot sequence. 

316 

317 * if ext=0, return the extrapolated value. 

318 * if ext=1, return 0 

319 * if ext=2, raise a ValueError 

320 * if ext=3, return the boundary value. 

321 

322 The default value is 0. 

323 

324 Returns 

325 ------- 

326 y : ndarray or list of ndarrays 

327 An array of values representing the spline function evaluated at 

328 the points in `x`. If `tck` was returned from `splprep`, then this 

329 is a list of arrays representing the curve in an N-D space. 

330 

331 Notes 

332 ----- 

333 Manipulating the tck-tuples directly is not recommended. In new code, 

334 prefer using `BSpline` objects. 

335 

336 See Also 

337 -------- 

338 splprep, splrep, sproot, spalde, splint 

339 bisplrep, bisplev 

340 BSpline 

341 

342 References 

343 ---------- 

344 .. [1] C. de Boor, "On calculating with b-splines", J. Approximation 

345 Theory, 6, p.50-62, 1972. 

346 .. [2] M. G. Cox, "The numerical evaluation of b-splines", J. Inst. Maths 

347 Applics, 10, p.134-149, 1972. 

348 .. [3] P. Dierckx, "Curve and surface fitting with splines", Monographs 

349 on Numerical Analysis, Oxford University Press, 1993. 

350 

351 """ 

352 if isinstance(tck, BSpline): 

353 if tck.c.ndim > 1: 

354 mesg = ("Calling splev() with BSpline objects with c.ndim > 1 is " 

355 "not recommended. Use BSpline.__call__(x) instead.") 

356 warnings.warn(mesg, DeprecationWarning) 

357 

358 # remap the out-of-bounds behavior 

359 try: 

360 extrapolate = {0: True, }[ext] 

361 except KeyError: 

362 raise ValueError("Extrapolation mode %s is not supported " 

363 "by BSpline." % ext) 

364 

365 return tck(x, der, extrapolate=extrapolate) 

366 else: 

367 return _impl.splev(x, tck, der, ext) 

368 

369 

370def splint(a, b, tck, full_output=0): 

371 """ 

372 Evaluate the definite integral of a B-spline between two given points. 

373 

374 Parameters 

375 ---------- 

376 a, b : float 

377 The end-points of the integration interval. 

378 tck : tuple or a BSpline instance 

379 If a tuple, then it should be a sequence of length 3, containing the 

380 vector of knots, the B-spline coefficients, and the degree of the 

381 spline (see `splev`). 

382 full_output : int, optional 

383 Non-zero to return optional output. 

384 

385 Returns 

386 ------- 

387 integral : float 

388 The resulting integral. 

389 wrk : ndarray 

390 An array containing the integrals of the normalized B-splines 

391 defined on the set of knots. 

392 (Only returned if `full_output` is non-zero) 

393 

394 Notes 

395 ----- 

396 `splint` silently assumes that the spline function is zero outside the data 

397 interval (`a`, `b`). 

398 

399 Manipulating the tck-tuples directly is not recommended. In new code, 

400 prefer using the `BSpline` objects. 

401 

402 See Also 

403 -------- 

404 splprep, splrep, sproot, spalde, splev 

405 bisplrep, bisplev 

406 BSpline 

407 

408 References 

409 ---------- 

410 .. [1] P.W. Gaffney, The calculation of indefinite integrals of b-splines", 

411 J. Inst. Maths Applics, 17, p.37-41, 1976. 

412 .. [2] P. Dierckx, "Curve and surface fitting with splines", Monographs 

413 on Numerical Analysis, Oxford University Press, 1993. 

414 

415 """ 

416 if isinstance(tck, BSpline): 

417 if tck.c.ndim > 1: 

418 mesg = ("Calling splint() with BSpline objects with c.ndim > 1 is " 

419 "not recommended. Use BSpline.integrate() instead.") 

420 warnings.warn(mesg, DeprecationWarning) 

421 

422 if full_output != 0: 

423 mesg = ("full_output = %s is not supported. Proceeding as if " 

424 "full_output = 0" % full_output) 

425 

426 return tck.integrate(a, b, extrapolate=False) 

427 else: 

428 return _impl.splint(a, b, tck, full_output) 

429 

430 

431def sproot(tck, mest=10): 

432 """ 

433 Find the roots of a cubic B-spline. 

434 

435 Given the knots (>=8) and coefficients of a cubic B-spline return the 

436 roots of the spline. 

437 

438 Parameters 

439 ---------- 

440 tck : tuple or a BSpline object 

441 If a tuple, then it should be a sequence of length 3, containing the 

442 vector of knots, the B-spline coefficients, and the degree of the 

443 spline. 

444 The number of knots must be >= 8, and the degree must be 3. 

445 The knots must be a montonically increasing sequence. 

446 mest : int, optional 

447 An estimate of the number of zeros (Default is 10). 

448 

449 Returns 

450 ------- 

451 zeros : ndarray 

452 An array giving the roots of the spline. 

453 

454 Notes 

455 ----- 

456 Manipulating the tck-tuples directly is not recommended. In new code, 

457 prefer using the `BSpline` objects. 

458 

459 See also 

460 -------- 

461 splprep, splrep, splint, spalde, splev 

462 bisplrep, bisplev 

463 BSpline 

464 

465 

466 References 

467 ---------- 

468 .. [1] C. de Boor, "On calculating with b-splines", J. Approximation 

469 Theory, 6, p.50-62, 1972. 

470 .. [2] M. G. Cox, "The numerical evaluation of b-splines", J. Inst. Maths 

471 Applics, 10, p.134-149, 1972. 

472 .. [3] P. Dierckx, "Curve and surface fitting with splines", Monographs 

473 on Numerical Analysis, Oxford University Press, 1993. 

474 

475 """ 

476 if isinstance(tck, BSpline): 

477 if tck.c.ndim > 1: 

478 mesg = ("Calling sproot() with BSpline objects with c.ndim > 1 is " 

479 "not recommended.") 

480 warnings.warn(mesg, DeprecationWarning) 

481 

482 t, c, k = tck.tck 

483 

484 # _impl.sproot expects the interpolation axis to be last, so roll it. 

485 # NB: This transpose is a no-op if c is 1D. 

486 sh = tuple(range(c.ndim)) 

487 c = c.transpose(sh[1:] + (0,)) 

488 return _impl.sproot((t, c, k), mest) 

489 else: 

490 return _impl.sproot(tck, mest) 

491 

492 

493def spalde(x, tck): 

494 """ 

495 Evaluate all derivatives of a B-spline. 

496 

497 Given the knots and coefficients of a cubic B-spline compute all 

498 derivatives up to order k at a point (or set of points). 

499 

500 Parameters 

501 ---------- 

502 x : array_like 

503 A point or a set of points at which to evaluate the derivatives. 

504 Note that ``t(k) <= x <= t(n-k+1)`` must hold for each `x`. 

505 tck : tuple 

506 A tuple ``(t, c, k)``, containing the vector of knots, the B-spline 

507 coefficients, and the degree of the spline (see `splev`). 

508 

509 Returns 

510 ------- 

511 results : {ndarray, list of ndarrays} 

512 An array (or a list of arrays) containing all derivatives 

513 up to order k inclusive for each point `x`. 

514 

515 See Also 

516 -------- 

517 splprep, splrep, splint, sproot, splev, bisplrep, bisplev, 

518 BSpline 

519 

520 References 

521 ---------- 

522 .. [1] C. de Boor: On calculating with b-splines, J. Approximation Theory 

523 6 (1972) 50-62. 

524 .. [2] M. G. Cox : The numerical evaluation of b-splines, J. Inst. Maths 

525 applics 10 (1972) 134-149. 

526 .. [3] P. Dierckx : Curve and surface fitting with splines, Monographs on 

527 Numerical Analysis, Oxford University Press, 1993. 

528 

529 """ 

530 if isinstance(tck, BSpline): 

531 raise TypeError("spalde does not accept BSpline instances.") 

532 else: 

533 return _impl.spalde(x, tck) 

534 

535 

536def insert(x, tck, m=1, per=0): 

537 """ 

538 Insert knots into a B-spline. 

539 

540 Given the knots and coefficients of a B-spline representation, create a 

541 new B-spline with a knot inserted `m` times at point `x`. 

542 This is a wrapper around the FORTRAN routine insert of FITPACK. 

543 

544 Parameters 

545 ---------- 

546 x (u) : array_like 

547 A 1-D point at which to insert a new knot(s). If `tck` was returned 

548 from ``splprep``, then the parameter values, u should be given. 

549 tck : a `BSpline` instance or a tuple 

550 If tuple, then it is expected to be a tuple (t,c,k) containing 

551 the vector of knots, the B-spline coefficients, and the degree of 

552 the spline. 

553 m : int, optional 

554 The number of times to insert the given knot (its multiplicity). 

555 Default is 1. 

556 per : int, optional 

557 If non-zero, the input spline is considered periodic. 

558 

559 Returns 

560 ------- 

561 BSpline instance or a tuple 

562 A new B-spline with knots t, coefficients c, and degree k. 

563 ``t(k+1) <= x <= t(n-k)``, where k is the degree of the spline. 

564 In case of a periodic spline (``per != 0``) there must be 

565 either at least k interior knots t(j) satisfying ``t(k+1)<t(j)<=x`` 

566 or at least k interior knots t(j) satisfying ``x<=t(j)<t(n-k)``. 

567 A tuple is returned iff the input argument `tck` is a tuple, otherwise 

568 a BSpline object is constructed and returned. 

569 

570 Notes 

571 ----- 

572 Based on algorithms from [1]_ and [2]_. 

573 

574 Manipulating the tck-tuples directly is not recommended. In new code, 

575 prefer using the `BSpline` objects. 

576 

577 References 

578 ---------- 

579 .. [1] W. Boehm, "Inserting new knots into b-spline curves.", 

580 Computer Aided Design, 12, p.199-201, 1980. 

581 .. [2] P. Dierckx, "Curve and surface fitting with splines, Monographs on 

582 Numerical Analysis", Oxford University Press, 1993. 

583 

584 """ 

585 if isinstance(tck, BSpline): 

586 

587 t, c, k = tck.tck 

588 

589 # FITPACK expects the interpolation axis to be last, so roll it over 

590 # NB: if c array is 1D, transposes are no-ops 

591 sh = tuple(range(c.ndim)) 

592 c = c.transpose(sh[1:] + (0,)) 

593 t_, c_, k_ = _impl.insert(x, (t, c, k), m, per) 

594 

595 # and roll the last axis back 

596 c_ = np.asarray(c_) 

597 c_ = c_.transpose((sh[-1],) + sh[:-1]) 

598 return BSpline(t_, c_, k_) 

599 else: 

600 return _impl.insert(x, tck, m, per) 

601 

602 

603def splder(tck, n=1): 

604 """ 

605 Compute the spline representation of the derivative of a given spline 

606 

607 Parameters 

608 ---------- 

609 tck : BSpline instance or a tuple of (t, c, k) 

610 Spline whose derivative to compute 

611 n : int, optional 

612 Order of derivative to evaluate. Default: 1 

613 

614 Returns 

615 ------- 

616 `BSpline` instance or tuple 

617 Spline of order k2=k-n representing the derivative 

618 of the input spline. 

619 A tuple is returned iff the input argument `tck` is a tuple, otherwise 

620 a BSpline object is constructed and returned. 

621 

622 Notes 

623 ----- 

624 

625 .. versionadded:: 0.13.0 

626 

627 See Also 

628 -------- 

629 splantider, splev, spalde 

630 BSpline 

631 

632 Examples 

633 -------- 

634 This can be used for finding maxima of a curve: 

635 

636 >>> from scipy.interpolate import splrep, splder, sproot 

637 >>> x = np.linspace(0, 10, 70) 

638 >>> y = np.sin(x) 

639 >>> spl = splrep(x, y, k=4) 

640 

641 Now, differentiate the spline and find the zeros of the 

642 derivative. (NB: `sproot` only works for order 3 splines, so we 

643 fit an order 4 spline): 

644 

645 >>> dspl = splder(spl) 

646 >>> sproot(dspl) / np.pi 

647 array([ 0.50000001, 1.5 , 2.49999998]) 

648 

649 This agrees well with roots :math:`\\pi/2 + n\\pi` of 

650 :math:`\\cos(x) = \\sin'(x)`. 

651 

652 """ 

653 if isinstance(tck, BSpline): 

654 return tck.derivative(n) 

655 else: 

656 return _impl.splder(tck, n) 

657 

658 

659def splantider(tck, n=1): 

660 """ 

661 Compute the spline for the antiderivative (integral) of a given spline. 

662 

663 Parameters 

664 ---------- 

665 tck : BSpline instance or a tuple of (t, c, k) 

666 Spline whose antiderivative to compute 

667 n : int, optional 

668 Order of antiderivative to evaluate. Default: 1 

669 

670 Returns 

671 ------- 

672 BSpline instance or a tuple of (t2, c2, k2) 

673 Spline of order k2=k+n representing the antiderivative of the input 

674 spline. 

675 A tuple is returned iff the input argument `tck` is a tuple, otherwise 

676 a BSpline object is constructed and returned. 

677 

678 See Also 

679 -------- 

680 splder, splev, spalde 

681 BSpline 

682 

683 Notes 

684 ----- 

685 The `splder` function is the inverse operation of this function. 

686 Namely, ``splder(splantider(tck))`` is identical to `tck`, modulo 

687 rounding error. 

688 

689 .. versionadded:: 0.13.0 

690 

691 Examples 

692 -------- 

693 >>> from scipy.interpolate import splrep, splder, splantider, splev 

694 >>> x = np.linspace(0, np.pi/2, 70) 

695 >>> y = 1 / np.sqrt(1 - 0.8*np.sin(x)**2) 

696 >>> spl = splrep(x, y) 

697 

698 The derivative is the inverse operation of the antiderivative, 

699 although some floating point error accumulates: 

700 

701 >>> splev(1.7, spl), splev(1.7, splder(splantider(spl))) 

702 (array(2.1565429877197317), array(2.1565429877201865)) 

703 

704 Antiderivative can be used to evaluate definite integrals: 

705 

706 >>> ispl = splantider(spl) 

707 >>> splev(np.pi/2, ispl) - splev(0, ispl) 

708 2.2572053588768486 

709 

710 This is indeed an approximation to the complete elliptic integral 

711 :math:`K(m) = \\int_0^{\\pi/2} [1 - m\\sin^2 x]^{-1/2} dx`: 

712 

713 >>> from scipy.special import ellipk 

714 >>> ellipk(0.8) 

715 2.2572053268208538 

716 

717 """ 

718 if isinstance(tck, BSpline): 

719 return tck.antiderivative(n) 

720 else: 

721 return _impl.splantider(tck, n)