external/io_esmf/ext_esmf_write_field.F90

   1
   2
   3 !TODO:  remove duplication between ext_esmf_read_field and
   4 !TODO:  ext_esmf_write_field
   5
   6 !TODO:  how to deal with time?  (via current ESMF_Clock)
   7 !TODO:  to begin, use it as an error check...
   8
   9
  10 !--- write_field
  11 SUBROUTINE ext_esmf_write_field ( DataHandle , DateStr , VarName , Field , FieldType , Comm , IOComm, &
  12                                   DomainDesc , MemoryOrder , Stagger , DimNames ,              &
  13                                   DomainStart , DomainEnd ,                                    &
  14                                   MemoryStart , MemoryEnd ,                                    &
  15                                   PatchStart , PatchEnd ,                                      &
  16                                   Status )
  17   USE module_ext_esmf
  18   IMPLICIT NONE
  19   INTEGER       ,INTENT(IN)    :: DataHandle
  20   CHARACTER*(*) ,intent(inout) :: DateStr
  21   CHARACTER*(*) ,intent(inout) :: VarName
  22   integer       ,intent(inout) :: FieldType
  23   integer       ,intent(inout) :: Comm
  24   integer       ,intent(inout) :: IOComm
  25   integer       ,intent(inout) :: DomainDesc
  26   character*(*) ,intent(inout) :: MemoryOrder
  27   character*(*) ,intent(inout) :: Stagger
  28   character*(*) ,intent(inout) :: DimNames(*)
  29   integer       ,intent(inout) :: DomainStart(*), DomainEnd(*)
  30   integer       ,intent(inout) :: MemoryStart(*), MemoryEnd(*)
  31   integer       ,intent(inout) :: PatchStart(*),  PatchEnd(*)
  32   REAL          ,INTENT(INOUT) :: Field(*)
  33   integer       ,intent(out)   :: Status
  34   ! Local declarations
  35   INTEGER :: ids,ide,jds,jde,kds,kde
  36   INTEGER :: ims,ime,jms,jme,kms,kme
  37   INTEGER :: ips,ipe,jps,jpe,kps,kpe
  38   TYPE(ESMF_State), POINTER :: exportstate
  39   TYPE(ESMF_Field) :: tmpField
  40   TYPE(ESMF_Array) :: tmpArray
  41   TYPE(ESMF_ArraySpec) :: arrayspec
  42 !  TYPE(ESMF_DataKind) :: esmf_kind
  43   INTEGER :: esmf_kind
  44   REAL(ESMF_KIND_R4), POINTER :: data_esmf_real_ptr(:,:)
  45   REAL(ESMF_KIND_R4), POINTER :: tmp_esmf_r4_ptr(:,:)
  46   INTEGER(ESMF_KIND_I4), POINTER :: data_esmf_int_ptr(:,:)
  47   INTEGER, PARAMETER :: esmf_rank = 2
  48   INTEGER :: DomainEndFull(esmf_rank), idefull, jdefull, ict, i, j
  49   INTEGER :: PatchEndFull(esmf_rank), ipefull, jpefull
  50   ! esmf_counts is redundant.  remove it as soon as ESMF_ArrayCreate no
  51   ! longer requires it
  52   INTEGER :: esmf_counts(esmf_rank)
  53   INTEGER :: rc, debug_level
  54   LOGICAL, EXTERNAL :: has_char
  55   character*256 mess
  56
  57   CALL get_wrf_debug_level( debug_level )
  58
  59   IF ( .NOT. int_valid_handle( DataHandle ) ) THEN
  60     CALL wrf_error_fatal("ext_esmf_write_field: invalid data handle" )
  61   ENDIF
  62   IF ( .NOT. int_handle_in_use( DataHandle ) ) THEN
  63     CALL wrf_error_fatal("ext_esmf_write_field: DataHandle not opened" )
  64   ENDIF
  65   IF ( .NOT. opened_for_write( DataHandle ) ) THEN
  66     CALL wrf_error_fatal("ext_esmf_write_field: DataHandle not opened for write" )
  67   ENDIF
  68
  69 write(mess,*)'ext_esmf_write_field ',DataHandle, TRIM(DateStr), TRIM(VarName)
  70 call wrf_debug( 300, TRIM(mess) )
  71
  72   IF      ( FieldType .EQ. WRF_REAL ) THEN
  73     esmf_kind = ESMF_KIND_R4
  74   ELSE IF ( FieldType .EQ. WRF_DOUBLE ) THEN
  75     CALL wrf_error_fatal( 'ext_esmf_write_field, WRF_DOUBLE not yet supported')
  76   ELSE IF ( FieldType .EQ. WRF_INTEGER ) THEN
  77     esmf_kind = ESMF_KIND_I4
  78 !TODO:  implement this (below)
  79     CALL wrf_error_fatal( 'ext_esmf_write_field, WRF_INTEGER not yet implemented')
  80   ELSE IF ( FieldType .EQ. WRF_LOGICAL ) THEN
  81     CALL wrf_error_fatal( 'ext_esmf_write_field, WRF_LOGICAL not yet supported')
  82   ENDIF
  83
  84   ims = MemoryStart(1) ; ime = MemoryEnd(1)
  85   jms = MemoryStart(2) ; jme = MemoryEnd(2)
  86   kms = MemoryStart(3) ; kme = MemoryEnd(3)
  87
  88   ips = PatchStart(1) ; ipe = PatchEnd(1)
  89   jps = PatchStart(2) ; jpe = PatchEnd(2)
  90   kps = PatchStart(3) ; kpe = PatchEnd(3)
  91
  92   ids = DomainStart(1) ; ide = DomainEnd(1)
  93   jds = DomainStart(2) ; jde = DomainEnd(2)
  94   kds = DomainStart(3) ; kde = DomainEnd(3)
  95
  96   ! For now, treat all arrays as 2D...
  97 !TODO:  Eventually, use ../io_netcdf subroutines Transpose() and reorder()
  98 !TODO:  (and etc.) to handle general array ranks and index orderings.
  99 !TODO:  Some copies of these exist in ../../frame/module_io.F.
 100 !TODO:  Then use ESMF_ArrayDataMap class to handle index mapping.
 101   IF ( kms /= kme ) THEN
 102     write(mess,*)'ext_esmf_write_field: rank > 2 not yet supported. field = ',TRIM(VarName)
 103     CALL wrf_error_fatal( mess )
 104   ENDIF
 105
 106 ! The non-staggered variables come in at one-less than
 107 ! domain dimensions, but io_esmf is currently hacked to use full
 108 ! domain spec, so adjust if not staggered.
 109 !TODO:  Remove EndFull hackery once ESMF can support staggered
 110 !TODO:  grids in regional models.  (This hack works around the current
 111 !TODO:  need to use only larger staggered dimensions for ESMF_Arrays.)
 112   CALL ioesmf_endfullhack( esmf_rank, DomainEnd, PatchEnd, Stagger, &
 113                            DomainEndFull, PatchEndFull )
 114   idefull = DomainEndFull(1)
 115   jdefull = DomainEndFull(2)
 116   ipefull = PatchEndFull(1)
 117   jpefull = PatchEndFull(2)
 118
 119 write(mess,*) ' ext_esmf_write_field: okay_to_write: ', DataHandle, okay_to_write(DataHandle)
 120 call wrf_debug( 300, TRIM(mess) )
 121
 122   ! case 1: the file is opened for write but not committed ("training")
 123   IF ( .NOT. okay_to_write( DataHandle ) )  THEN
 124
 125     ! Training:  build the ESMF export state
 126 write(mess,*) ' ext_esmf_write_field: TRAINING WRITE:  DataHandle = ', DataHandle
 127 call wrf_debug( 300, TRIM(mess) )
 128
 129     ! First, build the ESMF_Grid for this DataHandle, if it does not
 130     ! already exist
 131 write(0,*)__FILE__,__LINE__,'ext_esmf_write_field Stagger',TRIM(Stagger)
 132 write(0,*)__FILE__,__LINE__,'ext_esmf_write_field VarName',TRIM(VarName)
 133 write(0,*)__FILE__,__LINE__,'ext_esmf_write_field DomainEnd ', DomainEnd(1:esmf_rank)
 134 write(0,*)__FILE__,__LINE__,'ext_esmf_write_field PatchEnd  ', PatchEnd(1:esmf_rank)
 135     CALL ioesmf_create_grid( DataHandle, esmf_rank, MemoryOrder, Stagger,      &
 136                              DomainStart(1:esmf_rank), DomainEnd(1:esmf_rank), &
 137                              MemoryStart(1:esmf_rank), MemoryEnd(1:esmf_rank), &
 138                              PatchStart(1:esmf_rank), PatchEnd(1:esmf_rank) )
 139     ! Grab the current exportState and add to it...
 140     CALL ESMF_ExportStateGetCurrent( exportstate, rc )
 141     IF ( rc /= ESMF_SUCCESS ) THEN
 142       CALL wrf_error_fatal("ext_esmf_write_field, training:  ESMF_ExportStateGetCurrent failed" )
 143     ENDIF
 144 ! BEGIN DOESNOTWORK
 145 ! The following code does not work for reasons as-yet unknown.
 146 ! A likely suspect is lbounds and ubounds which fail in other interfaces in
 147 ! ESMF 2.2.0rp1 ...
 148     ! Build ESMF objects...
 149     ! Build an ESMF_ArraySpec.  The use of ESMF_ArraySpec and ESMF_Array
 150     ! objects allows some of the code that follows to be type-kind-independent.
 151
 152     ! Build an ESMF_Array
 153     ! Implementation note:  since we do not yet have full control over how
 154     ! ESMF chooses to lay out a "patch" within "memory", we must copy by
 155     ! hand.  (Reasons include lack of support in ESMF for asymmetric halos,
 156     ! addition of "extra" rows/columns to optimize alignment on some machines,
 157     ! handling of periodic boundary conditions, etc.)  Thus, there
 158     ! is no point in using larger "memory" sizes to build the array -- patch
 159     ! is fine.  Also, since we must copy anyway, might as well let ESMF manage
 160     ! the memory for simplicity.
 161 !TODO:  Once ESMF can match WRF memory-patch mapping, replace this with a more
 162 !TODO:  efficient solution that does not require a copy.
 163 !TODO:  esmf_counts is redundant.  Remove it as soon as ESMF_ArrayCreate no
 164 !TODO:  longer requires it.
 165 !    esmf_counts(1:esmf_rank) = DomainEndFull(1:esmf_rank) - &
 166 !                               DomainStart(1:esmf_rank) + 1
 167 !    tmpArray = ESMF_ArrayCreate(arrayspec, counts=esmf_counts,      &
 168 !                                lbounds=DomainStart(1:esmf_rank),   &
 169 !                                ubounds=DomainEndFull(1:esmf_rank), &
 170 !                                rc=rc)
 171 !    IF ( rc /= ESMF_SUCCESS ) THEN
 172 !      WRITE(mess,*) ' ext_esmf_write_field: ESMF_ArrayCreate failed, rc = ', rc
 173 !      CALL wrf_error_fatal( TRIM(mess) )
 174 !    ENDIF
 175     ! Determine grid staggering for this Field
 176 !    IF ( has_char( Stagger, 'x' ) .AND. has_char( Stagger, 'y' ) ) THEN
 177 !      CALL wrf_error_fatal( &
 178 !        "ext_esmf_write_field:  ESMF does not yet support XY staggering for C-grid" )
 179 !    ELSE IF ( has_char( Stagger, 'x' ) ) THEN
 180 !      horzrelloc=ESMF_CELL_WFACE
 181 !    ELSE IF ( has_char( Stagger, 'y' ) ) THEN
 182 !      horzrelloc=ESMF_CELL_SFACE
 183 !    ELSE
 184 !      horzrelloc=ESMF_CELL_CENTER
 185 !    ENDIF
 186     ! Build an ESMF_Field
 187     ! Note:  though it is counter-intuitive, ESMF uses
 188     ! shallow-copy-masquerading-as-reference to implement the
 189     ! pseudo-equivalent of POINTER assignment under-the-hood.  What this means
 190     ! here is that it is OK to pass deep object tmpArray into
 191     ! ESMF_FieldCreate() and then return from this subroutine.  Even though
 192     ! tmpArray goes out of scope, it is OK.  However, if tmpArray were to be
 193     ! modified after this call, the changes would not be guaranteed to always
 194     ! appear in tmpField.  It works that way now, but ESMF Core team has
 195     ! plans that may make it break in the future.  Build-it, attach-it,
 196     ! flush-it will work.  Build-it, attach-it, modify-it, flush-it may not
 197     ! always work.
 198     ! Note:  unique Field name is required by ESMF_StateAdd().
 199 !TODO:  use CF "standard_name" once the WRF Registry supports it
 200 !    tmpField = ESMF_FieldCreate( grid( DataHandle )%ptr, tmpArray,          &
 201 !                                 copyflag=ESMF_DATA_REF,                    &
 202 !                                 horzrelloc=horzrelloc, name=TRIM(VarName), &
 203 !                                 rc=rc )
 204 ! END DOESNOTWORK
 205 !TODO:  Compute horzrelloc from Stagger as above once ESMF supports staggering
 206 !    horzrelloc=ESMF_CELL_CENTER
 207 !TODO:  Add code for other data types here...
 208 !    ALLOCATE( tmp_esmf_r4_ptr(ips:ipefull,jps:jpefull) )
 209     ALLOCATE( tmp_esmf_r4_ptr(ips:ipe,jps:jpe) )
 210     CALL wrf_debug ( 100, 'ext_esmf_write_field: calling ESMF_FieldCreate' )
 211     tmpField = ESMF_FieldCreate(         &
 212                  grid( DataHandle )%ptr, &
 213                  tmp_esmf_r4_ptr,        &
 214                  datacopyflag=ESMF_DATACOPY_REFERENCE, &
 215                  staggerloc=ESMF_STAGGERLOC_CENTER,    &
 216                  name=TRIM(VarName),     &
 217                  rc=rc )
 218     IF ( rc /= ESMF_SUCCESS ) THEN
 219       WRITE(mess,*) ' ext_esmf_write_field: ESMF_FieldCreate failed, rc = ', rc
 220       CALL wrf_error_fatal( TRIM(mess) )
 221     ENDIF
 222     CALL wrf_debug ( 100, 'ext_esmf_write_field: back from ESMF_FieldCreate' )
 223     WRITE(mess,*) 'ext_esmf_write_field: tmp_esmf_r4_ptr(',        &
 224       LBOUND(tmp_esmf_r4_ptr,1),':',UBOUND(tmp_esmf_r4_ptr,1),',', &
 225       LBOUND(tmp_esmf_r4_ptr,2),':',UBOUND(tmp_esmf_r4_ptr,2),')'
 226     CALL wrf_debug ( 100 , TRIM(mess) )
 227     ! Add the Field to the export state...
 228 !TODO:  for now, just build ESMF_Fields and stuff them in
 229 !TODO:  later, use a single ESMF_Bundle
 230     CALL ESMF_StateAdd( exportstate, (/tmpField/), rc=rc ) ! 5.2.0r only accepts arrays; use array constructor
 231     IF ( rc /= ESMF_SUCCESS ) THEN
 232       CALL wrf_error_fatal("ext_esmf_write_field:  ESMF_StateAddfailed" )
 233     ENDIF
 234 write(mess,*) ' ext_esmf_write_field: END TRAINING WRITE:  DataHandle = ', DataHandle
 235 call wrf_debug( 300, TRIM(mess) )
 236
 237   ! case 2: opened for write and committed
 238   ELSE IF ( okay_to_write( DataHandle ) )  THEN
 239
 240 write(mess,*) ' ext_esmf_write_field: ACTUAL WRITE:  DataHandle = ', DataHandle
 241 call wrf_debug( 300, TRIM(mess) )
 242
 243     ! write:  insert data into the ESMF export state
 244     ! Grab the current exportState
 245     CALL ESMF_ExportStateGetCurrent( exportstate, rc )
 246     IF ( rc /= ESMF_SUCCESS ) THEN
 247       CALL wrf_error_fatal("ext_esmf_write_field:  ESMF_ExportStateGetCurrent failed" )
 248     ENDIF
 249     ! grab the Field
 250     CALL ESMF_StateGet( exportstate, itemName=TRIM(VarName), &
 251                              field=tmpfield, rc=rc )
 252     IF ( rc /= ESMF_SUCCESS ) THEN
 253       CALL wrf_error_fatal("ext_esmf_write_field:  ESMF_StateGetfailed" )
 254     ENDIF
 255
 256 CALL wrf_debug ( 100, 'ext_esmf_write_field '//TRIM(VarName)//':  calling ESMF_FieldPrint( tmpField ) 1' )
 257 IF ( 100 .LE. debug_level ) THEN
 258   CALL ESMF_FieldPrint( tmpField, rc=rc )
 259 ENDIF
 260 CALL wrf_debug ( 100, 'ext_esmf_write_field '//TRIM(VarName)//':  back from ESMF_FieldPrint( tmpField ) 1' )
 261
 262     ! grab a pointer to the export state data and copy data from Field
 263     IF      ( FieldType .EQ. WRF_REAL ) THEN
 264       CALL ESMF_FieldGet( tmpField, 0, data_esmf_real_ptr, rc=rc )
 265       IF ( rc /= ESMF_SUCCESS ) THEN
 266         CALL wrf_error_fatal("ext_esmf_write_field:  ESMF_FieldGetDataPointer(r4) failed" )
 267       ENDIF
 268       IF ( ( PatchStart(1)   /= LBOUND(data_esmf_real_ptr,1) ) .OR. &
 269            ( PatchEnd(1) /= UBOUND(data_esmf_real_ptr,1) ) .OR. &
 270            ( PatchStart(2)   /= LBOUND(data_esmf_real_ptr,2) ) .OR. &
 271            ( PatchEnd(2) /= UBOUND(data_esmf_real_ptr,2) ) ) THEN
 272         WRITE( mess,* ) 'ESMF_FieldGetDataPointer bounds mismatch',          &
 273           __FILE__ ,                                                         &
 274           ', line ',                                                         &
 275           __LINE__ ,                                                         &
 276           ', ips:ipe,jps:jpe = ',PatchStart(1),':',PatchEnd(1),',',      &
 277                                  PatchStart(2),':',PatchEnd(2),          &
 278           ', data_esmf_real_ptr(BOUNDS) = ',                                 &
 279           LBOUND(data_esmf_real_ptr,1),':',UBOUND(data_esmf_real_ptr,1),',', &
 280           LBOUND(data_esmf_real_ptr,2),':',UBOUND(data_esmf_real_ptr,2)
 281         CALL wrf_error_fatal ( TRIM(mess) )
 282       ENDIF
 283
 284 WRITE( mess,* ) 'DEBUG:  ext_esmf_write_field:  ips:ipe,jps:jpe = ', &
 285   ips,':',ipe,',',jps,':',jpe,                                       &
 286   ', data_esmf_real_ptr(BOUNDS) = ',                                 &
 287   LBOUND(data_esmf_real_ptr,1),':',UBOUND(data_esmf_real_ptr,1),',', &
 288   LBOUND(data_esmf_real_ptr,2),':',UBOUND(data_esmf_real_ptr,2)
 289 CALL wrf_debug( 300, TRIM(mess) )
 290
 291       CALL ioesmf_insert_data_real( Field, data_esmf_real_ptr,            &
 292                                     ips, ipe, jps, jpe, kps, kpe, &
 293                                     ims, ime, jms, jme, kms, kme )
 294     ELSE IF ( FieldType .EQ. WRF_INTEGER ) THEN
 295       CALL ESMF_FieldGet( tmpField, 0, data_esmf_int_ptr, rc=rc )
 296       IF ( rc /= ESMF_SUCCESS ) THEN
 297         CALL wrf_error_fatal("ext_esmf_write_field:  ESMF_FieldGetDataPointer(i4) failed" )
 298       ENDIF
 299       IF ( ( PatchStart(1)   /= LBOUND(data_esmf_int_ptr,1) ) .OR. &
 300            ( PatchEnd(1) /= UBOUND(data_esmf_int_ptr,1) ) .OR. &
 301            ( PatchStart(2)   /= LBOUND(data_esmf_int_ptr,2) ) .OR. &
 302            ( PatchEnd(2) /= UBOUND(data_esmf_int_ptr,2) ) ) THEN
 303         WRITE( mess,* ) 'ESMF_FieldGetDataPointer bounds mismatch',        &
 304           __FILE__ ,                                                       &
 305           ', line ',                                                       &
 306           __LINE__ ,                                                       &
 307           ', ips:ipe,jps:jpe = ',PatchStart(1),':',PatchEnd(1),',',    &
 308                                  PatchStart(2),':',PatchEnd(2),        &
 309           ', data_esmf_int_ptr(BOUNDS) = ',                                &
 310           LBOUND(data_esmf_int_ptr,1),':',UBOUND(data_esmf_int_ptr,1),',', &
 311           LBOUND(data_esmf_int_ptr,2),':',UBOUND(data_esmf_int_ptr,2)
 312         CALL wrf_error_fatal ( TRIM(mess) )
 313       ENDIF
 314       CALL ioesmf_insert_data_int( Field, data_esmf_int_ptr,             &
 315                                    ips, ipe, jps, jpe, kps, kpe, &
 316                                    ims, ime, jms, jme, kms, kme )
 317     ENDIF
 318 write(mess,*) ' ext_esmf_write_field: END ACTUAL WRITE:  DataHandle = ', DataHandle
 319 call wrf_debug( 300, TRIM(mess) )
 320
 321   ENDIF
 322
 323 CALL wrf_debug ( 100, 'ext_esmf_write_field '//TRIM(VarName)//':  calling ESMF_FieldPrint( tmpField ) 2' )
 324 IF ( 100 .LE. debug_level ) THEN
 325   CALL ESMF_FieldPrint( tmpField, rc=rc )
 326 ENDIF
 327 CALL wrf_debug ( 100, 'ext_esmf_write_field '//TRIM(VarName)//':  back from ESMF_FieldPrint( tmpField ) 2' )
 328
 329   Status = 0
 330
 331   RETURN
 332
 333 END SUBROUTINE ext_esmf_write_field
 334