1 subroutine profil(iopt,tx,nx,ty,ny,c,kx,ky,u,nu,cu,ier) 2c if iopt=0 subroutine profil calculates the b-spline coefficients of 3c the univariate spline f(y) = s(u,y) with s(x,y) a bivariate spline of 4c degrees kx and ky, given in the b-spline representation. 5c if iopt = 1 it calculates the b-spline coefficients of the univariate 6c spline g(x) = s(x,u) 7c 8c calling sequence: 9c call profil(iopt,tx,nx,ty,ny,c,kx,ky,u,nu,cu,ier) 10c 11c input parameters: 12c iopt : integer flag, specifying whether the profile f(y) (iopt=0) 13c or the profile g(x) (iopt=1) must be determined. 14c tx : real array, length nx, which contains the position of the 15c knots in the x-direction. 16c nx : integer, giving the total number of knots in the x-direction 17c ty : real array, length ny, which contains the position of the 18c knots in the y-direction. 19c ny : integer, giving the total number of knots in the y-direction 20c c : real array, length (nx-kx-1)*(ny-ky-1), which contains the 21c b-spline coefficients. 22c kx,ky : integer values, giving the degrees of the spline. 23c u : real value, specifying the requested profile. 24c tx(kx+1)<=u<=tx(nx-kx), if iopt=0. 25c ty(ky+1)<=u<=ty(ny-ky), if iopt=1. 26c nu : on entry nu must specify the dimension of the array cu. 27c nu >= ny if iopt=0, nu >= nx if iopt=1. 28c 29c output parameters: 30c cu : real array of dimension (nu). 31c on successful exit this array contains the b-spline 32c ier : integer error flag 33c ier=0 : normal return 34c ier=10: invalid input data (see restrictions) 35c 36c restrictions: 37c if iopt=0 : tx(kx+1) <= u <= tx(nx-kx), nu >=ny. 38c if iopt=1 : ty(ky+1) <= u <= ty(ny-ky), nu >=nx. 39c 40c other subroutines required: 41c fpbspl 42c 43c author : 44c p.dierckx 45c dept. computer science, k.u.leuven 46c celestijnenlaan 200a, b-3001 heverlee, belgium. 47c e-mail : Paul.Dierckx@cs.kuleuven.ac.be 48c 49c latest update : march 1987 50c 51c ..scalar arguments.. 52 integer iopt,nx,ny,kx,ky,nu,ier 53 real*8 u 54c ..array arguments.. 55 real*8 tx(nx),ty(ny),c((nx-kx-1)*(ny-ky-1)),cu(nu) 56c ..local scalars.. 57 integer i,j,kx1,ky1,l,l1,m,m0,nkx1,nky1 58 real*8 sum 59c ..local array 60 real*8 h(6) 61c .. 62c before starting computations a data check is made. if the input data 63c are invalid control is immediately repassed to the calling program. 64 kx1 = kx+1 65 ky1 = ky+1 66 nkx1 = nx-kx1 67 nky1 = ny-ky1 68 ier = 10 69 if(iopt.ne.0) go to 200 70 if(nu.lt.ny) go to 300 71 if(u.lt.tx(kx1) .or. u.gt.tx(nkx1+1)) go to 300 72c the b-splinecoefficients of f(y) = s(u,y). 73 ier = 0 74 l = kx1 75 l1 = l+1 76 110 if(u.lt.tx(l1) .or. l.eq.nkx1) go to 120 77 l = l1 78 l1 = l+1 79 go to 110 80 120 call fpbspl(tx,nx,kx,u,l,h) 81 m0 = (l-kx1)*nky1+1 82 do 140 i=1,nky1 83 m = m0 84 sum = 0. 85 do 130 j=1,kx1 86 sum = sum+h(j)*c(m) 87 m = m+nky1 88 130 continue 89 cu(i) = sum 90 m0 = m0+1 91 140 continue 92 go to 300 93 200 if(nu.lt.nx) go to 300 94 if(u.lt.ty(ky1) .or. u.gt.ty(nky1+1)) go to 300 95c the b-splinecoefficients of g(x) = s(x,u). 96 ier = 0 97 l = ky1 98 l1 = l+1 99 210 if(u.lt.ty(l1) .or. l.eq.nky1) go to 220 100 l = l1 101 l1 = l+1 102 go to 210 103 220 call fpbspl(ty,ny,ky,u,l,h) 104 m0 = l-ky 105 do 240 i=1,nkx1 106 m = m0 107 sum = 0. 108 do 230 j=1,ky1 109 sum = sum+h(j)*c(m) 110 m = m+1 111 230 continue 112 cu(i) = sum 113 m0 = m0+nky1 114 240 continue 115 300 return 116 end 117 118