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