1 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
2 | ! |
---|
3 | ! this program is a short driver that tests the remappings using |
---|
4 | ! a simple analytic field. the results are written in netCDF |
---|
5 | ! format. |
---|
6 | ! |
---|
7 | !----------------------------------------------------------------------- |
---|
8 | ! |
---|
9 | ! CVS:$Id: scrip_test_repeat.f,v 1.3 2000/04/19 21:56:26 pwjones Exp $ |
---|
10 | ! |
---|
11 | ! Copyright (c) 1997, 1998 the Regents of the University of |
---|
12 | ! California. |
---|
13 | ! |
---|
14 | ! This software and ancillary information (herein called software) |
---|
15 | ! called SCRIP is made available under the terms described here. |
---|
16 | ! The software has been approved for release with associated |
---|
17 | ! LA-CC Number 98-45. |
---|
18 | ! |
---|
19 | ! Unless otherwise indicated, this software has been authored |
---|
20 | ! by an employee or employees of the University of California, |
---|
21 | ! operator of the Los Alamos National Laboratory under Contract |
---|
22 | ! No. W-7405-ENG-36 with the U.S. Department of Energy. The U.S. |
---|
23 | ! Government has rights to use, reproduce, and distribute this |
---|
24 | ! software. The public may copy and use this software without |
---|
25 | ! charge, provided that this Notice and any statement of authorship |
---|
26 | ! are reproduced on all copies. Neither the Government nor the |
---|
27 | ! University makes any warranty, express or implied, or assumes |
---|
28 | ! any liability or responsibility for the use of this software. |
---|
29 | ! |
---|
30 | ! If software is modified to produce derivative works, such modified |
---|
31 | ! software should be clearly marked, so as not to confuse it with |
---|
32 | ! the version available from Los Alamos National Laboratory. |
---|
33 | ! |
---|
34 | !*********************************************************************** |
---|
35 | |
---|
36 | program remap_test_repeat |
---|
37 | |
---|
38 | !----------------------------------------------------------------------- |
---|
39 | |
---|
40 | use kinds_mod ! defines common data types |
---|
41 | use constants ! defines common constants |
---|
42 | use iounits ! I/O unit manager |
---|
43 | use netcdf_mod ! netcdf I/O stuff |
---|
44 | use grids ! module containing grid info |
---|
45 | use remap_vars ! module containing remapping info |
---|
46 | |
---|
47 | implicit none |
---|
48 | |
---|
49 | !----------------------------------------------------------------------- |
---|
50 | ! |
---|
51 | ! interface for remap routine |
---|
52 | ! |
---|
53 | !----------------------------------------------------------------------- |
---|
54 | |
---|
55 | interface |
---|
56 | subroutine remap(dst_array, map_wts, dst_add, src_add, |
---|
57 | & src_array, src_grad1, src_grad2, src_grad3) |
---|
58 | |
---|
59 | use kinds_mod |
---|
60 | use constants |
---|
61 | |
---|
62 | implicit none |
---|
63 | |
---|
64 | integer (kind=int_kind), dimension(:), intent(in) :: |
---|
65 | & dst_add, src_add |
---|
66 | |
---|
67 | real (kind=dbl_kind), dimension(:,:), intent(in) :: map_wts |
---|
68 | |
---|
69 | real (kind=dbl_kind), dimension(:), intent(in) :: src_array |
---|
70 | |
---|
71 | real (kind=dbl_kind), dimension(:), intent(in), optional :: |
---|
72 | & src_grad1, src_grad2, src_grad3 |
---|
73 | |
---|
74 | real (kind=dbl_kind), dimension(:), intent(inout) :: dst_array |
---|
75 | |
---|
76 | end subroutine remap |
---|
77 | end interface |
---|
78 | |
---|
79 | !----------------------------------------------------------------------- |
---|
80 | ! |
---|
81 | ! input namelist variables |
---|
82 | ! |
---|
83 | !----------------------------------------------------------------------- |
---|
84 | |
---|
85 | integer (kind=int_kind) :: |
---|
86 | & field_choice, ! choice of field to be interpolated |
---|
87 | & num_repeats ! number of times to repeat remappings |
---|
88 | |
---|
89 | character (char_len) :: |
---|
90 | & interp_file1, ! filename containing remap data (map1) |
---|
91 | & interp_file2, ! filename containing remap data (map2) |
---|
92 | & output_file ! filename containing output test data |
---|
93 | |
---|
94 | namelist /remap_inputs/ field_choice, num_repeats, |
---|
95 | & interp_file1, interp_file2, output_file |
---|
96 | |
---|
97 | !----------------------------------------------------------------------- |
---|
98 | ! |
---|
99 | ! local variables |
---|
100 | ! |
---|
101 | !----------------------------------------------------------------------- |
---|
102 | |
---|
103 | character (char_len) :: |
---|
104 | & map_name1, ! name for mapping from grid1 to grid2 |
---|
105 | & map_name2 ! name for mapping from grid2 to grid1 |
---|
106 | |
---|
107 | integer (kind=int_kind) :: ! netCDF ids for files and arrays |
---|
108 | & n, ncstat, nc_outfile_id, |
---|
109 | & nc_srcgrdcntrlat_id, nc_srcgrdcntrlon_id, |
---|
110 | & nc_dstgrdcntrlat_id, nc_dstgrdcntrlon_id, |
---|
111 | & nc_srcgrdrank_id, nc_dstgrdrank_id, |
---|
112 | & nc_srcgrdimask_id, nc_dstgrdimask_id, |
---|
113 | & nc_srcgrdarea_id, nc_dstgrdarea_id, |
---|
114 | & nc_srcgrdfrac_id, nc_dstgrdfrac_id, |
---|
115 | & nc_srcarray_id, nc_srcgradlat_id, nc_srcgradlon_id, |
---|
116 | & nc_dstarray1_id, nc_dstarray2_id, |
---|
117 | & nc_dsterror1_id, nc_dsterror2_id |
---|
118 | |
---|
119 | integer (kind=int_kind), dimension(:), allocatable :: |
---|
120 | & nc_grid1size_id, nc_grid2size_id |
---|
121 | |
---|
122 | !----------------------------------------------------------------------- |
---|
123 | |
---|
124 | character (char_len) :: |
---|
125 | & dim_name ! netCDF dimension name |
---|
126 | |
---|
127 | integer (kind=int_kind) :: i,j,n, |
---|
128 | & iunit ! unit number for namelist file |
---|
129 | |
---|
130 | integer (kind=int_kind), dimension(:), allocatable :: |
---|
131 | & grid1_imask, grid2_imask |
---|
132 | |
---|
133 | real (kind=dbl_kind) :: |
---|
134 | & length ! length scale for cosine hill test field |
---|
135 | |
---|
136 | real (kind=dbl_kind), dimension(:), allocatable :: |
---|
137 | & grid1_array, |
---|
138 | & grid1_tmp, |
---|
139 | & grad1_lat, |
---|
140 | & grad1_lon, |
---|
141 | & grid2_array, |
---|
142 | & grid2_err, |
---|
143 | & grid2_tmp, |
---|
144 | & grad2_lat, |
---|
145 | & grad2_lon |
---|
146 | |
---|
147 | !----------------------------------------------------------------------- |
---|
148 | ! |
---|
149 | ! read namelist for file and mapping info |
---|
150 | ! |
---|
151 | !----------------------------------------------------------------------- |
---|
152 | |
---|
153 | call get_unit(iunit) |
---|
154 | open(iunit, file='repeat_test_in', status='old', form='formatted') |
---|
155 | read(iunit, nml=remap_inputs) |
---|
156 | call release_unit(iunit) |
---|
157 | write(*,nml=remap_inputs) |
---|
158 | |
---|
159 | !----------------------------------------------------------------------- |
---|
160 | ! |
---|
161 | ! read remapping data |
---|
162 | ! |
---|
163 | !----------------------------------------------------------------------- |
---|
164 | |
---|
165 | num_maps = 2 |
---|
166 | |
---|
167 | call read_remap(map_name1, map_name2, interp_file1, interp_file2) |
---|
168 | |
---|
169 | !----------------------------------------------------------------------- |
---|
170 | ! |
---|
171 | ! allocate arrays |
---|
172 | ! |
---|
173 | !----------------------------------------------------------------------- |
---|
174 | |
---|
175 | allocate (grid1_array (grid1_size), |
---|
176 | & grid1_tmp (grid1_size), |
---|
177 | & grad1_lat (grid1_size), |
---|
178 | & grad1_lon (grid1_size), |
---|
179 | & grid1_imask (grid1_size), |
---|
180 | & grid2_array (grid2_size), |
---|
181 | & grid2_err (grid2_size), |
---|
182 | & grid2_tmp (grid2_size), |
---|
183 | & grad2_lat (grid2_size), |
---|
184 | & grad2_lon (grid2_size), |
---|
185 | & grid2_imask (grid2_size)) |
---|
186 | |
---|
187 | where (grid1_mask) |
---|
188 | grid1_imask = 1 |
---|
189 | elsewhere |
---|
190 | grid1_imask = 0 |
---|
191 | endwhere |
---|
192 | where (grid2_mask) |
---|
193 | grid2_imask = 1 |
---|
194 | elsewhere |
---|
195 | grid2_imask = 0 |
---|
196 | endwhere |
---|
197 | |
---|
198 | !----------------------------------------------------------------------- |
---|
199 | ! |
---|
200 | ! setup a NetCDF file for output |
---|
201 | ! |
---|
202 | !----------------------------------------------------------------------- |
---|
203 | |
---|
204 | !*** |
---|
205 | !*** create netCDF dataset |
---|
206 | !*** |
---|
207 | |
---|
208 | ncstat = nf_create (output_file, NF_CLOBBER, nc_outfile_id) |
---|
209 | call netcdf_error_handler(ncstat) |
---|
210 | |
---|
211 | ncstat = nf_put_att_text (nc_outfile_id, NF_GLOBAL, 'title', |
---|
212 | & len_trim(map_name1), map_name1) |
---|
213 | call netcdf_error_handler(ncstat) |
---|
214 | |
---|
215 | !*** |
---|
216 | !*** define grid size dimensions |
---|
217 | !*** |
---|
218 | |
---|
219 | allocate( nc_grid1size_id(grid1_rank), |
---|
220 | & nc_grid2size_id(grid2_rank)) |
---|
221 | |
---|
222 | do n=1,grid1_rank |
---|
223 | write(dim_name,1000) 'grid1_dim',n |
---|
224 | ncstat = nf_def_dim (nc_outfile_id, dim_name, |
---|
225 | & grid1_dims(n), nc_grid1size_id(n)) |
---|
226 | call netcdf_error_handler(ncstat) |
---|
227 | end do |
---|
228 | |
---|
229 | do n=1,grid2_rank |
---|
230 | write(dim_name,1000) 'grid2_dim',n |
---|
231 | ncstat = nf_def_dim (nc_outfile_id, dim_name, |
---|
232 | & grid2_dims(n), nc_grid2size_id(n)) |
---|
233 | call netcdf_error_handler(ncstat) |
---|
234 | end do |
---|
235 | 1000 format(a9,i1) |
---|
236 | |
---|
237 | !*** |
---|
238 | !*** define grid center latitude array |
---|
239 | !*** |
---|
240 | |
---|
241 | ncstat = nf_def_var (nc_outfile_id, 'src_grid_center_lat', |
---|
242 | & NF_DOUBLE, grid1_rank, nc_grid1size_id, |
---|
243 | & nc_srcgrdcntrlat_id) |
---|
244 | call netcdf_error_handler(ncstat) |
---|
245 | |
---|
246 | ncstat = nf_put_att_text (nc_outfile_id, nc_srcgrdcntrlat_id, |
---|
247 | & 'units', 7, 'radians') |
---|
248 | call netcdf_error_handler(ncstat) |
---|
249 | |
---|
250 | ncstat = nf_def_var (nc_outfile_id, 'dst_grid_center_lat', |
---|
251 | & NF_DOUBLE, grid2_rank, nc_grid2size_id, |
---|
252 | & nc_dstgrdcntrlat_id) |
---|
253 | call netcdf_error_handler(ncstat) |
---|
254 | |
---|
255 | ncstat = nf_put_att_text (nc_outfile_id, nc_dstgrdcntrlat_id, |
---|
256 | & 'units', 7, 'radians') |
---|
257 | call netcdf_error_handler(ncstat) |
---|
258 | |
---|
259 | !*** |
---|
260 | !*** define grid center longitude array |
---|
261 | !*** |
---|
262 | |
---|
263 | ncstat = nf_def_var (nc_outfile_id, 'src_grid_center_lon', |
---|
264 | & NF_DOUBLE, grid1_rank, nc_grid1size_id, |
---|
265 | & nc_srcgrdcntrlon_id) |
---|
266 | call netcdf_error_handler(ncstat) |
---|
267 | |
---|
268 | ncstat = nf_put_att_text (nc_outfile_id, nc_srcgrdcntrlon_id, |
---|
269 | & 'units', 7, 'radians') |
---|
270 | call netcdf_error_handler(ncstat) |
---|
271 | |
---|
272 | ncstat = nf_def_var (nc_outfile_id, 'dst_grid_center_lon', |
---|
273 | & NF_DOUBLE, grid2_rank, nc_grid2size_id, |
---|
274 | & nc_dstgrdcntrlon_id) |
---|
275 | call netcdf_error_handler(ncstat) |
---|
276 | |
---|
277 | ncstat = nf_put_att_text (nc_outfile_id, nc_dstgrdcntrlon_id, |
---|
278 | & 'units', 7, 'radians') |
---|
279 | call netcdf_error_handler(ncstat) |
---|
280 | |
---|
281 | !*** |
---|
282 | !*** define grid mask |
---|
283 | !*** |
---|
284 | |
---|
285 | ncstat = nf_def_var (nc_outfile_id, 'src_grid_imask', NF_INT, |
---|
286 | & grid1_rank, nc_grid1size_id, nc_srcgrdimask_id) |
---|
287 | call netcdf_error_handler(ncstat) |
---|
288 | |
---|
289 | ncstat = nf_put_att_text (nc_outfile_id, nc_srcgrdimask_id, |
---|
290 | & 'units', 8, 'unitless') |
---|
291 | call netcdf_error_handler(ncstat) |
---|
292 | |
---|
293 | ncstat = nf_def_var (nc_outfile_id, 'dst_grid_imask', NF_INT, |
---|
294 | & grid2_rank, nc_grid2size_id, nc_dstgrdimask_id) |
---|
295 | call netcdf_error_handler(ncstat) |
---|
296 | |
---|
297 | ncstat = nf_put_att_text (nc_outfile_id, nc_dstgrdimask_id, |
---|
298 | & 'units', 8, 'unitless') |
---|
299 | call netcdf_error_handler(ncstat) |
---|
300 | |
---|
301 | !*** |
---|
302 | !*** define grid area arrays |
---|
303 | !*** |
---|
304 | |
---|
305 | ncstat = nf_def_var (nc_outfile_id, 'src_grid_area', |
---|
306 | & NF_DOUBLE, grid1_rank, nc_grid1size_id, |
---|
307 | & nc_srcgrdarea_id) |
---|
308 | call netcdf_error_handler(ncstat) |
---|
309 | |
---|
310 | ncstat = nf_def_var (nc_outfile_id, 'dst_grid_area', |
---|
311 | & NF_DOUBLE, grid2_rank, nc_grid2size_id, |
---|
312 | & nc_dstgrdarea_id) |
---|
313 | call netcdf_error_handler(ncstat) |
---|
314 | |
---|
315 | !*** |
---|
316 | !*** define grid fraction arrays |
---|
317 | !*** |
---|
318 | |
---|
319 | ncstat = nf_def_var (nc_outfile_id, 'src_grid_frac', |
---|
320 | & NF_DOUBLE, grid1_rank, nc_grid1size_id, |
---|
321 | & nc_srcgrdfrac_id) |
---|
322 | call netcdf_error_handler(ncstat) |
---|
323 | |
---|
324 | ncstat = nf_def_var (nc_outfile_id, 'dst_grid_frac', |
---|
325 | & NF_DOUBLE, grid2_rank, nc_grid2size_id, |
---|
326 | & nc_dstgrdfrac_id) |
---|
327 | call netcdf_error_handler(ncstat) |
---|
328 | |
---|
329 | !*** |
---|
330 | !*** define source array |
---|
331 | !*** |
---|
332 | |
---|
333 | ncstat = nf_def_var (nc_outfile_id, 'src_array', |
---|
334 | & NF_DOUBLE, grid1_rank, nc_grid1size_id, |
---|
335 | & nc_srcarray_id) |
---|
336 | call netcdf_error_handler(ncstat) |
---|
337 | |
---|
338 | !*** |
---|
339 | !*** define gradient arrays |
---|
340 | !*** |
---|
341 | |
---|
342 | ncstat = nf_def_var (nc_outfile_id, 'src_grad_lat', |
---|
343 | & NF_DOUBLE, grid1_rank, nc_grid1size_id, |
---|
344 | & nc_srcgradlat_id) |
---|
345 | call netcdf_error_handler(ncstat) |
---|
346 | |
---|
347 | ncstat = nf_def_var (nc_outfile_id, 'src_grad_lon', |
---|
348 | & NF_DOUBLE, grid1_rank, nc_grid1size_id, |
---|
349 | & nc_srcgradlon_id) |
---|
350 | call netcdf_error_handler(ncstat) |
---|
351 | |
---|
352 | !*** |
---|
353 | !*** define destination arrays |
---|
354 | !*** |
---|
355 | |
---|
356 | ncstat = nf_def_var (nc_outfile_id, 'dst_array1', |
---|
357 | & NF_DOUBLE, grid2_rank, nc_grid2size_id, |
---|
358 | & nc_dstarray1_id) |
---|
359 | call netcdf_error_handler(ncstat) |
---|
360 | |
---|
361 | ncstat = nf_def_var (nc_outfile_id, 'dst_array2', |
---|
362 | & NF_DOUBLE, grid2_rank, nc_grid2size_id, |
---|
363 | & nc_dstarray2_id) |
---|
364 | call netcdf_error_handler(ncstat) |
---|
365 | |
---|
366 | !*** |
---|
367 | !*** define error arrays |
---|
368 | !*** |
---|
369 | |
---|
370 | ncstat = nf_def_var (nc_outfile_id, 'dst_error1', |
---|
371 | & NF_DOUBLE, grid2_rank, nc_grid2size_id, |
---|
372 | & nc_dsterror1_id) |
---|
373 | call netcdf_error_handler(ncstat) |
---|
374 | |
---|
375 | ncstat = nf_def_var (nc_outfile_id, 'dst_error2', |
---|
376 | & NF_DOUBLE, grid2_rank, nc_grid2size_id, |
---|
377 | & nc_dsterror2_id) |
---|
378 | call netcdf_error_handler(ncstat) |
---|
379 | |
---|
380 | !*** |
---|
381 | !*** end definition stage |
---|
382 | !*** |
---|
383 | |
---|
384 | ncstat = nf_enddef(nc_outfile_id) |
---|
385 | call netcdf_error_handler(ncstat) |
---|
386 | |
---|
387 | !----------------------------------------------------------------------- |
---|
388 | ! |
---|
389 | ! write some grid info |
---|
390 | ! |
---|
391 | !----------------------------------------------------------------------- |
---|
392 | |
---|
393 | !*** |
---|
394 | !*** write grid center latitude array |
---|
395 | !*** |
---|
396 | |
---|
397 | ncstat = nf_put_var_double(nc_outfile_id, nc_srcgrdcntrlat_id, |
---|
398 | & grid1_center_lat) |
---|
399 | call netcdf_error_handler(ncstat) |
---|
400 | |
---|
401 | ncstat = nf_put_var_double(nc_outfile_id, nc_dstgrdcntrlat_id, |
---|
402 | & grid2_center_lat) |
---|
403 | call netcdf_error_handler(ncstat) |
---|
404 | |
---|
405 | !*** |
---|
406 | !*** write grid center longitude array |
---|
407 | !*** |
---|
408 | |
---|
409 | ncstat = nf_put_var_double(nc_outfile_id, nc_srcgrdcntrlon_id, |
---|
410 | & grid1_center_lon) |
---|
411 | call netcdf_error_handler(ncstat) |
---|
412 | |
---|
413 | ncstat = nf_put_var_double(nc_outfile_id, nc_dstgrdcntrlon_id, |
---|
414 | & grid2_center_lon) |
---|
415 | call netcdf_error_handler(ncstat) |
---|
416 | |
---|
417 | !*** |
---|
418 | !*** write grid mask |
---|
419 | !*** |
---|
420 | |
---|
421 | ncstat = nf_put_var_int(nc_outfile_id, nc_srcgrdimask_id, |
---|
422 | & grid1_imask) |
---|
423 | call netcdf_error_handler(ncstat) |
---|
424 | |
---|
425 | ncstat = nf_put_var_int(nc_outfile_id, nc_dstgrdimask_id, |
---|
426 | & grid2_imask) |
---|
427 | call netcdf_error_handler(ncstat) |
---|
428 | |
---|
429 | !*** |
---|
430 | !*** define grid area arrays |
---|
431 | !*** |
---|
432 | |
---|
433 | ncstat = nf_put_var_double(nc_outfile_id, nc_srcgrdarea_id, |
---|
434 | & grid1_area) |
---|
435 | call netcdf_error_handler(ncstat) |
---|
436 | |
---|
437 | ncstat = nf_put_var_double(nc_outfile_id, nc_dstgrdarea_id, |
---|
438 | & grid2_area) |
---|
439 | call netcdf_error_handler(ncstat) |
---|
440 | |
---|
441 | !*** |
---|
442 | !*** define grid fraction arrays |
---|
443 | !*** |
---|
444 | |
---|
445 | ncstat = nf_put_var_double(nc_outfile_id, nc_srcgrdfrac_id, |
---|
446 | & grid1_frac) |
---|
447 | call netcdf_error_handler(ncstat) |
---|
448 | |
---|
449 | ncstat = nf_put_var_double(nc_outfile_id, nc_dstgrdfrac_id, |
---|
450 | & grid2_frac) |
---|
451 | call netcdf_error_handler(ncstat) |
---|
452 | |
---|
453 | !----------------------------------------------------------------------- |
---|
454 | ! |
---|
455 | ! set up fields for test cases based on user choice |
---|
456 | ! |
---|
457 | !----------------------------------------------------------------------- |
---|
458 | |
---|
459 | select case (field_choice) |
---|
460 | case(1) !*** cosine hill at lon=pi and lat=0 |
---|
461 | |
---|
462 | length = 0.1*pi2 |
---|
463 | |
---|
464 | grid1_array = cos(grid1_center_lat)*cos(grid1_center_lon) |
---|
465 | grid2_array = cos(grid2_center_lat)*cos(grid2_center_lon) |
---|
466 | |
---|
467 | grid1_tmp = acos(-grid1_array)/length |
---|
468 | grid2_tmp = acos(-grid2_array)/length |
---|
469 | |
---|
470 | where (grid1_tmp <= one) |
---|
471 | grad1_lat = (pi/length)*sin(pi*grid1_tmp)* |
---|
472 | & sin(grid1_center_lat)*cos(grid1_center_lon)/ |
---|
473 | & sqrt(one-grid1_array**2) |
---|
474 | grad1_lon = (pi/length)*sin(pi*grid1_tmp)* |
---|
475 | & sin(grid1_center_lon)/ |
---|
476 | & sqrt(one-grid1_array**2) |
---|
477 | grid1_array = two + cos(pi*grid1_tmp) |
---|
478 | elsewhere |
---|
479 | grid1_array = one |
---|
480 | grad1_lat = zero |
---|
481 | grad1_lon = zero |
---|
482 | endwhere |
---|
483 | |
---|
484 | where (grid2_tmp <= one) |
---|
485 | grad2_lat = (pi/length)*sin(pi*grid2_tmp)* |
---|
486 | & sin(grid2_center_lat)*cos(grid2_center_lon)/ |
---|
487 | & sqrt(one-grid2_array**2) |
---|
488 | grad2_lon = (pi/length)*sin(pi*grid2_tmp)* |
---|
489 | & sin(grid2_center_lon)/ |
---|
490 | & sqrt(one-grid2_array**2) |
---|
491 | grid2_array = two + cos(pi*grid2_tmp) |
---|
492 | elsewhere |
---|
493 | grid2_array = one |
---|
494 | grad2_lat = zero |
---|
495 | grad2_lon = zero |
---|
496 | endwhere |
---|
497 | |
---|
498 | where (.not. grid1_mask) |
---|
499 | grid1_array = zero |
---|
500 | grad1_lat = zero |
---|
501 | grad1_lon = zero |
---|
502 | end where |
---|
503 | |
---|
504 | where (.not. grid2_mask) |
---|
505 | grid2_array = zero |
---|
506 | grad2_lat = zero |
---|
507 | grad2_lon = zero |
---|
508 | end where |
---|
509 | |
---|
510 | case(2) !*** pseudo-spherical harmonic l=2,m=2 |
---|
511 | |
---|
512 | where (grid1_mask) |
---|
513 | grid1_array = two + cos(grid1_center_lat)**2* |
---|
514 | & cos(two*grid1_center_lon) |
---|
515 | grad1_lat = -sin(two*grid1_center_lat)* |
---|
516 | & cos(two*grid1_center_lon) |
---|
517 | grad1_lon = -two*cos(grid1_center_lat)* |
---|
518 | & sin(two*grid1_center_lon) |
---|
519 | elsewhere |
---|
520 | grid1_array = zero |
---|
521 | grad1_lat = zero |
---|
522 | grad1_lon = zero |
---|
523 | end where |
---|
524 | |
---|
525 | where (grid2_mask) |
---|
526 | grid2_array = two + cos(grid2_center_lat)**2* |
---|
527 | & cos(two*grid2_center_lon) |
---|
528 | grad2_lat = -sin(two*grid2_center_lat)* |
---|
529 | & cos(two*grid2_center_lon) |
---|
530 | grad2_lon = -two*cos(grid2_center_lat)* |
---|
531 | & sin(two*grid2_center_lon) |
---|
532 | elsewhere |
---|
533 | grid2_array = zero |
---|
534 | grad2_lat = zero |
---|
535 | grad2_lon = zero |
---|
536 | end where |
---|
537 | |
---|
538 | case(3) !*** pseudo-spherical harmonic l=32, m=16 |
---|
539 | |
---|
540 | where (grid1_mask) |
---|
541 | grid1_array = two + sin(two*grid1_center_lat)**16* |
---|
542 | & cos(16.*grid1_center_lon) |
---|
543 | grad1_lat = 32.*sin(two*grid1_center_lat)**15* |
---|
544 | & cos(two*grid1_center_lat)* |
---|
545 | & cos(16.*grid1_center_lon) |
---|
546 | grad1_lon = -32.*sin(two*grid1_center_lat)**15* |
---|
547 | & sin(grid1_center_lat)* |
---|
548 | & sin(16.*grid1_center_lon) |
---|
549 | elsewhere |
---|
550 | grid1_array = zero |
---|
551 | grad1_lat = zero |
---|
552 | grad1_lon = zero |
---|
553 | end where |
---|
554 | |
---|
555 | where (grid2_mask) |
---|
556 | grid2_array = two + sin(two*grid2_center_lat)**16* |
---|
557 | & cos(16.*grid2_center_lon) |
---|
558 | grad2_lat = 32.*sin(two*grid2_center_lat)**15* |
---|
559 | & cos(two*grid2_center_lat)* |
---|
560 | & cos(16.*grid2_center_lon) |
---|
561 | grad2_lon = -32.*sin(two*grid2_center_lat)**15* |
---|
562 | & sin(grid2_center_lat)* |
---|
563 | & sin(16.*grid2_center_lon) |
---|
564 | elsewhere |
---|
565 | grid2_array = zero |
---|
566 | grad2_lat = zero |
---|
567 | grad2_lon = zero |
---|
568 | end where |
---|
569 | |
---|
570 | case default |
---|
571 | |
---|
572 | stop 'Bad choice for field to interpolate' |
---|
573 | |
---|
574 | end select |
---|
575 | |
---|
576 | !----------------------------------------------------------------------- |
---|
577 | ! |
---|
578 | ! test repeated first-order maps between grid1 and grid2 |
---|
579 | ! |
---|
580 | !----------------------------------------------------------------------- |
---|
581 | |
---|
582 | call remap(grid2_tmp, wts_map1, grid2_add_map1, grid1_add_map1, |
---|
583 | & grid1_array) |
---|
584 | do n=1,num_repeats |
---|
585 | call remap(grid1_tmp, wts_map2, grid1_add_map2, grid2_add_map2, |
---|
586 | & grid2_tmp) |
---|
587 | call remap(grid2_tmp, wts_map1, grid2_add_map1, grid1_add_map1, |
---|
588 | & grid1_tmp) |
---|
589 | end do |
---|
590 | |
---|
591 | where (grid2_frac > .999) |
---|
592 | grid2_err = (grid2_tmp - grid2_array)/grid2_array |
---|
593 | elsewhere |
---|
594 | grid2_err = zero |
---|
595 | end where |
---|
596 | |
---|
597 | print *,'First order mapping from grid1 to grid2:' |
---|
598 | print *,'----------------------------------------' |
---|
599 | print *,'Grid1 min,max: ',minval(grid1_array),maxval(grid1_array) |
---|
600 | print *,'Grid2 min,max: ',minval(grid2_tmp ),maxval(grid2_tmp ) |
---|
601 | print *,' Err2 min,max: ',minval(grid2_err),maxval(grid2_err) |
---|
602 | print *,' Err2 mean: ',sum(abs(grid2_err))/ |
---|
603 | & count(grid2_frac > .001) |
---|
604 | |
---|
605 | !*** |
---|
606 | !*** Conservation Test |
---|
607 | !*** |
---|
608 | |
---|
609 | print *,'Conservation:' |
---|
610 | print *,'Grid1 Integral = ',sum(grid1_array*grid1_area*grid1_frac) |
---|
611 | print *,'Grid2 Integral = ',sum(grid2_tmp *grid2_area*grid2_frac) |
---|
612 | |
---|
613 | !----------------------------------------------------------------------- |
---|
614 | ! |
---|
615 | ! write results to NetCDF file |
---|
616 | ! |
---|
617 | !----------------------------------------------------------------------- |
---|
618 | |
---|
619 | ncstat = nf_put_var_double(nc_outfile_id, nc_srcarray_id, |
---|
620 | & grid1_array) |
---|
621 | call netcdf_error_handler(ncstat) |
---|
622 | |
---|
623 | ncstat = nf_put_var_double(nc_outfile_id, nc_dstarray1_id, |
---|
624 | & grid2_tmp ) |
---|
625 | call netcdf_error_handler(ncstat) |
---|
626 | |
---|
627 | ncstat = nf_put_var_double(nc_outfile_id, nc_dsterror1_id, |
---|
628 | & grid2_err) |
---|
629 | call netcdf_error_handler(ncstat) |
---|
630 | |
---|
631 | !----------------------------------------------------------------------- |
---|
632 | ! |
---|
633 | ! test repeated second-order mapppings between grid1 and grid2 |
---|
634 | ! |
---|
635 | !----------------------------------------------------------------------- |
---|
636 | |
---|
637 | if (num_wts > 1) then |
---|
638 | |
---|
639 | call remap(grid2_tmp , wts_map1, grid2_add_map1, grid1_add_map1, |
---|
640 | & grid1_array, src_grad1=grad1_lat, |
---|
641 | & src_grad2=grad1_lon) |
---|
642 | |
---|
643 | do n=1,num_repeats |
---|
644 | call remap(grid1_tmp, wts_map2, grid1_add_map2, grid2_add_map2, |
---|
645 | & grid2_tmp, src_grad1=grad2_lat, |
---|
646 | & src_grad2=grad2_lon) |
---|
647 | |
---|
648 | call remap(grid2_tmp, wts_map1, grid2_add_map1, grid1_add_map1, |
---|
649 | & grid1_tmp, src_grad1=grad1_lat, |
---|
650 | & src_grad2=grad1_lon) |
---|
651 | |
---|
652 | end do |
---|
653 | |
---|
654 | where (grid2_frac > .999) |
---|
655 | grid2_err = (grid2_tmp - grid2_array)/grid2_array |
---|
656 | elsewhere |
---|
657 | grid2_err = zero |
---|
658 | end where |
---|
659 | |
---|
660 | print *,'Second order mapping from grid1 to grid2:' |
---|
661 | print *,'-----------------------------------------' |
---|
662 | print *,'Grid1 min,max: ',minval(grid1_array),maxval(grid1_array) |
---|
663 | print *,'Grid2 min,max: ',minval(grid2_tmp ),maxval(grid2_tmp ) |
---|
664 | print *,' Err2 min,max: ',minval(grid2_err),maxval(grid2_err) |
---|
665 | print *,' Err2 mean: ',sum(abs(grid2_err))/ |
---|
666 | & count(grid2_frac > .001) |
---|
667 | |
---|
668 | !*** |
---|
669 | !*** Conservation Test |
---|
670 | !*** |
---|
671 | |
---|
672 | print *,'Conservation:' |
---|
673 | print *,'Grid1 Integral = ',sum(grid1_array*grid1_area*grid1_frac) |
---|
674 | print *,'Grid2 Integral = ',sum(grid2_tmp *grid2_area*grid2_frac) |
---|
675 | |
---|
676 | !----------------------------------------------------------------------- |
---|
677 | ! |
---|
678 | ! write results to NetCDF file |
---|
679 | ! |
---|
680 | !----------------------------------------------------------------------- |
---|
681 | |
---|
682 | ncstat = nf_put_var_double(nc_outfile_id, nc_srcgradlon_id, |
---|
683 | & grad1_lon) |
---|
684 | call netcdf_error_handler(ncstat) |
---|
685 | |
---|
686 | ncstat = nf_put_var_double(nc_outfile_id, nc_dstarray2_id, |
---|
687 | & grid2_tmp ) |
---|
688 | call netcdf_error_handler(ncstat) |
---|
689 | |
---|
690 | ncstat = nf_put_var_double(nc_outfile_id, nc_dsterror2_id, |
---|
691 | & grid2_err) |
---|
692 | call netcdf_error_handler(ncstat) |
---|
693 | |
---|
694 | endif |
---|
695 | |
---|
696 | !----------------------------------------------------------------------- |
---|
697 | ! |
---|
698 | ! close netCDF file |
---|
699 | ! |
---|
700 | !----------------------------------------------------------------------- |
---|
701 | |
---|
702 | ncstat = nf_close(nc_outfile_id) |
---|
703 | call netcdf_error_handler(ncstat) |
---|
704 | |
---|
705 | !----------------------------------------------------------------------- |
---|
706 | |
---|
707 | end program remap_test_repeat |
---|
708 | |
---|
709 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|