V3FIT
ezspline_load.f90
1 !/////
2 ! R8 !
3 !/////
4 !
5 ! 1-D
6 !
7 
8 ! DMC modifications -- spl_name optional argument.
9 ! a single file can contain more than one spline (or not). If it does
10 ! contain more than one spline, then each one must have a distinct name,
11 ! which must be specified in the call to tell which one to read this time.
12 
13 ! fullsv -- ezspline_save now has a "full_save" option which means, save
14 ! the spline including coefficients. ezspline_load is modified to look
15 ! for and read fullsv data if it is present; if not, the spline coefficients
16 ! are recalculated.
17 
18  subroutine ezspline_load1_r8(spline_o, filename, ier, spl_name)
19  use ezspline_obj
20  use ezcdf
21  implicit none
22  type(EZspline1_r8) :: spline_o
23  character*(*) :: filename
24  ! ier:
25  ! 21=could not load spline object from file filename
26  ! 22=loaded spline object from file filename but failed at coefficient set-up
27  integer, intent(out) :: ier
28 
29  character*(*), intent(in), OPTIONAL :: spl_name ! specify name of spline
30  ! to be read-- in case file contains more than one.
31 
32  integer ncid, ifail
33  integer dimlens(3)
34  character(2), parameter :: real8='R8'
35  character(3), parameter :: int='INT'
36  integer n1, BCS1(2)
37  real(ezspline_r8), dimension(:), allocatable :: f
38 
39  character*4 :: data_type
40 
41  character*32 zpre
42  logical :: fullsv
43 
44  if(present(spl_name)) then
45  zpre = spl_name
46  else
47  zpre = ' '
48  endif
49 
50  ier = 0
51  call cdfopn(ncid, filename, 'r') ! no error flag??
52  if(ncid==0) then
53  ier=43
54  return
55  endif
56 
57  ! check if n1 is present; check if spline object has correct rank
58  ! n2 should NOT be present.
59 
60  call cdfinqvar(ncid, trim(zpre)//'n1', dimlens, data_type, ifail)
61  if(ifail.ne.0) then
62  ier = 21
63  return
64  endif
65  call cdfinqvar(ncid, trim(zpre)//'n2', dimlens, data_type, ifail)
66  if(ifail.eq.0) then
67  ier = 20
68  return
69  endif
70 
71  ! check if fullsv was in effect -- if so, x1pkg will be found.
72 
73  fullsv=.false.
74  call cdfinqvar(ncid, trim(zpre)//'x1pkg', dimlens, data_type, ifail)
75  if(ifail.eq.0) then
76  fullsv=.true.
77  endif
78 
79  if( spline_o%nguard /= 123456789 ) call ezspline_preinit(spline_o)
80 
81  if( ezspline_allocated(spline_o) ) call ezspline_free1_r8(spline_o, ier)
82 
83  call cdfgetvar(ncid, trim(zpre)//'n1', n1, ifail)
84  call cdfgetvar(ncid, trim(zpre)//'isLinear', spline_o%isLinear, ifail)
85  if(ifail.ne.0) spline_o%isLinear=0
86 
87  if(spline_o%isLinear.eq.1) then
88  call ezlinear_init1_r8(spline_o, n1, ier)
89  else
90  bcs1 = (/0, 0/)
91  call ezspline_init1_r8(spline_o, n1, bcs1, ier)
92  endif
93  if(ier.ne.0) return
94 
95  call cdfgetvar(ncid, trim(zpre)//'klookup1', spline_o%klookup1, ifail)
96  if(ifail.ne.0) spline_o%klookup1=-3 ! the old default
97 
98  call cdfgetvar(ncid, trim(zpre)//'isHermite', spline_o%isHermite, ifail)
99  call cdfgetvar(ncid, trim(zpre)//'isReady', spline_o%isReady, ifail)
100  call cdfgetvar(ncid, trim(zpre)//'ibctype1', spline_o%ibctype1, ifail)
101  call cdfgetvar(ncid, trim(zpre)//'x1', spline_o%x1, ifail)
102  call cdfgetvar(ncid, trim(zpre)//'bcval1min', spline_o%bcval1min, ifail)
103  call cdfgetvar(ncid, trim(zpre)//'bcval1max', spline_o%bcval1max, ifail)
104  if(ifail/=0) then
105  ier=44
106  return
107  endif
108 
109  if(.not.fullsv) then
110  allocate(f(n1))
111  if(spline_o%isReady==1) then
112  call cdfgetvar(ncid, trim(zpre)//'f', f, ifail)
113  else
114  ier = 92
115  endif
116  else
117  call cdfgetvar(ncid, trim(zpre)//'x1min', spline_o%x1min, ifail)
118  call cdfgetvar(ncid, trim(zpre)//'x1max', spline_o%x1max, ifail)
119  call cdfgetvar(ncid, trim(zpre)//'ilin1', spline_o%ilin1, ifail)
120  call cdfgetvar(ncid, trim(zpre)//'x1pkg', spline_o%x1pkg, ifail)
121  call cdfgetvar(ncid, trim(zpre)//'fspl', spline_o%fspl, ifail)
122  endif
123 
124  if(ifail /=0) then
125  ier = 21
126  return
127  endif
128 
129  call cdfcls(ncid)
130 
131  if(.not.fullsv) then
132  call ezspline_setup1_r8x(spline_o, f, ier)
133  deallocate(f)
134  if(ifail /=0) ier = 22
135  endif
136 
137  end subroutine ezspline_load1_r8
138 
139 
140 
141 
142 
143 !!
144 !! 2-D
145 !!
146 
147 
148  subroutine ezspline_load2_r8(spline_o, filename, ier, spl_name)
149  use ezspline_obj
150  use ezcdf
151  implicit none
152  type(EZspline2_r8) :: spline_o
153  character*(*) :: filename
154  ! ier:
155  ! 20=attempt to load spline object with incorrect rank
156  ! 21=could not load spline object from file filename
157  ! 22=loaded spline object from file filename but failed at coefficient set-up
158  integer, intent(out) :: ier
159 
160  character*(*), intent(in), OPTIONAL :: spl_name ! specify name of spline
161  ! to be read-- in case file contains more than one.
162 
163  integer ncid, ifail, in0, in1, in2
164  integer dimlens(3)
165  character(2), parameter :: real8='R8'
166  character(3), parameter :: int='INT'
167  integer n1, n2, BCS1(2), BCS2(2), hspline(2)
168  real(ezspline_r8), dimension(:,:), allocatable :: f
169 
170  character*4 :: data_type
171 
172  character*32 zpre
173  logical :: fullsv
174 
175  if(present(spl_name)) then
176  zpre = spl_name
177  else
178  zpre = ' '
179  endif
180 
181  ier = 0
182  call cdfopn(ncid, filename, 'r') ! no error flag??
183  if(ncid==0) then
184  ier=43
185  return
186  endif
187 
188  ! check if n1 is present; check if spline object has correct rank
189  ! n2 should be present but not n3...
190 
191  call cdfinqvar(ncid, trim(zpre)//'n1', dimlens, data_type, ifail)
192  if(ifail.ne.0) then
193  ier = 21
194  return
195  endif
196  call cdfinqvar(ncid, trim(zpre)//'n2', dimlens, data_type, ifail)
197  if(ifail.ne.0) then
198  ier = 20
199  return
200  endif
201  call cdfinqvar(ncid, trim(zpre)//'n3', dimlens, data_type, ifail)
202  if(ifail.eq.0) then
203  ier = 20
204  return
205  endif
206 
207  ! check if fullsv was in effect -- if so, x1pkg will be found.
208 
209  fullsv=.false.
210  call cdfinqvar(ncid, trim(zpre)//'x1pkg', dimlens, data_type, ifail)
211  if(ifail.eq.0) then
212  fullsv=.true.
213  endif
214 
215  if(spline_o%nguard /= 123456789 ) call ezspline_preinit(spline_o)
216 
217  if( ezspline_allocated(spline_o) ) call ezspline_free2_r8(spline_o, ier)
218 
219  call cdfgetvar(ncid, trim(zpre)//'n1', n1, ifail)
220  call cdfgetvar(ncid, trim(zpre)//'n2', n2, ifail)
221  call cdfgetvar(ncid, trim(zpre)//'isLinear', spline_o%isLinear, ifail)
222  if(ifail.ne.0) spline_o%isLinear=0
223  call cdfgetvar(ncid, trim(zpre)//'isHybrid', spline_o%isHybrid, ifail)
224  if(ifail.ne.0) spline_o%isHybrid=0
225  spline_o%hspline = 0
226 
227  if(spline_o%isLinear.eq.1) then
228  call ezlinear_init2_r8(spline_o, n1, n2, ier)
229  else if(spline_o%isHybrid.eq.1) then
230  call cdfgetvar(ncid, trim(zpre)//'hspline', hspline, ifail)
231  call ezhybrid_init2_r8x(spline_o, n1, n2, hspline, ier)
232  else
233  bcs1 = (/0, 0/); bcs2 = (/0, 0/)
234  call ezspline_init2_r8(spline_o, n1, n2, bcs1, bcs2, ier)
235  endif
236  if(ier.ne.0) return
237 
238  call cdfgetvar(ncid, trim(zpre)//'klookup1', spline_o%klookup1, ifail)
239  if(ifail.ne.0) spline_o%klookup1=-3 ! the old default
240  call cdfgetvar(ncid, trim(zpre)//'klookup2', spline_o%klookup2, ifail)
241  if(ifail.ne.0) spline_o%klookup2=-3 ! the old default
242 
243  call cdfgetvar(ncid, trim(zpre)//'isHermite', spline_o%isHermite, ifail)
244  call cdfgetvar(ncid, trim(zpre)//'isReady', spline_o%isReady, ifail)
245  call cdfgetvar(ncid, trim(zpre)//'ibctype1', spline_o%ibctype1, ifail)
246  call cdfgetvar(ncid, trim(zpre)//'ibctype2', spline_o%ibctype2, ifail)
247  call cdfgetvar(ncid, trim(zpre)//'x1', spline_o%x1, ifail)
248  call cdfgetvar(ncid, trim(zpre)//'x2', spline_o%x2, ifail)
249  call cdfgetvar(ncid, trim(zpre)//'bcval1min', spline_o%bcval1min, ifail)
250  call cdfgetvar(ncid, trim(zpre)//'bcval1max', spline_o%bcval1max, ifail)
251  call cdfgetvar(ncid, trim(zpre)//'bcval2min', spline_o%bcval2min, ifail)
252  call cdfgetvar(ncid, trim(zpre)//'bcval2max', spline_o%bcval2max, ifail)
253  if(ifail/=0) then
254  ier=44
255  return
256  endif
257 
258  in0=size(spline_o%fspl,1)
259  in1=size(spline_o%fspl,2)
260  in2=size(spline_o%fspl,3)
261 
262  if(.not.fullsv) then
263  allocate(f(in1,in2))
264  if(spline_o%isReady==1) then
265  call cdfgetvar(ncid, trim(zpre)//'f', f, ifail)
266  else
267  ier = 92
268  endif
269  else
270  call cdfgetvar(ncid, trim(zpre)//'x1min', spline_o%x1min, ifail)
271  call cdfgetvar(ncid, trim(zpre)//'x1max', spline_o%x1max, ifail)
272  call cdfgetvar(ncid, trim(zpre)//'ilin1', spline_o%ilin1, ifail)
273  call cdfgetvar(ncid, trim(zpre)//'x2min', spline_o%x2min, ifail)
274  call cdfgetvar(ncid, trim(zpre)//'x2max', spline_o%x2max, ifail)
275  call cdfgetvar(ncid, trim(zpre)//'ilin2', spline_o%ilin2, ifail)
276  call cdfgetvar(ncid, trim(zpre)//'x1pkg', spline_o%x1pkg, ifail)
277  call cdfgetvar(ncid, trim(zpre)//'x2pkg', spline_o%x2pkg, ifail)
278  call cdfgetvar(ncid, trim(zpre)//'fspl', spline_o%fspl, ifail)
279  endif
280 
281  if(ifail /=0) then
282  ier = 21
283  return
284  endif
285 
286  call cdfcls(ncid)
287 
288  if(.not.fullsv) then
289  call ezspline_setup2_r8x(spline_o, f, ier)
290  deallocate(f)
291  if(ifail /=0) ier = 22
292  endif
293 
294  end subroutine ezspline_load2_r8
295 
296 
297 !!!
298 !!! 3-D
299 !!!
300 
301 
302  subroutine ezspline_load3_r8(spline_o, filename, ier, spl_name)
303  use ezspline_obj
304  use ezcdf
305  implicit none
306  type(EZspline3_r8) :: spline_o
307  character*(*) :: filename
308  ! ier:
309  ! 21=could not load spline object from file filename
310  ! 22=loaded spline object from file filename but failed at coefficient set-up
311  integer, intent(out) :: ier
312 
313  character*(*), intent(in), OPTIONAL :: spl_name ! specify name of spline
314  ! to be read-- in case file contains more than one.
315 
316  integer ncid, ifail, in0, in1, in2, in3
317  integer dimlens(3)
318  character(2), parameter :: real8='R8'
319  character(3), parameter :: int='INT'
320  integer n1, n2, n3, BCS1(2), BCS2(2), BCS3(2), hspline(3)
321  real(ezspline_r8), dimension(:,:,:), allocatable :: f
322 
323  character*4 :: data_type
324 
325  character*32 zpre
326  logical :: fullsv
327 
328  if(present(spl_name)) then
329  zpre = spl_name
330  else
331  zpre = ' '
332  endif
333 
334  ier = 0
335  call cdfopn(ncid, filename, 'r') ! no error flag??
336  if(ncid==0) then
337  ier=43
338  return
339  endif
340 
341  ! check if n1 is present; check if spline object has correct rank
342  ! n2 and n3 should also be present.
343 
344  call cdfinqvar(ncid, trim(zpre)//'n1', dimlens, data_type, ifail)
345  if(ifail.ne.0) then
346  ier = 21
347  return
348  endif
349  call cdfinqvar(ncid, trim(zpre)//'n2', dimlens, data_type, ifail)
350  if(ifail.ne.0) then
351  ier = 20
352  return
353  endif
354  call cdfinqvar(ncid, trim(zpre)//'n3', dimlens, data_type, ifail)
355  if(ifail.ne.0) then
356  ier = 20
357  return
358  endif
359 
360  ! check if fullsv was in effect -- if so, x1pkg will be found.
361 
362  fullsv=.false.
363  call cdfinqvar(ncid, trim(zpre)//'x1pkg', dimlens, data_type, ifail)
364  if(ifail.eq.0) then
365  fullsv=.true.
366  endif
367 
368  if(spline_o%nguard /= 123456789 ) call ezspline_preinit(spline_o)
369 
370  if( ezspline_allocated(spline_o) ) call ezspline_free3_r8(spline_o, ier)
371 
372  call cdfgetvar(ncid, trim(zpre)//'n1', n1, ifail)
373  call cdfgetvar(ncid, trim(zpre)//'n2', n2, ifail)
374  call cdfgetvar(ncid, trim(zpre)//'n3', n3, ifail)
375  call cdfgetvar(ncid, trim(zpre)//'isLinear', spline_o%isLinear, ifail)
376  if(ifail.ne.0) spline_o%isLinear=0
377  call cdfgetvar(ncid, trim(zpre)//'isHybrid', spline_o%isHybrid, ifail)
378  if(ifail.ne.0) spline_o%isHybrid=0
379  spline_o%hspline = 0
380 
381  if(spline_o%isLinear.eq.1) then
382  call ezlinear_init3_r8(spline_o, n1, n2, n3, ier)
383  else if(spline_o%isHybrid.eq.1) then
384  call cdfgetvar(ncid, trim(zpre)//'hspline', hspline, ifail)
385  call ezhybrid_init3_r8x(spline_o, n1, n2, n3, hspline, ier)
386  else
387  bcs1 = (/0, 0/); bcs2 = (/0, 0/); bcs3 = (/0, 0/)
388  call ezspline_init3_r8(spline_o, n1, n2, n3, bcs1, bcs2, bcs3, ier)
389  endif
390  if(ier.ne.0) return
391 
392  call cdfgetvar(ncid, trim(zpre)//'klookup1', spline_o%klookup1, ifail)
393  if(ifail.ne.0) spline_o%klookup1=-3 ! the old default
394  call cdfgetvar(ncid, trim(zpre)//'klookup2', spline_o%klookup2, ifail)
395  if(ifail.ne.0) spline_o%klookup2=-3 ! the old default
396  call cdfgetvar(ncid, trim(zpre)//'klookup3', spline_o%klookup3, ifail)
397  if(ifail.ne.0) spline_o%klookup3=-3 ! the old default
398 
399  call cdfgetvar(ncid, trim(zpre)//'isHermite', spline_o%isHermite, ifail)
400  call cdfgetvar(ncid, trim(zpre)//'isReady', spline_o%isReady, ifail)
401  call cdfgetvar(ncid, trim(zpre)//'ibctype1', spline_o%ibctype1, ifail)
402  call cdfgetvar(ncid, trim(zpre)//'ibctype2', spline_o%ibctype2, ifail)
403  call cdfgetvar(ncid, trim(zpre)//'ibctype3', spline_o%ibctype3, ifail)
404  call cdfgetvar(ncid, trim(zpre)//'x1', spline_o%x1, ifail)
405  call cdfgetvar(ncid, trim(zpre)//'x2', spline_o%x2, ifail)
406  call cdfgetvar(ncid, trim(zpre)//'x3', spline_o%x3, ifail)
407  call cdfgetvar(ncid, trim(zpre)//'bcval1min', spline_o%bcval1min, ifail)
408  call cdfgetvar(ncid, trim(zpre)//'bcval1max', spline_o%bcval1max, ifail)
409  call cdfgetvar(ncid, trim(zpre)//'bcval2min', spline_o%bcval2min, ifail)
410  call cdfgetvar(ncid, trim(zpre)//'bcval2max', spline_o%bcval2max, ifail)
411  call cdfgetvar(ncid, trim(zpre)//'bcval3min', spline_o%bcval3min, ifail)
412  call cdfgetvar(ncid, trim(zpre)//'bcval3max', spline_o%bcval3max, ifail)
413  if(ifail/=0) then
414  ier=44
415  return
416  endif
417 
418  in0=size(spline_o%fspl,1)
419  in1=size(spline_o%fspl,2)
420  in2=size(spline_o%fspl,3)
421  in3=size(spline_o%fspl,4)
422 
423  if(.not.fullsv) then
424  allocate(f(in1,in2,in3))
425  if(spline_o%isReady==1) then
426  call cdfgetvar(ncid, trim(zpre)//'f', f, ifail)
427  else
428  ier = 92
429  endif
430  else
431  call cdfgetvar(ncid, trim(zpre)//'x1min', spline_o%x1min, ifail)
432  call cdfgetvar(ncid, trim(zpre)//'x1max', spline_o%x1max, ifail)
433  call cdfgetvar(ncid, trim(zpre)//'ilin1', spline_o%ilin1, ifail)
434  call cdfgetvar(ncid, trim(zpre)//'x2min', spline_o%x2min, ifail)
435  call cdfgetvar(ncid, trim(zpre)//'x2max', spline_o%x2max, ifail)
436  call cdfgetvar(ncid, trim(zpre)//'ilin2', spline_o%ilin2, ifail)
437  call cdfgetvar(ncid, trim(zpre)//'x3min', spline_o%x3min, ifail)
438  call cdfgetvar(ncid, trim(zpre)//'x3max', spline_o%x3max, ifail)
439  call cdfgetvar(ncid, trim(zpre)//'ilin3', spline_o%ilin3, ifail)
440  call cdfgetvar(ncid, trim(zpre)//'x1pkg', spline_o%x1pkg, ifail)
441  call cdfgetvar(ncid, trim(zpre)//'x2pkg', spline_o%x2pkg, ifail)
442  call cdfgetvar(ncid, trim(zpre)//'x3pkg', spline_o%x3pkg, ifail)
443  call ezspline_cdfget3(ncid, trim(zpre)//'fspl', spline_o%fspl, &
444  in0*in1, in2, in3, ifail)
445  endif
446 
447  if(ifail /=0) then
448  ier = 21
449  return
450  endif
451 
452  call cdfcls(ncid)
453 
454  if(.not.fullsv) then
455  call ezspline_setup3_r8x(spline_o, f, ier)
456  deallocate(f)
457  if(ifail /=0) ier = 22
458  endif
459 
460  end subroutine ezspline_load3_r8
461 !/////
462 ! R4 !
463 !/////
464 !
465 ! 1-D
466 !
467 
468 
469  subroutine ezspline_load1_r4(spline_o, filename, ier, spl_name)
470  use ezspline_obj
471  use ezcdf
472  implicit none
473  type(EZspline1_r4) :: spline_o
474  character*(*) :: filename
475  ! ier:
476  ! 21=could not load spline object from file filename
477  ! 22=loaded spline object from file filename but failed at coefficient set-up
478  integer, intent(out) :: ier
479 
480  character*(*), intent(in), OPTIONAL :: spl_name ! specify name of spline
481  ! to be read-- in case file contains more than one.
482 
483  integer ncid, ifail
484  integer dimlens(3)
485  character(2), parameter :: real8='R8'
486  character(3), parameter :: int='INT'
487  integer n1, BCS1(2)
488  real(ezspline_r4), dimension(:), allocatable :: f
489 
490  character*4 :: data_type
491 
492  character*32 zpre
493  logical :: fullsv
494 
495  if(present(spl_name)) then
496  zpre = spl_name
497  else
498  zpre = ' '
499  endif
500 
501  ier = 0
502  call cdfopn(ncid, filename, 'r') ! no error flag??
503  if(ncid==0) then
504  ier=43
505  return
506  endif
507 
508  ! check if n1 is present; check if spline object has correct rank
509  ! n2 should NOT be present.
510 
511  call cdfinqvar(ncid, trim(zpre)//'n1', dimlens, data_type, ifail)
512  if(ifail.ne.0) then
513  ier = 21
514  return
515  endif
516  call cdfinqvar(ncid, trim(zpre)//'n2', dimlens, data_type, ifail)
517  if(ifail.eq.0) then
518  ier = 20
519  return
520  endif
521 
522  ! check if fullsv was in effect -- if so, x1pkg will be found.
523 
524  fullsv=.false.
525  call cdfinqvar(ncid, trim(zpre)//'x1pkg', dimlens, data_type, ifail)
526  if(ifail.eq.0) then
527  fullsv=.true.
528  endif
529 
530  if(spline_o%nguard /= 123456789 ) call ezspline_preinit(spline_o)
531 
532  if( ezspline_allocated(spline_o) ) call ezspline_free1_r4(spline_o, ier)
533 
534  call cdfgetvar(ncid, trim(zpre)//'n1', n1, ifail)
535  call cdfgetvar(ncid, trim(zpre)//'isLinear', spline_o%isLinear, ifail)
536  if(ifail.ne.0) spline_o%isLinear=0
537 
538  if(spline_o%isLinear.eq.1) then
539  call ezlinear_init1_r4(spline_o, n1, ier)
540  else
541  bcs1 = (/0, 0/)
542  call ezspline_init1_r4(spline_o, n1, bcs1, ier)
543  endif
544  if(ier.ne.0) return
545 
546  call cdfgetvar(ncid, trim(zpre)//'klookup1', spline_o%klookup1, ifail)
547  if(ifail.ne.0) spline_o%klookup1=-3 ! the old default
548 
549  call cdfgetvar(ncid, trim(zpre)//'isHermite', spline_o%isHermite, ifail)
550  call cdfgetvar(ncid, trim(zpre)//'isReady', spline_o%isReady, ifail)
551  call cdfgetvar(ncid, trim(zpre)//'ibctype1', spline_o%ibctype1, ifail)
552  call cdfgetvar(ncid, trim(zpre)//'x1', spline_o%x1, ifail)
553  call cdfgetvar(ncid, trim(zpre)//'bcval1min', spline_o%bcval1min, ifail)
554  call cdfgetvar(ncid, trim(zpre)//'bcval1max', spline_o%bcval1max, ifail)
555  if(ifail/=0) then
556  ier=44
557  return
558  endif
559 
560  if(.not.fullsv) then
561  allocate(f(n1))
562  if(spline_o%isReady==1) then
563  call cdfgetvar(ncid, trim(zpre)//'f', f, ifail)
564  else
565  ier = 92
566  endif
567  else
568  call cdfgetvar(ncid, trim(zpre)//'x1min', spline_o%x1min, ifail)
569  call cdfgetvar(ncid, trim(zpre)//'x1max', spline_o%x1max, ifail)
570  call cdfgetvar(ncid, trim(zpre)//'ilin1', spline_o%ilin1, ifail)
571  call cdfgetvar(ncid, trim(zpre)//'x1pkg', spline_o%x1pkg, ifail)
572  call cdfgetvar(ncid, trim(zpre)//'fspl', spline_o%fspl, ifail)
573  endif
574 
575  if(ifail /=0) then
576  ier = 21
577  return
578  endif
579 
580  call cdfcls(ncid)
581 
582  if(.not.fullsv) then
583  call ezspline_setup1_r4x(spline_o, f, ier)
584  deallocate(f)
585  if(ifail /=0) ier = 22
586  endif
587 
588  end subroutine ezspline_load1_r4
589 
590 
591 
592 
593 
594 !!
595 !! 2-D
596 !!
597 
598 
599  subroutine ezspline_load2_r4(spline_o, filename, ier, spl_name)
600  use ezspline_obj
601  use ezcdf
602  implicit none
603  type(EZspline2_r4) :: spline_o
604  character*(*) :: filename
605  ! ier:
606  ! 20=attempt to load spline object with incorrect rank
607  ! 21=could not load spline object from file filename
608  ! 22=loaded spline object from file filename but failed at coefficient set-up
609  integer, intent(out) :: ier
610 
611  character*(*), intent(in), OPTIONAL :: spl_name ! specify name of spline
612  ! to be read-- in case file contains more than one.
613 
614  integer ncid, ifail, in0, in1, in2
615  integer dimlens(3)
616  character(2), parameter :: real8='R8'
617  character(3), parameter :: int='INT'
618  integer n1, n2, BCS1(2), BCS2(2), hspline(2)
619  real(ezspline_r4), dimension(:,:), allocatable :: f
620 
621  character*4 :: data_type
622 
623  character*32 zpre
624  logical :: fullsv
625 
626  if(present(spl_name)) then
627  zpre = spl_name
628  else
629  zpre = ' '
630  endif
631 
632  ier = 0
633  call cdfopn(ncid, filename, 'r') ! no error flag??
634  if(ncid==0) then
635  ier=43
636  return
637  endif
638 
639  ! check if n1 is present; check if spline object has correct rank
640  ! n2 should be present but not n3...
641 
642  call cdfinqvar(ncid, trim(zpre)//'n1', dimlens, data_type, ifail)
643  if(ifail.ne.0) then
644  ier = 21
645  return
646  endif
647  call cdfinqvar(ncid, trim(zpre)//'n2', dimlens, data_type, ifail)
648  if(ifail.ne.0) then
649  ier = 20
650  return
651  endif
652  call cdfinqvar(ncid, trim(zpre)//'n3', dimlens, data_type, ifail)
653  if(ifail.eq.0) then
654  ier = 20
655  return
656  endif
657 
658  ! check if fullsv was in effect -- if so, x1pkg will be found.
659 
660  fullsv=.false.
661  call cdfinqvar(ncid, trim(zpre)//'x1pkg', dimlens, data_type, ifail)
662  if(ifail.eq.0) then
663  fullsv=.true.
664  endif
665 
666  if(spline_o%nguard /= 123456789 ) call ezspline_preinit(spline_o)
667 
668  if( ezspline_allocated(spline_o) ) call ezspline_free2_r4(spline_o, ier)
669 
670  call cdfgetvar(ncid, trim(zpre)//'n1', n1, ifail)
671  call cdfgetvar(ncid, trim(zpre)//'n2', n2, ifail)
672  call cdfgetvar(ncid, trim(zpre)//'isLinear', spline_o%isLinear, ifail)
673  if(ifail.ne.0) spline_o%isLinear=0
674  call cdfgetvar(ncid, trim(zpre)//'isHybrid', spline_o%isHybrid, ifail)
675  if(ifail.ne.0) spline_o%isHybrid=0
676  spline_o%hspline = 0
677 
678  if(spline_o%isLinear.eq.1) then
679  call ezlinear_init2_r4(spline_o, n1, n2, ier)
680  else if(spline_o%isHybrid.eq.1) then
681  call cdfgetvar(ncid, trim(zpre)//'hspline', hspline, ifail)
682  call ezhybrid_init2_r4x(spline_o, n1, n2, hspline, ier)
683  else
684  bcs1 = (/0, 0/); bcs2 = (/0, 0/)
685  call ezspline_init2_r4(spline_o, n1, n2, bcs1, bcs2, ier)
686  endif
687  if(ier.ne.0) return
688 
689  call cdfgetvar(ncid, trim(zpre)//'klookup1', spline_o%klookup1, ifail)
690  if(ifail.ne.0) spline_o%klookup1=-3 ! the old default
691  call cdfgetvar(ncid, trim(zpre)//'klookup2', spline_o%klookup2, ifail)
692  if(ifail.ne.0) spline_o%klookup2=-3 ! the old default
693 
694  call cdfgetvar(ncid, trim(zpre)//'isHermite', spline_o%isHermite, ifail)
695  call cdfgetvar(ncid, trim(zpre)//'isReady', spline_o%isReady, ifail)
696  call cdfgetvar(ncid, trim(zpre)//'ibctype1', spline_o%ibctype1, ifail)
697  call cdfgetvar(ncid, trim(zpre)//'ibctype2', spline_o%ibctype2, ifail)
698  call cdfgetvar(ncid, trim(zpre)//'x1', spline_o%x1, ifail)
699  call cdfgetvar(ncid, trim(zpre)//'x2', spline_o%x2, ifail)
700  call cdfgetvar(ncid, trim(zpre)//'bcval1min', spline_o%bcval1min, ifail)
701  call cdfgetvar(ncid, trim(zpre)//'bcval1max', spline_o%bcval1max, ifail)
702  call cdfgetvar(ncid, trim(zpre)//'bcval2min', spline_o%bcval2min, ifail)
703  call cdfgetvar(ncid, trim(zpre)//'bcval2max', spline_o%bcval2max, ifail)
704  if(ifail/=0) then
705  ier=44
706  return
707  endif
708 
709  in0=size(spline_o%fspl,1)
710  in1=size(spline_o%fspl,2)
711  in2=size(spline_o%fspl,3)
712 
713  if(.not.fullsv) then
714  allocate(f(in1,in2))
715  if(spline_o%isReady==1) then
716  call cdfgetvar(ncid, trim(zpre)//'f', f, ifail)
717  else
718  ier = 92
719  endif
720  else
721  call cdfgetvar(ncid, trim(zpre)//'x1min', spline_o%x1min, ifail)
722  call cdfgetvar(ncid, trim(zpre)//'x1max', spline_o%x1max, ifail)
723  call cdfgetvar(ncid, trim(zpre)//'ilin1', spline_o%ilin1, ifail)
724  call cdfgetvar(ncid, trim(zpre)//'x2min', spline_o%x2min, ifail)
725  call cdfgetvar(ncid, trim(zpre)//'x2max', spline_o%x2max, ifail)
726  call cdfgetvar(ncid, trim(zpre)//'ilin2', spline_o%ilin2, ifail)
727  call cdfgetvar(ncid, trim(zpre)//'x1pkg', spline_o%x1pkg, ifail)
728  call cdfgetvar(ncid, trim(zpre)//'x2pkg', spline_o%x2pkg, ifail)
729  call cdfgetvar(ncid, trim(zpre)//'fspl', spline_o%fspl, ifail)
730  endif
731 
732  if(ifail /=0) then
733  ier = 21
734  return
735  endif
736 
737  call cdfcls(ncid)
738 
739  if(.not.fullsv) then
740  call ezspline_setup2_r4x(spline_o, f, ier)
741  deallocate(f)
742  if(ifail /=0) ier = 22
743  endif
744 
745  end subroutine ezspline_load2_r4
746 
747 
748 !!!
749 !!! 3-D
750 !!!
751 
752 
753  subroutine ezspline_load3_r4(spline_o, filename, ier, spl_name)
754  use ezspline_obj
755  use ezcdf
756  implicit none
757  type(EZspline3_r4) :: spline_o
758  character*(*) :: filename
759  ! ier:
760  ! 21=could not load spline object from file filename
761  ! 22=loaded spline object from file filename but failed at coefficient set-up
762  integer, intent(out) :: ier
763 
764  character*(*), intent(in), OPTIONAL :: spl_name ! specify name of spline
765  ! to be read-- in case file contains more than one.
766 
767  integer ncid, ifail, in0, in1, in2, in3
768  integer dimlens(3)
769  character(2), parameter :: real8='R8'
770  character(3), parameter :: int='INT'
771  integer n1, n2, n3, BCS1(2), BCS2(2), BCS3(2), hspline(3)
772  real(ezspline_r4), dimension(:,:,:), allocatable :: f
773 
774  character*4 :: data_type
775 
776  character*32 zpre
777  logical :: fullsv
778 
779  if(present(spl_name)) then
780  zpre = spl_name
781  else
782  zpre = ' '
783  endif
784 
785  ier = 0
786  call cdfopn(ncid, filename, 'r') ! no error flag??
787  if(ncid==0) then
788  ier=43
789  return
790  endif
791 
792  ! check if n1 is present; check if spline object has correct rank
793  ! n2 and n3 should also be present.
794 
795  call cdfinqvar(ncid, trim(zpre)//'n1', dimlens, data_type, ifail)
796  if(ifail.ne.0) then
797  ier = 21
798  return
799  endif
800  call cdfinqvar(ncid, trim(zpre)//'n2', dimlens, data_type, ifail)
801  if(ifail.ne.0) then
802  ier = 20
803  return
804  endif
805  call cdfinqvar(ncid, trim(zpre)//'n3', dimlens, data_type, ifail)
806  if(ifail.ne.0) then
807  ier = 20
808  return
809  endif
810 
811  ! check if fullsv was in effect -- if so, x1pkg will be found.
812 
813  fullsv=.false.
814  call cdfinqvar(ncid, trim(zpre)//'x1pkg', dimlens, data_type, ifail)
815  if(ifail.eq.0) then
816  fullsv=.true.
817  endif
818 
819  if(spline_o%nguard /= 123456789 ) call ezspline_preinit(spline_o)
820 
821  if( ezspline_allocated(spline_o) ) call ezspline_free3_r4(spline_o, ier)
822 
823  call cdfgetvar(ncid, trim(zpre)//'n1', n1, ifail)
824  call cdfgetvar(ncid, trim(zpre)//'n2', n2, ifail)
825  call cdfgetvar(ncid, trim(zpre)//'n3', n3, ifail)
826  call cdfgetvar(ncid, trim(zpre)//'isLinear', spline_o%isLinear, ifail)
827  if(ifail.ne.0) spline_o%isLinear=0
828  call cdfgetvar(ncid, trim(zpre)//'isHybrid', spline_o%isHybrid, ifail)
829  if(ifail.ne.0) spline_o%isHybrid=0
830  spline_o%hspline = 0
831 
832  if(spline_o%isLinear.eq.1) then
833  call ezlinear_init3_r4(spline_o, n1, n2, n3, ier)
834  else if(spline_o%isHybrid.eq.1) then
835  call cdfgetvar(ncid, trim(zpre)//'hspline', hspline, ifail)
836  call ezhybrid_init3_r4x(spline_o, n1, n2, n3, hspline, ier)
837  else
838  bcs1 = (/0, 0/); bcs2 = (/0, 0/); bcs3 = (/0, 0/)
839  call ezspline_init3_r4(spline_o, n1, n2, n3, bcs1, bcs2, bcs3, ier)
840  endif
841  if(ier.ne.0) return
842 
843  call cdfgetvar(ncid, trim(zpre)//'klookup1', spline_o%klookup1, ifail)
844  if(ifail.ne.0) spline_o%klookup1=-3 ! the old default
845  call cdfgetvar(ncid, trim(zpre)//'klookup2', spline_o%klookup2, ifail)
846  if(ifail.ne.0) spline_o%klookup2=-3 ! the old default
847  call cdfgetvar(ncid, trim(zpre)//'klookup3', spline_o%klookup3, ifail)
848  if(ifail.ne.0) spline_o%klookup3=-3 ! the old default
849 
850  call cdfgetvar(ncid, trim(zpre)//'isHermite', spline_o%isHermite, ifail)
851  call cdfgetvar(ncid, trim(zpre)//'isReady', spline_o%isReady, ifail)
852  call cdfgetvar(ncid, trim(zpre)//'ibctype1', spline_o%ibctype1, ifail)
853  call cdfgetvar(ncid, trim(zpre)//'ibctype2', spline_o%ibctype2, ifail)
854  call cdfgetvar(ncid, trim(zpre)//'ibctype3', spline_o%ibctype3, ifail)
855  call cdfgetvar(ncid, trim(zpre)//'x1', spline_o%x1, ifail)
856  call cdfgetvar(ncid, trim(zpre)//'x2', spline_o%x2, ifail)
857  call cdfgetvar(ncid, trim(zpre)//'x3', spline_o%x3, ifail)
858  call cdfgetvar(ncid, trim(zpre)//'bcval1min', spline_o%bcval1min, ifail)
859  call cdfgetvar(ncid, trim(zpre)//'bcval1max', spline_o%bcval1max, ifail)
860  call cdfgetvar(ncid, trim(zpre)//'bcval2min', spline_o%bcval2min, ifail)
861  call cdfgetvar(ncid, trim(zpre)//'bcval2max', spline_o%bcval2max, ifail)
862  call cdfgetvar(ncid, trim(zpre)//'bcval3min', spline_o%bcval3min, ifail)
863  call cdfgetvar(ncid, trim(zpre)//'bcval3max', spline_o%bcval3max, ifail)
864  if(ifail/=0) then
865  ier=44
866  return
867  endif
868 
869  in0=size(spline_o%fspl,1)
870  in1=size(spline_o%fspl,2)
871  in2=size(spline_o%fspl,3)
872  in3=size(spline_o%fspl,4)
873 
874  if(.not.fullsv) then
875  allocate(f(in1,in2,in3))
876  if(spline_o%isReady==1) then
877  call cdfgetvar(ncid, trim(zpre)//'f', f, ifail)
878  else
879  ier = 92
880  endif
881  else
882  call cdfgetvar(ncid, trim(zpre)//'x1min', spline_o%x1min, ifail)
883  call cdfgetvar(ncid, trim(zpre)//'x1max', spline_o%x1max, ifail)
884  call cdfgetvar(ncid, trim(zpre)//'ilin1', spline_o%ilin1, ifail)
885  call cdfgetvar(ncid, trim(zpre)//'x2min', spline_o%x2min, ifail)
886  call cdfgetvar(ncid, trim(zpre)//'x2max', spline_o%x2max, ifail)
887  call cdfgetvar(ncid, trim(zpre)//'ilin2', spline_o%ilin2, ifail)
888  call cdfgetvar(ncid, trim(zpre)//'x3min', spline_o%x3min, ifail)
889  call cdfgetvar(ncid, trim(zpre)//'x3max', spline_o%x3max, ifail)
890  call cdfgetvar(ncid, trim(zpre)//'ilin3', spline_o%ilin3, ifail)
891  call cdfgetvar(ncid, trim(zpre)//'x1pkg', spline_o%x1pkg, ifail)
892  call cdfgetvar(ncid, trim(zpre)//'x2pkg', spline_o%x2pkg, ifail)
893  call cdfgetvar(ncid, trim(zpre)//'x3pkg', spline_o%x3pkg, ifail)
894  call ezspline_cdfget3(ncid, trim(zpre)//'fspl', spline_o%fspl, &
895  in0*in1, in2, in3, ifail)
896  endif
897 
898  if(ifail /=0) then
899  ier = 21
900  return
901  endif
902 
903  call cdfcls(ncid)
904 
905  if(.not.fullsv) then
906  call ezspline_setup3_r4x(spline_o, f, ier)
907  deallocate(f)
908  if(ifail /=0) ier = 22
909  endif
910 
911  end subroutine ezspline_load3_r4
ezspline_obj::ezspline_allocated
Definition: ezspline_obj.f90:19
ezspline_obj::ezspline_preinit
Definition: ezspline_obj.f90:3