source: CPL/oasis3/trunk/src/lib/scrip/src/remap_write.f @ 1677

Last change on this file since 1677 was 1677, checked in by aclsce, 12 years ago

Imported oasis3 (tag ipslcm5a) from cvs server to svn server (igcmg project).

File size: 50.0 KB
Line 
1C****
2C                    ************************
3C                    *     OASIS MODULE     *
4C                    *     ------------     *
5C                    ************************
6C**** 
7C***********************************************************************
8C     This module belongs to the SCRIP library. It is modified to run
9C     within OASIS. 
10C     Modifications:
11C       - Added CASE for GAUSWGT
12C
13C     Modified by            D. Declat,  CERFACS              27.06.2002
14C***********************************************************************
15!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16!
17!     This module contains routines for writing the remapping data to
18!     a file.  Before writing the data for each mapping, the links are
19!     sorted by destination grid address.
20!
21!-----------------------------------------------------------------------
22!
23!     CVS:$Id: remap_write.f,v 1.1.1.1 2005/03/23 16:01:12 adm Exp $
24!
25!     Copyright (c) 1997, 1998 the Regents of the University of
26!       California.
27!
28!     This software and ancillary information (herein called software)
29!     called SCRIP is made available under the terms described here. 
30!     The software has been approved for release with associated
31!     LA-CC Number 98-45.
32!
33!     Unless otherwise indicated, this software has been authored
34!     by an employee or employees of the University of California,
35!     operator of the Los Alamos National Laboratory under Contract
36!     No. W-7405-ENG-36 with the U.S. Department of Energy.  The U.S.
37!     Government has rights to use, reproduce, and distribute this
38!     software.  The public may copy and use this software without
39!     charge, provided that this Notice and any statement of authorship
40!     are reproduced on all copies.  Neither the Government nor the
41!     University makes any warranty, express or implied, or assumes
42!     any liability or responsibility for the use of this software.
43!
44!     If software is modified to produce derivative works, such modified
45!     software should be clearly marked, so as not to confuse it with
46!     the version available from Los Alamos National Laboratory.
47!
48!***********************************************************************
49
50      module remap_write
51
52!-----------------------------------------------------------------------
53
54      use kinds_mod     ! defines common data types
55      use constants     ! defines common scalar constants
56      use grids         ! module containing grid information
57      use remap_vars    ! module containing remap information
58      use netcdf_mod    ! module with netCDF stuff
59
60      implicit none
61
62!-----------------------------------------------------------------------
63!
64!     module variables
65!
66!-----------------------------------------------------------------------
67
68      character(char_len), private :: 
69     &   map_method       ! character string for map_type
70     &,  normalize_opt    ! character string for normalization option
71     &,  history          ! character string for history information
72     &,  convention       ! character string for output convention
73
74      character(8), private :: 
75     &   cdate            ! character date string
76
77      integer (kind=int_kind), dimension(:), allocatable, private ::
78     &   src_mask_int     ! integer masks to determine
79     &,  dst_mask_int     ! cells that participate in map
80
81!-----------------------------------------------------------------------
82!
83!     various netCDF identifiers used by output routines
84!
85!-----------------------------------------------------------------------
86
87      integer (kind=int_kind), private ::
88     &   ncstat               ! error flag for netCDF calls
89     &,  nc_file_id           ! id for netCDF file
90     &,  nc_srcgrdsize_id     ! id for source grid size
91     &,  nc_dstgrdsize_id     ! id for destination grid size
92     &,  nc_srcgrdcorn_id     ! id for number of source grid corners
93     &,  nc_dstgrdcorn_id     ! id for number of dest grid corners
94     &,  nc_srcgrdrank_id     ! id for source grid rank
95     &,  nc_dstgrdrank_id     ! id for dest grid rank
96     &,  nc_numlinks_id       ! id for number of links in mapping
97     &,  nc_numwgts_id        ! id for number of weights for mapping
98     &,  nc_srcgrddims_id     ! id for source grid dimensions
99     &,  nc_dstgrddims_id     ! id for dest grid dimensions
100     &,  nc_srcgrdcntrlat_id  ! id for source grid center latitude
101     &,  nc_dstgrdcntrlat_id  ! id for dest grid center latitude
102     &,  nc_srcgrdcntrlon_id  ! id for source grid center longitude
103     &,  nc_dstgrdcntrlon_id  ! id for dest grid center longitude
104     &,  nc_srcgrdimask_id    ! id for source grid mask
105     &,  nc_dstgrdimask_id    ! id for dest grid mask
106     &,  nc_srcgrdcrnrlat_id  ! id for latitude of source grid corners
107     &,  nc_srcgrdcrnrlon_id  ! id for longitude of source grid corners
108     &,  nc_dstgrdcrnrlat_id  ! id for latitude of dest grid corners
109     &,  nc_dstgrdcrnrlon_id  ! id for longitude of dest grid corners
110     &,  nc_srcgrdarea_id     ! id for area of source grid cells
111     &,  nc_dstgrdarea_id     ! id for area of dest grid cells
112     &,  nc_srcgrdfrac_id     ! id for area fraction on source grid
113     &,  nc_dstgrdfrac_id     ! id for area fraction on dest grid
114     &,  nc_srcadd_id         ! id for map source address
115     &,  nc_dstadd_id         ! id for map destination address
116     &,  nc_rmpmatrix_id      ! id for remapping matrix
117
118      integer (kind=int_kind), dimension(2), private ::
119     &   nc_dims2_id  ! netCDF ids for 2d array dims
120
121!***********************************************************************
122
123      contains
124
125!***********************************************************************
126
127      subroutine write_remap(map1_name, map2_name, 
128     &                       interp_file1, interp_file2, output_opt)
129
130!-----------------------------------------------------------------------
131!
132!     calls correct output routine based on output format choice
133!
134!-----------------------------------------------------------------------
135
136!-----------------------------------------------------------------------
137!
138!     input variables
139!
140!-----------------------------------------------------------------------
141
142      character(char_len), intent(in) ::
143     &            map1_name,    ! name for mapping grid1 to grid2
144     &            map2_name,    ! name for mapping grid2 to grid1
145     &            interp_file1, ! filename for map1 remap data
146     &            interp_file2, ! filename for map2 remap data
147     &            output_opt    ! option for output conventions
148
149!-----------------------------------------------------------------------
150!
151!     local variables
152!
153!-----------------------------------------------------------------------
154!-----------------------------------------------------------------------
155!
156!     define some common variables to be used in all routines
157!
158!-----------------------------------------------------------------------
159      select case(norm_opt)
160      case (norm_opt_none)
161        normalize_opt = 'none'
162      case (norm_opt_frcarea)
163        normalize_opt = 'fracarea'
164      case (norm_opt_dstarea)
165        normalize_opt = 'destarea'
166      case (norm_opt_nonorm)
167        normalize_opt = 'no norm'
168      end select
169      select case(map_type)
170      case(map_type_conserv)
171        map_method = 'Conservative remapping'
172      case(map_type_bilinear)
173        map_method = 'Bilinear remapping'
174      case(map_type_distwgt)
175        map_method = 'Distance weighted avg of nearest neighbors'
176      case(map_type_gauswgt)
177        map_method = 'Gaussian weighted avg of nearest neighbors'
178      case(map_type_bicubic)
179        map_method = 'Bicubic remapping'
180      case default
181        stop 'Invalid Map Type'
182      end select
183
184      call date_and_time(date=cdate)
185      write (history,1000) cdate(5:6),cdate(7:8),cdate(1:4)
186 1000 format('Created: ',a2,'-',a2,'-',a4)
187
188!-----------------------------------------------------------------------
189!
190!     call appropriate output routine
191!
192!-----------------------------------------------------------------------
193
194      select case(output_opt)
195      case ('scrip')
196        call write_remap_scrip(map1_name, interp_file1, 1)
197      case ('ncar-csm')
198        call write_remap_csm  (map1_name, interp_file1, 1)
199      case default
200        stop 'unknown output file convention'
201      end select
202
203!-----------------------------------------------------------------------
204!
205!     call appropriate output routine for second mapping if required
206!
207!-----------------------------------------------------------------------
208
209      if (num_maps > 1) then
210        select case(output_opt)
211        case ('scrip')
212          call write_remap_scrip(map2_name, interp_file2, 2)
213        case ('ncar-csm')
214          call write_remap_csm  (map2_name, interp_file2, 2)
215        case default
216          stop 'unknown output file convention'
217        end select
218      endif
219
220!-----------------------------------------------------------------------
221
222      end subroutine write_remap
223
224!***********************************************************************
225
226      subroutine write_remap_scrip(map_name, interp_file, direction)
227
228!-----------------------------------------------------------------------
229!
230!     writes remap data to a netCDF file using SCRIP conventions
231!
232!-----------------------------------------------------------------------
233
234!-----------------------------------------------------------------------
235!
236!     input variables
237!
238!-----------------------------------------------------------------------
239
240      character(char_len), intent(in) ::
241     &            map_name     ! name for mapping
242     &,           interp_file  ! filename for remap data
243
244      integer (kind=int_kind), intent(in) ::
245     &  direction              ! direction of map (1=grid1 to grid2
246                               !                   2=grid2 to grid1)
247
248!-----------------------------------------------------------------------
249!
250!     local variables
251!
252!-----------------------------------------------------------------------
253
254      character(char_len) ::
255     &  grid1_ctmp        ! character temp for grid1 names
256     &, grid2_ctmp        ! character temp for grid2 names
257
258      integer (kind=int_kind) ::
259     &  itmp1             ! integer temp
260     &, itmp2             ! integer temp
261     &, itmp3             ! integer temp
262     &, itmp4             ! integer temp
263
264!-----------------------------------------------------------------------
265!
266!     create netCDF file for mapping and define some global attributes
267!
268!-----------------------------------------------------------------------
269      ncstat = nf_create (interp_file, NF_CLOBBER, nc_file_id)
270      call netcdf_error_handler(ncstat)
271      !***
272      !*** map name
273      !***
274      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, 'title',
275     &                          len_trim(map_name), map_name)
276      call netcdf_error_handler(ncstat)
277
278      !***
279      !*** normalization option
280      !***
281      ncstat = nf_put_att_text(nc_file_id, NF_GLOBAL, 'normalization',
282     &                         len_trim(normalize_opt), normalize_opt)
283      call netcdf_error_handler(ncstat)
284
285      !***
286      !*** map method
287      !***
288      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, 'map_method',
289     &                          len_trim(map_method), map_method)
290      call netcdf_error_handler(ncstat)
291
292      !***
293      !*** history
294      !***
295      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, 'history',
296     &                          len_trim(history), history)
297      call netcdf_error_handler(ncstat)
298
299      !***
300      !*** file convention
301      !***
302      convention = 'SCRIP'
303      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, 'conventions',
304     &                          len_trim(convention), convention)
305      call netcdf_error_handler(ncstat)
306
307      !***
308      !*** source and destination grid names
309      !***
310
311      if (direction == 1) then
312        grid1_ctmp = 'source_grid'
313        grid2_ctmp = 'dest_grid'
314      else
315        grid1_ctmp = 'dest_grid'
316        grid2_ctmp = 'source_grid'
317      endif
318
319      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, trim(grid1_ctmp),
320     &                          len_trim(grid1_name), grid1_name)
321      call netcdf_error_handler(ncstat)
322      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, trim(grid2_ctmp),
323     &                          len_trim(grid2_name), grid2_name)
324      call netcdf_error_handler(ncstat)
325
326!-----------------------------------------------------------------------
327!
328!     prepare netCDF dimension info
329!
330!-----------------------------------------------------------------------
331
332      !***
333      !*** define grid size dimensions
334      !***
335
336      if (direction == 1) then
337        itmp1 = grid1_size
338        itmp2 = grid2_size
339      else
340        itmp1 = grid2_size
341        itmp2 = grid1_size
342      endif
343
344      ncstat = nf_def_dim (nc_file_id, 'src_grid_size', itmp1, 
345     &                     nc_srcgrdsize_id)
346      call netcdf_error_handler(ncstat)
347
348      ncstat = nf_def_dim (nc_file_id, 'dst_grid_size', itmp2, 
349     &                     nc_dstgrdsize_id)
350      call netcdf_error_handler(ncstat)
351
352      !***
353      !*** define grid corner dimension
354      !***
355
356      if (direction == 1) then
357        itmp1 = grid1_corners
358        itmp2 = grid2_corners
359      else
360        itmp1 = grid2_corners
361        itmp2 = grid1_corners
362      endif
363
364      ncstat = nf_def_dim (nc_file_id, 'src_grid_corners', 
365     &                     itmp1, nc_srcgrdcorn_id)
366      call netcdf_error_handler(ncstat)
367
368      ncstat = nf_def_dim (nc_file_id, 'dst_grid_corners', 
369     &                     itmp2, nc_dstgrdcorn_id)
370      call netcdf_error_handler(ncstat)
371      !***
372      !*** define grid rank dimension
373      !***
374
375      if (direction == 1) then
376        itmp1 = grid1_rank
377        itmp2 = grid2_rank
378      else
379        itmp1 = grid2_rank
380        itmp2 = grid1_rank
381      endif
382
383      ncstat = nf_def_dim (nc_file_id, 'src_grid_rank', 
384     &                     itmp1, nc_srcgrdrank_id)
385      call netcdf_error_handler(ncstat)
386
387      ncstat = nf_def_dim (nc_file_id, 'dst_grid_rank', 
388     &                     itmp2, nc_dstgrdrank_id)
389      call netcdf_error_handler(ncstat)
390      !***
391      !*** define map size dimensions
392      !***
393
394      if (direction == 1) then
395        itmp1 = num_links_map1
396      else
397        itmp1 = num_links_map2
398      endif
399
400      ncstat = nf_def_dim (nc_file_id, 'num_links', 
401     &                     itmp1, nc_numlinks_id)
402      call netcdf_error_handler(ncstat)
403      ncstat = nf_def_dim (nc_file_id, 'num_wgts', 
404     &                     num_wts, nc_numwgts_id)
405      call netcdf_error_handler(ncstat)
406      !***
407      !*** define grid dimensions
408      !***
409
410      ncstat = nf_def_var (nc_file_id, 'src_grid_dims', NF_INT,
411     &                     1, nc_srcgrdrank_id, nc_srcgrddims_id)
412      call netcdf_error_handler(ncstat)
413      ncstat = nf_def_var (nc_file_id, 'dst_grid_dims', NF_INT,
414     &                     1, nc_dstgrdrank_id, nc_dstgrddims_id)
415      call netcdf_error_handler(ncstat)
416!-----------------------------------------------------------------------
417!
418!     define all arrays for netCDF descriptors
419!
420!-----------------------------------------------------------------------
421
422      !***
423      !*** define grid center latitude array
424      !***
425
426      ncstat = nf_def_var (nc_file_id, 'src_grid_center_lat', 
427     &                     NF_DOUBLE, 1, nc_srcgrdsize_id, 
428     &                     nc_srcgrdcntrlat_id)
429      call netcdf_error_handler(ncstat)
430      ncstat = nf_def_var (nc_file_id, 'dst_grid_center_lat', 
431     &                     NF_DOUBLE, 1, nc_dstgrdsize_id, 
432     &                     nc_dstgrdcntrlat_id)
433      call netcdf_error_handler(ncstat)
434      !***
435      !*** define grid center longitude array
436      !***
437
438      ncstat = nf_def_var (nc_file_id, 'src_grid_center_lon', 
439     &                     NF_DOUBLE, 1, nc_srcgrdsize_id, 
440     &                     nc_srcgrdcntrlon_id)
441      call netcdf_error_handler(ncstat)
442      ncstat = nf_def_var (nc_file_id, 'dst_grid_center_lon', 
443     &                     NF_DOUBLE, 1, nc_dstgrdsize_id, 
444     &                     nc_dstgrdcntrlon_id)
445      call netcdf_error_handler(ncstat)
446      !***
447      !*** define grid corner lat/lon arrays
448      !***
449
450      nc_dims2_id(1) = nc_srcgrdcorn_id
451      nc_dims2_id(2) = nc_srcgrdsize_id
452
453      ncstat = nf_def_var (nc_file_id, 'src_grid_corner_lat', 
454     &                     NF_DOUBLE, 2, nc_dims2_id, 
455     &                     nc_srcgrdcrnrlat_id)
456      call netcdf_error_handler(ncstat)
457      ncstat = nf_def_var (nc_file_id, 'src_grid_corner_lon', 
458     &                     NF_DOUBLE, 2, nc_dims2_id, 
459     &                     nc_srcgrdcrnrlon_id)
460      call netcdf_error_handler(ncstat)
461
462      nc_dims2_id(1) = nc_dstgrdcorn_id
463      nc_dims2_id(2) = nc_dstgrdsize_id
464
465      ncstat = nf_def_var (nc_file_id, 'dst_grid_corner_lat', 
466     &                     NF_DOUBLE, 2, nc_dims2_id, 
467     &                     nc_dstgrdcrnrlat_id)
468      call netcdf_error_handler(ncstat)
469      ncstat = nf_def_var (nc_file_id, 'dst_grid_corner_lon', 
470     &                     NF_DOUBLE, 2, nc_dims2_id, 
471     &                     nc_dstgrdcrnrlon_id)
472      call netcdf_error_handler(ncstat)
473
474      !***
475      !*** define units for all coordinate arrays
476      !***
477
478      if (direction == 1) then
479        grid1_ctmp = grid1_units
480        grid2_ctmp = grid2_units
481      else
482        grid1_ctmp = grid2_units
483        grid2_ctmp = grid1_units
484      endif
485
486      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdcntrlat_id, 
487     &                          'units', 7, grid1_ctmp)
488      call netcdf_error_handler(ncstat)
489      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdcntrlat_id, 
490     &                          'units', 7, grid2_ctmp)
491      call netcdf_error_handler(ncstat)
492
493      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdcntrlon_id, 
494     &                          'units', 7, grid1_ctmp)
495      call netcdf_error_handler(ncstat)
496
497      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdcntrlon_id, 
498     &                          'units', 7, grid2_ctmp)
499      call netcdf_error_handler(ncstat)
500
501      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdcrnrlat_id, 
502     &                          'units', 7, grid1_ctmp)
503      call netcdf_error_handler(ncstat)
504      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdcrnrlon_id, 
505     &                          'units', 7, grid1_ctmp)
506      call netcdf_error_handler(ncstat)
507
508      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdcrnrlat_id, 
509     &                          'units', 7, grid2_ctmp)
510      call netcdf_error_handler(ncstat)
511
512      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdcrnrlon_id, 
513     &                          'units', 7, grid2_ctmp)
514      call netcdf_error_handler(ncstat)
515
516      !***
517      !*** define grid mask
518      !***
519
520      ncstat = nf_def_var (nc_file_id, 'src_grid_imask', NF_INT,
521     &                     1, nc_srcgrdsize_id, nc_srcgrdimask_id)
522      call netcdf_error_handler(ncstat)
523
524      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdimask_id, 
525     &                          'units', 8, 'unitless')
526      call netcdf_error_handler(ncstat)
527      ncstat = nf_def_var (nc_file_id, 'dst_grid_imask', NF_INT,
528     &                     1, nc_dstgrdsize_id, nc_dstgrdimask_id)
529      call netcdf_error_handler(ncstat)
530
531      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdimask_id, 
532     &                          'units', 8, 'unitless')
533      call netcdf_error_handler(ncstat)
534      !***
535      !*** define grid area arrays
536      !***
537
538      ncstat = nf_def_var (nc_file_id, 'src_grid_area', 
539     &                     NF_DOUBLE, 1, nc_srcgrdsize_id, 
540     &                     nc_srcgrdarea_id)
541      call netcdf_error_handler(ncstat)
542      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdarea_id, 
543     &                          'units', 14, 'square radians')
544      call netcdf_error_handler(ncstat)
545
546      ncstat = nf_def_var (nc_file_id, 'dst_grid_area', 
547     &                     NF_DOUBLE, 1, nc_dstgrdsize_id, 
548     &                     nc_dstgrdarea_id)
549      call netcdf_error_handler(ncstat)
550
551      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdarea_id, 
552     &                          'units', 14, 'square radians')
553      call netcdf_error_handler(ncstat)
554
555      !***
556      !*** define grid fraction arrays
557      !***
558
559      ncstat = nf_def_var (nc_file_id, 'src_grid_frac', 
560     &                     NF_DOUBLE, 1, nc_srcgrdsize_id, 
561     &                     nc_srcgrdfrac_id)
562      call netcdf_error_handler(ncstat)
563      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdfrac_id, 
564     &                          'units', 8, 'unitless')
565      call netcdf_error_handler(ncstat)
566
567      ncstat = nf_def_var (nc_file_id, 'dst_grid_frac', 
568     &                     NF_DOUBLE, 1, nc_dstgrdsize_id, 
569     &                     nc_dstgrdfrac_id)
570      call netcdf_error_handler(ncstat)
571
572      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdfrac_id, 
573     &                          'units', 8, 'unitless')
574      call netcdf_error_handler(ncstat)
575      !***
576      !*** define mapping arrays
577      !***
578
579      ncstat = nf_def_var (nc_file_id, 'src_address', 
580     &                     NF_INT, 1, nc_numlinks_id, 
581     &                     nc_srcadd_id)
582      call netcdf_error_handler(ncstat)
583
584      ncstat = nf_def_var (nc_file_id, 'dst_address', 
585     &                     NF_INT, 1, nc_numlinks_id, 
586     &                     nc_dstadd_id)
587      call netcdf_error_handler(ncstat)
588
589      nc_dims2_id(1) = nc_numwgts_id
590      nc_dims2_id(2) = nc_numlinks_id
591
592      ncstat = nf_def_var (nc_file_id, 'remap_matrix', 
593     &                     NF_DOUBLE, 2, nc_dims2_id, 
594     &                     nc_rmpmatrix_id)
595      call netcdf_error_handler(ncstat)
596
597      !***
598      !*** end definition stage
599      !***
600
601      ncstat = nf_enddef(nc_file_id)
602      call netcdf_error_handler(ncstat)
603!-----------------------------------------------------------------------
604!
605!     compute integer masks
606!
607!-----------------------------------------------------------------------
608
609      if (direction == 1) then
610        allocate (src_mask_int(grid1_size),
611     &            dst_mask_int(grid2_size))
612
613        where (grid2_mask)
614          dst_mask_int = 1
615        elsewhere
616          dst_mask_int = 0
617        endwhere
618
619        where (grid1_mask)
620          src_mask_int = 1
621        elsewhere
622          src_mask_int = 0
623        endwhere
624      else
625        allocate (src_mask_int(grid2_size),
626     &            dst_mask_int(grid1_size))
627
628        where (grid1_mask)
629          dst_mask_int = 1
630        elsewhere
631          dst_mask_int = 0
632        endwhere
633
634        where (grid2_mask)
635          src_mask_int = 1
636        elsewhere
637          src_mask_int = 0
638        endwhere
639      endif
640
641!-----------------------------------------------------------------------
642!
643!     change units of lat/lon coordinates if input units different
644!     from radians
645!
646!-----------------------------------------------------------------------
647
648      if (grid1_units(1:7) == 'degrees' .and. direction == 1) then
649        grid1_center_lat = grid1_center_lat/deg2rad
650        grid1_center_lon = grid1_center_lon/deg2rad
651        grid1_corner_lat = grid1_corner_lat/deg2rad
652        grid1_corner_lon = grid1_corner_lon/deg2rad
653      endif
654
655      if (grid2_units(1:7) == 'degrees' .and. direction == 1) then
656        grid2_center_lat = grid2_center_lat/deg2rad
657        grid2_center_lon = grid2_center_lon/deg2rad
658        grid2_corner_lat = grid2_corner_lat/deg2rad
659        grid2_corner_lon = grid2_corner_lon/deg2rad
660      endif
661
662!-----------------------------------------------------------------------
663!
664!     write mapping data
665!
666!-----------------------------------------------------------------------
667
668      if (direction == 1) then
669        itmp1 = nc_srcgrddims_id
670        itmp2 = nc_dstgrddims_id
671      else
672        itmp2 = nc_srcgrddims_id
673        itmp1 = nc_dstgrddims_id
674      endif
675
676      ncstat = nf_put_var_int(nc_file_id, itmp1, grid1_dims)
677      call netcdf_error_handler(ncstat)
678
679      ncstat = nf_put_var_int(nc_file_id, itmp2, grid2_dims)
680      call netcdf_error_handler(ncstat)
681
682      ncstat = nf_put_var_int(nc_file_id, nc_srcgrdimask_id, 
683     &                        src_mask_int)
684      call netcdf_error_handler(ncstat)
685      ncstat = nf_put_var_int(nc_file_id, nc_dstgrdimask_id,
686     &                        dst_mask_int)
687      call netcdf_error_handler(ncstat)
688
689      deallocate(src_mask_int, dst_mask_int)
690
691      if (direction == 1) then
692        itmp1 = nc_srcgrdcntrlat_id
693        itmp2 = nc_srcgrdcntrlon_id
694        itmp3 = nc_srcgrdcrnrlat_id
695        itmp4 = nc_srcgrdcrnrlon_id
696      else
697        itmp1 = nc_dstgrdcntrlat_id
698        itmp2 = nc_dstgrdcntrlon_id
699        itmp3 = nc_dstgrdcrnrlat_id
700        itmp4 = nc_dstgrdcrnrlon_id
701      endif
702
703      ncstat = nf_put_var_double(nc_file_id, itmp1, grid1_center_lat)
704      call netcdf_error_handler(ncstat)
705
706      ncstat = nf_put_var_double(nc_file_id, itmp2, grid1_center_lon)
707      call netcdf_error_handler(ncstat)
708      if (.not. luse_grid_centers) then
709          ncstat = nf_put_var_double(nc_file_id, itmp3, 
710     $        grid1_corner_lat)
711          call netcdf_error_handler(ncstat)
712
713          ncstat = nf_put_var_double(nc_file_id, itmp4, 
714     $        grid1_corner_lon)
715          call netcdf_error_handler(ncstat)
716      endif
717
718      if (direction == 1) then
719        itmp1 = nc_dstgrdcntrlat_id
720        itmp2 = nc_dstgrdcntrlon_id
721        itmp3 = nc_dstgrdcrnrlat_id
722        itmp4 = nc_dstgrdcrnrlon_id
723      else
724        itmp1 = nc_srcgrdcntrlat_id
725        itmp2 = nc_srcgrdcntrlon_id
726        itmp3 = nc_srcgrdcrnrlat_id
727        itmp4 = nc_srcgrdcrnrlon_id
728      endif
729
730      ncstat = nf_put_var_double(nc_file_id, itmp1, grid2_center_lat)
731      call netcdf_error_handler(ncstat)
732
733      ncstat = nf_put_var_double(nc_file_id, itmp2, grid2_center_lon)
734      call netcdf_error_handler(ncstat)
735
736      if (.not. luse_grid_centers) then
737          ncstat = nf_put_var_double(nc_file_id, itmp3, 
738     $        grid2_corner_lat)
739          call netcdf_error_handler(ncstat)
740
741          ncstat = nf_put_var_double(nc_file_id, itmp4, 
742     $        grid2_corner_lon)
743          call netcdf_error_handler(ncstat)
744      endif
745      if (direction == 1) then
746        itmp1 = nc_srcgrdarea_id
747        itmp2 = nc_srcgrdfrac_id
748        itmp3 = nc_dstgrdarea_id
749        itmp4 = nc_dstgrdfrac_id
750      else
751        itmp1 = nc_dstgrdarea_id
752        itmp2 = nc_dstgrdfrac_id
753        itmp3 = nc_srcgrdarea_id
754        itmp4 = nc_srcgrdfrac_id
755      endif
756
757      if (luse_grid1_area) then
758        ncstat = nf_put_var_double(nc_file_id, itmp1, grid1_area_in)
759      else
760        ncstat = nf_put_var_double(nc_file_id, itmp1, grid1_area)
761      endif
762      call netcdf_error_handler(ncstat)
763
764      ncstat = nf_put_var_double(nc_file_id, itmp2, grid1_frac)
765      call netcdf_error_handler(ncstat)
766
767      if (luse_grid2_area) then
768        ncstat = nf_put_var_double(nc_file_id, itmp3, grid2_area_in)
769      else
770        ncstat = nf_put_var_double(nc_file_id, itmp3, grid2_area)
771      endif
772      call netcdf_error_handler(ncstat)
773
774      ncstat = nf_put_var_double(nc_file_id, itmp4, grid2_frac)
775      call netcdf_error_handler(ncstat)
776      if (direction == 1) then
777        ncstat = nf_put_var_int(nc_file_id, nc_srcadd_id, 
778     &                          grid1_add_map1)
779        call netcdf_error_handler(ncstat)
780
781        ncstat = nf_put_var_int(nc_file_id, nc_dstadd_id, 
782     &                          grid2_add_map1)
783        call netcdf_error_handler(ncstat)
784
785        ncstat = nf_put_var_double(nc_file_id, nc_rmpmatrix_id, 
786     &                             wts_map1)
787        call netcdf_error_handler(ncstat)
788      else
789        ncstat = nf_put_var_int(nc_file_id, nc_srcadd_id, 
790     &                          grid2_add_map2)
791        call netcdf_error_handler(ncstat)
792
793        ncstat = nf_put_var_int(nc_file_id, nc_dstadd_id, 
794     &                          grid1_add_map2)
795        call netcdf_error_handler(ncstat)
796
797        ncstat = nf_put_var_double(nc_file_id, nc_rmpmatrix_id, 
798     &                             wts_map2)
799        call netcdf_error_handler(ncstat)
800      endif
801
802      ncstat = nf_close(nc_file_id)
803      call netcdf_error_handler(ncstat)
804!-----------------------------------------------------------------------
805
806      end subroutine write_remap_scrip
807
808!***********************************************************************
809
810      subroutine write_remap_csm(map_name, interp_file, direction)
811
812!-----------------------------------------------------------------------
813!
814!     writes remap data to a netCDF file using NCAR-CSM conventions
815!
816!-----------------------------------------------------------------------
817
818!-----------------------------------------------------------------------
819!
820!     input variables
821!
822!-----------------------------------------------------------------------
823
824      character(char_len), intent(in) ::
825     &            map_name     ! name for mapping
826     &,           interp_file  ! filename for remap data
827
828      integer (kind=int_kind), intent(in) ::
829     &  direction              ! direction of map (1=grid1 to grid2
830                               !                   2=grid2 to grid1)
831
832!-----------------------------------------------------------------------
833!
834!     local variables
835!
836!-----------------------------------------------------------------------
837
838      character(char_len) ::
839     &  grid1_ctmp        ! character temp for grid1 names
840     &, grid2_ctmp        ! character temp for grid2 names
841
842      integer (kind=int_kind) ::
843     &  itmp1             ! integer temp
844     &, itmp2             ! integer temp
845     &, itmp3             ! integer temp
846     &, itmp4             ! integer temp
847     &, nc_numwgts1_id    ! extra netCDF id for additional weights
848     &, nc_src_isize_id   ! extra netCDF id for ni_a
849     &, nc_src_jsize_id   ! extra netCDF id for nj_a
850     &, nc_dst_isize_id   ! extra netCDF id for ni_b
851     &, nc_dst_jsize_id   ! extra netCDF id for nj_b
852     &, nc_rmpmatrix2_id  ! extra netCDF id for high-order remap matrix
853
854      real (kind=dbl_kind), dimension(:),allocatable ::
855     &  wts1              ! CSM wants single array for 1st-order wts
856
857      real (kind=dbl_kind), dimension(:,:),allocatable ::
858     &  wts2              ! write remaining weights in different array
859
860!-----------------------------------------------------------------------
861!
862!     create netCDF file for mapping and define some global attributes
863!
864!-----------------------------------------------------------------------
865
866      ncstat = nf_create (interp_file, NF_CLOBBER, nc_file_id)
867      call netcdf_error_handler(ncstat)
868
869      !***
870      !*** map name
871      !***
872      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, 'title',
873     &                          len_trim(map_name), map_name)
874      call netcdf_error_handler(ncstat)
875
876      !***
877      !*** normalization option
878      !***
879      ncstat = nf_put_att_text(nc_file_id, NF_GLOBAL, 'normalization',
880     &                         len_trim(normalize_opt), normalize_opt)
881      call netcdf_error_handler(ncstat)
882
883      !***
884      !*** map method
885      !***
886      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, 'map_method',
887     &                          len_trim(map_method), map_method)
888      call netcdf_error_handler(ncstat)
889
890      !***
891      !*** history
892      !***
893      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, 'history',
894     &                          len_trim(history), history)
895      call netcdf_error_handler(ncstat)
896
897      !***
898      !*** file convention
899      !***
900      convention = 'NCAR-CSM'
901      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, 'conventions',
902     &                          len_trim(convention), convention)
903      call netcdf_error_handler(ncstat)
904
905      !***
906      !*** source and destination grid names
907      !***
908
909      if (direction == 1) then
910        grid1_ctmp = 'domain_a'
911        grid2_ctmp = 'domain_b'
912      else
913        grid1_ctmp = 'domain_b'
914        grid2_ctmp = 'domain_a'
915      endif
916
917      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, trim(grid1_ctmp),
918     &                          len_trim(grid1_name), grid1_name)
919      call netcdf_error_handler(ncstat)
920
921      ncstat = nf_put_att_text (nc_file_id, NF_GLOBAL, trim(grid2_ctmp),
922     &                          len_trim(grid2_name), grid2_name)
923      call netcdf_error_handler(ncstat)
924
925!-----------------------------------------------------------------------
926!
927!     prepare netCDF dimension info
928!
929!-----------------------------------------------------------------------
930
931      !***
932      !*** define grid size dimensions
933      !***
934
935      if (direction == 1) then
936        itmp1 = grid1_size
937        itmp2 = grid2_size
938      else
939        itmp1 = grid2_size
940        itmp2 = grid1_size
941      endif
942
943      ncstat = nf_def_dim (nc_file_id, 'n_a', itmp1, nc_srcgrdsize_id)
944      call netcdf_error_handler(ncstat)
945
946      ncstat = nf_def_dim (nc_file_id, 'n_b', itmp2, nc_dstgrdsize_id)
947      call netcdf_error_handler(ncstat)
948
949      !***
950      !*** define grid corner dimension
951      !***
952
953      if (direction == 1) then
954        itmp1 = grid1_corners
955        itmp2 = grid2_corners
956      else
957        itmp1 = grid2_corners
958        itmp2 = grid1_corners
959      endif
960
961      ncstat = nf_def_dim (nc_file_id, 'nv_a', itmp1, nc_srcgrdcorn_id)
962      call netcdf_error_handler(ncstat)
963
964      ncstat = nf_def_dim (nc_file_id, 'nv_b', itmp2, nc_dstgrdcorn_id)
965      call netcdf_error_handler(ncstat)
966
967      !***
968      !*** define grid rank dimension
969      !***
970
971      if (direction == 1) then
972        itmp1 = grid1_rank
973        itmp2 = grid2_rank
974      else
975        itmp1 = grid2_rank
976        itmp2 = grid1_rank
977      endif
978
979      ncstat = nf_def_dim (nc_file_id, 'src_grid_rank', 
980     &                     itmp1, nc_srcgrdrank_id)
981      call netcdf_error_handler(ncstat)
982
983      ncstat = nf_def_dim (nc_file_id, 'dst_grid_rank', 
984     &                     itmp2, nc_dstgrdrank_id)
985      call netcdf_error_handler(ncstat)
986
987      !***
988      !*** define first two dims as if 2-d cartesian domain
989      !***
990
991      if (direction == 1) then
992        itmp1 = grid1_dims(1)
993        if (grid1_rank > 1) then
994          itmp2 = grid1_dims(2)
995        else
996          itmp2 = 0
997        endif
998        itmp3 = grid2_dims(1)
999        if (grid2_rank > 1) then
1000          itmp4 = grid2_dims(2)
1001        else
1002          itmp4 = 0
1003        endif
1004      else
1005        itmp1 = grid2_dims(1)
1006        if (grid2_rank > 1) then
1007          itmp2 = grid2_dims(2)
1008        else
1009          itmp2 = 0
1010        endif
1011        itmp3 = grid1_dims(1)
1012        if (grid1_rank > 1) then
1013          itmp4 = grid1_dims(2)
1014        else
1015          itmp4 = 0
1016        endif
1017      endif
1018
1019      ncstat = nf_def_dim (nc_file_id, 'ni_a', itmp1, nc_src_isize_id)
1020      call netcdf_error_handler(ncstat)
1021
1022      ncstat = nf_def_dim (nc_file_id, 'nj_a', itmp2, nc_src_jsize_id)
1023      call netcdf_error_handler(ncstat)
1024
1025      ncstat = nf_def_dim (nc_file_id, 'ni_b', itmp3, nc_dst_isize_id)
1026      call netcdf_error_handler(ncstat)
1027
1028      ncstat = nf_def_dim (nc_file_id, 'nj_b', itmp4, nc_dst_jsize_id)
1029      call netcdf_error_handler(ncstat)
1030
1031      !***
1032      !*** define map size dimensions
1033      !***
1034
1035      if (direction == 1) then
1036        itmp1 = num_links_map1
1037      else
1038        itmp1 = num_links_map2
1039      endif
1040
1041      ncstat = nf_def_dim (nc_file_id, 'n_s', itmp1, nc_numlinks_id)
1042      call netcdf_error_handler(ncstat)
1043
1044      ncstat = nf_def_dim (nc_file_id, 'num_wgts', 
1045     &                     num_wts, nc_numwgts_id)
1046      call netcdf_error_handler(ncstat)
1047
1048      if (num_wts > 1) then
1049        ncstat = nf_def_dim (nc_file_id, 'num_wgts1', 
1050     &                       num_wts-1, nc_numwgts1_id)
1051        call netcdf_error_handler(ncstat)
1052      endif
1053
1054      !***
1055      !*** define grid dimensions
1056      !***
1057
1058      ncstat = nf_def_var (nc_file_id, 'src_grid_dims', NF_INT,
1059     &                     1, nc_srcgrdrank_id, nc_srcgrddims_id)
1060      call netcdf_error_handler(ncstat)
1061
1062      ncstat = nf_def_var (nc_file_id, 'dst_grid_dims', NF_INT,
1063     &                     1, nc_dstgrdrank_id, nc_dstgrddims_id)
1064      call netcdf_error_handler(ncstat)
1065
1066!-----------------------------------------------------------------------
1067!
1068!     define all arrays for netCDF descriptors
1069!
1070!-----------------------------------------------------------------------
1071
1072      !***
1073      !*** define grid center latitude array
1074      !***
1075
1076      ncstat = nf_def_var (nc_file_id, 'yc_a',
1077     &                     NF_DOUBLE, 1, nc_srcgrdsize_id, 
1078     &                     nc_srcgrdcntrlat_id)
1079      call netcdf_error_handler(ncstat)
1080
1081      ncstat = nf_def_var (nc_file_id, 'yc_b', 
1082     &                     NF_DOUBLE, 1, nc_dstgrdsize_id, 
1083     &                     nc_dstgrdcntrlat_id)
1084      call netcdf_error_handler(ncstat)
1085
1086      !***
1087      !*** define grid center longitude array
1088      !***
1089
1090      ncstat = nf_def_var (nc_file_id, 'xc_a', 
1091     &                     NF_DOUBLE, 1, nc_srcgrdsize_id, 
1092     &                     nc_srcgrdcntrlon_id)
1093      call netcdf_error_handler(ncstat)
1094
1095      ncstat = nf_def_var (nc_file_id, 'xc_b', 
1096     &                     NF_DOUBLE, 1, nc_dstgrdsize_id, 
1097     &                     nc_dstgrdcntrlon_id)
1098      call netcdf_error_handler(ncstat)
1099
1100      !***
1101      !*** define grid corner lat/lon arrays
1102      !***
1103
1104      nc_dims2_id(1) = nc_srcgrdcorn_id
1105      nc_dims2_id(2) = nc_srcgrdsize_id
1106
1107      ncstat = nf_def_var (nc_file_id, 'yv_a', 
1108     &                     NF_DOUBLE, 2, nc_dims2_id, 
1109     &                     nc_srcgrdcrnrlat_id)
1110      call netcdf_error_handler(ncstat)
1111
1112      ncstat = nf_def_var (nc_file_id, 'xv_a', 
1113     &                     NF_DOUBLE, 2, nc_dims2_id, 
1114     &                     nc_srcgrdcrnrlon_id)
1115      call netcdf_error_handler(ncstat)
1116
1117      nc_dims2_id(1) = nc_dstgrdcorn_id
1118      nc_dims2_id(2) = nc_dstgrdsize_id
1119
1120      ncstat = nf_def_var (nc_file_id, 'yv_b', 
1121     &                     NF_DOUBLE, 2, nc_dims2_id, 
1122     &                     nc_dstgrdcrnrlat_id)
1123      call netcdf_error_handler(ncstat)
1124
1125      ncstat = nf_def_var (nc_file_id, 'xv_b', 
1126     &                     NF_DOUBLE, 2, nc_dims2_id, 
1127     &                     nc_dstgrdcrnrlon_id)
1128      call netcdf_error_handler(ncstat)
1129
1130      !***
1131      !*** CSM wants all in degrees
1132      !***
1133
1134      grid1_units = 'degrees'
1135      grid2_units = 'degrees'
1136
1137      if (direction == 1) then
1138        grid1_ctmp = grid1_units
1139        grid2_ctmp = grid2_units
1140      else
1141        grid1_ctmp = grid2_units
1142        grid2_ctmp = grid1_units
1143      endif
1144
1145      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdcntrlat_id, 
1146     &                          'units', 7, grid1_ctmp)
1147      call netcdf_error_handler(ncstat)
1148
1149      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdcntrlat_id, 
1150     &                          'units', 7, grid2_ctmp)
1151      call netcdf_error_handler(ncstat)
1152
1153      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdcntrlon_id, 
1154     &                          'units', 7, grid1_ctmp)
1155      call netcdf_error_handler(ncstat)
1156
1157      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdcntrlon_id, 
1158     &                          'units', 7, grid2_ctmp)
1159      call netcdf_error_handler(ncstat)
1160
1161      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdcrnrlat_id, 
1162     &                          'units', 7, grid1_ctmp)
1163      call netcdf_error_handler(ncstat)
1164
1165      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdcrnrlon_id, 
1166     &                          'units', 7, grid1_ctmp)
1167      call netcdf_error_handler(ncstat)
1168
1169      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdcrnrlat_id, 
1170     &                          'units', 7, grid2_ctmp)
1171      call netcdf_error_handler(ncstat)
1172
1173      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdcrnrlon_id, 
1174     &                          'units', 7, grid2_ctmp)
1175      call netcdf_error_handler(ncstat)
1176
1177      !***
1178      !*** define grid mask
1179      !***
1180
1181      ncstat = nf_def_var (nc_file_id, 'mask_a', NF_INT,
1182     &                     1, nc_srcgrdsize_id, nc_srcgrdimask_id)
1183      call netcdf_error_handler(ncstat)
1184
1185      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdimask_id, 
1186     &                          'units', 8, 'unitless')
1187      call netcdf_error_handler(ncstat)
1188
1189      ncstat = nf_def_var (nc_file_id, 'mask_b', NF_INT,
1190     &                     1, nc_dstgrdsize_id, nc_dstgrdimask_id)
1191      call netcdf_error_handler(ncstat)
1192
1193      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdimask_id, 
1194     &                          'units', 8, 'unitless')
1195      call netcdf_error_handler(ncstat)
1196
1197      !***
1198      !*** define grid area arrays
1199      !***
1200
1201      ncstat = nf_def_var (nc_file_id, 'area_a', 
1202     &                     NF_DOUBLE, 1, nc_srcgrdsize_id, 
1203     &                     nc_srcgrdarea_id)
1204      call netcdf_error_handler(ncstat)
1205
1206      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdarea_id, 
1207     &                          'units', 14, 'square radians')
1208      call netcdf_error_handler(ncstat)
1209
1210      ncstat = nf_def_var (nc_file_id, 'area_b', 
1211     &                     NF_DOUBLE, 1, nc_dstgrdsize_id, 
1212     &                     nc_dstgrdarea_id)
1213      call netcdf_error_handler(ncstat)
1214
1215      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdarea_id, 
1216     &                          'units', 14, 'square radians')
1217      call netcdf_error_handler(ncstat)
1218
1219      !***
1220      !*** define grid fraction arrays
1221      !***
1222
1223      ncstat = nf_def_var (nc_file_id, 'frac_a', 
1224     &                     NF_DOUBLE, 1, nc_srcgrdsize_id, 
1225     &                     nc_srcgrdfrac_id)
1226      call netcdf_error_handler(ncstat)
1227
1228      ncstat = nf_put_att_text (nc_file_id, nc_srcgrdfrac_id, 
1229     &                          'units', 8, 'unitless')
1230      call netcdf_error_handler(ncstat)
1231
1232      ncstat = nf_def_var (nc_file_id, 'frac_b', 
1233     &                     NF_DOUBLE, 1, nc_dstgrdsize_id, 
1234     &                     nc_dstgrdfrac_id)
1235      call netcdf_error_handler(ncstat)
1236
1237      ncstat = nf_put_att_text (nc_file_id, nc_dstgrdfrac_id, 
1238     &                          'units', 8, 'unitless')
1239      call netcdf_error_handler(ncstat)
1240
1241      !***
1242      !*** define mapping arrays
1243      !***
1244
1245      ncstat = nf_def_var (nc_file_id, 'col', 
1246     &                     NF_INT, 1, nc_numlinks_id, 
1247     &                     nc_srcadd_id)
1248      call netcdf_error_handler(ncstat)
1249
1250      ncstat = nf_def_var (nc_file_id, 'row', 
1251     &                     NF_INT, 1, nc_numlinks_id, 
1252     &                     nc_dstadd_id)
1253      call netcdf_error_handler(ncstat)
1254
1255      ncstat = nf_def_var (nc_file_id, 'S', 
1256     &                     NF_DOUBLE, 1, nc_numlinks_id, 
1257     &                     nc_rmpmatrix_id)
1258      call netcdf_error_handler(ncstat)
1259
1260      if (num_wts > 1) then
1261        nc_dims2_id(1) = nc_numwgts1_id
1262        nc_dims2_id(2) = nc_numlinks_id
1263
1264        ncstat = nf_def_var (nc_file_id, 'S2', 
1265     &                     NF_DOUBLE, 2, nc_dims2_id, 
1266     &                     nc_rmpmatrix2_id)
1267        call netcdf_error_handler(ncstat)
1268      endif
1269
1270      !***
1271      !*** end definition stage
1272      !***
1273
1274      ncstat = nf_enddef(nc_file_id)
1275      call netcdf_error_handler(ncstat)
1276
1277!-----------------------------------------------------------------------
1278!
1279!     compute integer masks
1280!
1281!-----------------------------------------------------------------------
1282
1283      if (direction == 1) then
1284        allocate (src_mask_int(grid1_size),
1285     &            dst_mask_int(grid2_size))
1286
1287        where (grid2_mask)
1288          dst_mask_int = 1
1289        elsewhere
1290          dst_mask_int = 0
1291        endwhere
1292
1293        where (grid1_mask)
1294          src_mask_int = 1
1295        elsewhere
1296          src_mask_int = 0
1297        endwhere
1298      else
1299        allocate (src_mask_int(grid2_size),
1300     &            dst_mask_int(grid1_size))
1301
1302        where (grid1_mask)
1303          dst_mask_int = 1
1304        elsewhere
1305          dst_mask_int = 0
1306        endwhere
1307
1308        where (grid2_mask)
1309          src_mask_int = 1
1310        elsewhere
1311          src_mask_int = 0
1312        endwhere
1313      endif
1314
1315!-----------------------------------------------------------------------
1316!
1317!     change units of lat/lon coordinates if input units different
1318!     from radians. if this is the second mapping, the conversion has
1319!     alread been done.
1320!
1321!-----------------------------------------------------------------------
1322
1323      if (grid1_units(1:7) == 'degrees' .and. direction == 1) then
1324        grid1_center_lat = grid1_center_lat/deg2rad
1325        grid1_center_lon = grid1_center_lon/deg2rad
1326        grid1_corner_lat = grid1_corner_lat/deg2rad
1327        grid1_corner_lon = grid1_corner_lon/deg2rad
1328      endif
1329
1330      if (grid2_units(1:7) == 'degrees' .and. direction == 1) then
1331        grid2_center_lat = grid2_center_lat/deg2rad
1332        grid2_center_lon = grid2_center_lon/deg2rad
1333        grid2_corner_lat = grid2_corner_lat/deg2rad
1334        grid2_corner_lon = grid2_corner_lon/deg2rad
1335      endif
1336
1337!-----------------------------------------------------------------------
1338!
1339!     write mapping data
1340!
1341!-----------------------------------------------------------------------
1342
1343      if (direction == 1) then
1344        itmp1 = nc_srcgrddims_id
1345        itmp2 = nc_dstgrddims_id
1346      else
1347        itmp2 = nc_srcgrddims_id
1348        itmp1 = nc_dstgrddims_id
1349      endif
1350
1351      ncstat = nf_put_var_int(nc_file_id, itmp1, grid1_dims)
1352      call netcdf_error_handler(ncstat)
1353
1354      ncstat = nf_put_var_int(nc_file_id, itmp2, grid2_dims)
1355      call netcdf_error_handler(ncstat)
1356
1357      ncstat = nf_put_var_int(nc_file_id, nc_srcgrdimask_id, 
1358     &                        src_mask_int)
1359      call netcdf_error_handler(ncstat)
1360
1361      ncstat = nf_put_var_int(nc_file_id, nc_dstgrdimask_id,
1362     &                        dst_mask_int)
1363      call netcdf_error_handler(ncstat)
1364
1365      deallocate(src_mask_int, dst_mask_int)
1366
1367      if (direction == 1) then
1368        itmp1 = nc_srcgrdcntrlat_id
1369        itmp2 = nc_srcgrdcntrlon_id
1370        itmp3 = nc_srcgrdcrnrlat_id
1371        itmp4 = nc_srcgrdcrnrlon_id
1372      else
1373        itmp1 = nc_dstgrdcntrlat_id
1374        itmp2 = nc_dstgrdcntrlon_id
1375        itmp3 = nc_dstgrdcrnrlat_id
1376        itmp4 = nc_dstgrdcrnrlon_id
1377      endif
1378
1379      ncstat = nf_put_var_double(nc_file_id, itmp1, grid1_center_lat)
1380      call netcdf_error_handler(ncstat)
1381
1382      ncstat = nf_put_var_double(nc_file_id, itmp2, grid1_center_lon)
1383      call netcdf_error_handler(ncstat)
1384
1385      ncstat = nf_put_var_double(nc_file_id, itmp3, grid1_corner_lat)
1386      call netcdf_error_handler(ncstat)
1387
1388      ncstat = nf_put_var_double(nc_file_id, itmp4, grid1_corner_lon)
1389      call netcdf_error_handler(ncstat)
1390
1391      if (direction == 1) then
1392        itmp1 = nc_dstgrdcntrlat_id
1393        itmp2 = nc_dstgrdcntrlon_id
1394        itmp3 = nc_dstgrdcrnrlat_id
1395        itmp4 = nc_dstgrdcrnrlon_id
1396      else
1397        itmp1 = nc_srcgrdcntrlat_id
1398        itmp2 = nc_srcgrdcntrlon_id
1399        itmp3 = nc_srcgrdcrnrlat_id
1400        itmp4 = nc_srcgrdcrnrlon_id
1401      endif
1402
1403      ncstat = nf_put_var_double(nc_file_id, itmp1, grid2_center_lat)
1404      call netcdf_error_handler(ncstat)
1405
1406      ncstat = nf_put_var_double(nc_file_id, itmp2, grid2_center_lon)
1407      call netcdf_error_handler(ncstat)
1408
1409      ncstat = nf_put_var_double(nc_file_id, itmp3, grid2_corner_lat)
1410      call netcdf_error_handler(ncstat)
1411
1412      ncstat = nf_put_var_double(nc_file_id, itmp4, grid2_corner_lon)
1413      call netcdf_error_handler(ncstat)
1414
1415      if (direction == 1) then
1416        itmp1 = nc_srcgrdarea_id
1417        itmp2 = nc_srcgrdfrac_id
1418        itmp3 = nc_dstgrdarea_id
1419        itmp4 = nc_dstgrdfrac_id
1420      else
1421        itmp1 = nc_dstgrdarea_id
1422        itmp2 = nc_dstgrdfrac_id
1423        itmp3 = nc_srcgrdarea_id
1424        itmp4 = nc_srcgrdfrac_id
1425      endif
1426
1427      if (luse_grid1_area) then
1428        ncstat = nf_put_var_double(nc_file_id, itmp1, grid1_area_in)
1429      else
1430        ncstat = nf_put_var_double(nc_file_id, itmp1, grid1_area)
1431      endif
1432      call netcdf_error_handler(ncstat)
1433
1434      ncstat = nf_put_var_double(nc_file_id, itmp2, grid1_frac)
1435      call netcdf_error_handler(ncstat)
1436
1437      if (luse_grid2_area) then
1438        ncstat = nf_put_var_double(nc_file_id, itmp3, grid2_area)
1439      else
1440        ncstat = nf_put_var_double(nc_file_id, itmp3, grid2_area)
1441      endif
1442      call netcdf_error_handler(ncstat)
1443
1444      ncstat = nf_put_var_double(nc_file_id, itmp4, grid2_frac)
1445      call netcdf_error_handler(ncstat)
1446
1447      if (direction == 1) then
1448        ncstat = nf_put_var_int(nc_file_id, nc_srcadd_id, 
1449     &                          grid1_add_map1)
1450        call netcdf_error_handler(ncstat)
1451
1452        ncstat = nf_put_var_int(nc_file_id, nc_dstadd_id, 
1453     &                          grid2_add_map1)
1454        call netcdf_error_handler(ncstat)
1455
1456        if (num_wts == 1) then
1457          ncstat = nf_put_var_double(nc_file_id, nc_rmpmatrix_id, 
1458     &                               wts_map1)
1459          call netcdf_error_handler(ncstat)
1460        else
1461          allocate(wts1(num_links_map1),wts2(num_wts-1,num_links_map1))
1462
1463          wts1 = wts_map1(1,:)
1464          wts2 = wts_map1(2:,:)
1465
1466          ncstat = nf_put_var_double(nc_file_id, nc_rmpmatrix_id, 
1467     &                               wts1)
1468          call netcdf_error_handler(ncstat)
1469          ncstat = nf_put_var_double(nc_file_id, nc_rmpmatrix2_id, 
1470     &                               wts2)
1471          call netcdf_error_handler(ncstat)
1472          deallocate(wts1,wts2)
1473        endif
1474      else
1475        ncstat = nf_put_var_int(nc_file_id, nc_srcadd_id, 
1476     &                          grid2_add_map2)
1477        call netcdf_error_handler(ncstat)
1478
1479        ncstat = nf_put_var_int(nc_file_id, nc_dstadd_id, 
1480     &                          grid1_add_map2)
1481        call netcdf_error_handler(ncstat)
1482
1483        if (num_wts == 1) then
1484          ncstat = nf_put_var_double(nc_file_id, nc_rmpmatrix_id, 
1485     &                               wts_map2)
1486          call netcdf_error_handler(ncstat)
1487        else
1488          allocate(wts1(num_links_map2),wts2(num_wts-1,num_links_map2))
1489
1490          wts1 = wts_map2(1,:)
1491          wts2 = wts_map2(2:,:)
1492
1493          ncstat = nf_put_var_double(nc_file_id, nc_rmpmatrix_id, 
1494     &                               wts1)
1495          call netcdf_error_handler(ncstat)
1496          ncstat = nf_put_var_double(nc_file_id, nc_rmpmatrix2_id, 
1497     &                               wts2)
1498          call netcdf_error_handler(ncstat)
1499          deallocate(wts1,wts2)
1500        endif
1501      endif
1502
1503      ncstat = nf_close(nc_file_id)
1504      call netcdf_error_handler(ncstat)
1505
1506!-----------------------------------------------------------------------
1507
1508      end subroutine write_remap_csm
1509
1510!***********************************************************************
1511
1512      end module remap_write
1513
1514!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1515
Note: See TracBrowser for help on using the repository browser.