XMDF  2.2
Tests.f90

Tests.f90 provides Fortran code examples for many XMDF APIs.

1 MODULE testdatasets
2 
3 USE xmdf
4 
5 CHARACTER(LEN=*), PARAMETER :: datasets_location = 'Datasets'
6 CHARACTER(LEN=*), PARAMETER :: scalar_a_location = 'Scalars/ScalarA'
7 CHARACTER(LEN=*), PARAMETER :: scalar_a_location_full = 'Datasets/Scalars/ScalarA'
8 CHARACTER(LEN=*), PARAMETER :: scalar_b_location = 'Scalars/ScalarB'
9 CHARACTER(LEN=*), PARAMETER :: vector2d_a_location = 'Vectors/Vector2D_A'
10 CHARACTER(LEN=*), PARAMETER :: vector2d_b_location = 'Vectors/Vector2D_B'
11 
12 CONTAINS
13 ! ---------------------------------------------------------------------------
14 ! FUNCTION tdEditScalarAValues
15 ! PURPOSE
16 ! NOTES
17 ! ---------------------------------------------------------------------------
18 SUBROUTINE td_edit_scalar_a_values(a_Filename, a_Compression, error)
19  CHARACTER(LEN=*), INTENT(IN) :: a_filename
20  INTEGER, INTENT(IN) :: a_compression
21  INTEGER, INTENT(OUT) :: error
22  INTEGER xfileid, xscalarid
23  INTEGER, PARAMETER :: editnumvalues = 3
24  INTEGER edittimestep
25  INTEGER, DIMENSION(editNumValues) :: indices
26  REAL, DIMENSION(editNumValues) :: new_values
27 
28  CALL td_write_scalar_a(a_filename, a_compression, error)
29  if (error < 0) then
30  return
31  endif
32 
33  ! open the file and edit the values
34  CALL xf_open_file(a_filename, .false., xfileid, error)
35  if (error < 0) then
36  return
37  endif
38 
39  CALL xf_open_group(xfileid, scalar_a_location_full, xscalarid, error);
40  if (error < 0) then
41  CALL xf_close_file(xfileid, error)
42  return
43  endif
44 
45  ! Edit values in timestep 1, make index 1 = 4, index 5 = 40,
46  ! and index 10 = 400
47  edittimestep = 1
48  indices(1) = 1;
49  indices(2) = 5;
50  indices(3) = 10;
51  new_values(1) = 4.0;
52  new_values(2) = 40.0;
53  new_values(3) = 400.0;
54 
55  CALL xf_change_scalar_values_timestep_float(xscalarid, edittimestep, editnumvalues, &
56  indices, new_values, error)
57  if (error < 0) then
58  CALL xf_close_group(xscalarid, error)
59  CALL xf_close_file(xfileid, error)
60  return
61  endif
62 
63 ! RDJ - Technically we shouldn't have to close and reopen a file that we are
64 ! editing but it seems to fix a crash that the tests are having.
65 ! CALL XF_CLOSE_GROUP(xScalarId, error)
66 ! CALL XF_CLOSE_FILE(xFileId, error)
67 ! ! open the file and edit the values
68 ! CALL XF_OPEN_FILE(a_Filename, .FALSE., xFileId, error)
69 ! if (error < 0) then
70 ! return
71 ! endif
72 !
73 ! CALL XF_OPEN_GROUP(xFileId, SCALAR_A_LOCATION_FULL, xScalarId, error);
74 ! if (error < 0) then
75 ! CALL XF_CLOSE_FILE(xFileId, error)
76 ! return
77 ! endif
78 
79  ! Edit values in timestep 2, make index 2 = 6, index 3 = 60, and
80  ! index 9 = 600
81  edittimestep = 2
82  indices(1) = 2
83  indices(2) = 3
84  indices(3) = 9
85  new_values(1) = -6.0
86  new_values(2) = 60.0
87  new_values(3) = 6000.0
88 
89  CALL xf_change_scalar_values_timestep_float(xscalarid, edittimestep, editnumvalues, &
90  indices, new_values, error)
91  if (error < 0) then
92  CALL xf_close_group(xscalarid, error)
93  CALL xf_close_file(xfileid, error)
94  return
95  endif
96 
97  CALL xf_close_group(xscalarid, error)
98  CALL xf_close_file(xfileid, error)
99 
100  return
101 END SUBROUTINE
102 ! tdEditScalarAValues
103 
104 ! ---------------------------------------------------------------------------
105 ! FUNCTION tdReadScalarAIndices
106 ! PURPOSE
107 ! NOTES
108 ! ---------------------------------------------------------------------------
109 SUBROUTINE td_read_scalar_a_indices (a_Filename, a_nIndices, a_indices, error)
110  CHARACTER(LEN=*), INTENT(IN) :: a_filename
111  INTEGER, INTENT(IN) :: a_nindices
112  INTEGER, DIMENSION(*), INTENT(IN) :: a_indices
113  INTEGER, INTENT(OUT) :: error
114  INTEGER xfileid, xdsetsid, xscalaraid
115  INTEGER ntimesteps
116  REAL, ALLOCATABLE, DIMENSION(:) :: fvalues
117  INTEGER nvalues
118  INTEGER id, i, j
119 
120  ! open the file
121  CALL xf_open_file(a_filename, .true., xfileid, error)
122  if (error < 0) then
123  return
124  endif
125 
126  ! open the dataset group
127  CALL xf_open_group(xfileid, datasets_location, xdsetsid, error)
128  if (error >= 0) then
129  CALL xf_open_group(xdsetsid, scalar_a_location, xscalaraid, error)
130  endif
131  if (error < 0) then
132  return
133  endif
134 
135  ! Find out the number of timesteps in the file
136  CALL xf_get_dataset_num_times(xscalaraid, ntimesteps, error)
137  if (error < 0) then
138  return
139  endif
140  if (ntimesteps < 1) then
141  error = -1
142  return
143  endif
144 
145  ! Read the values for the index
146  nvalues = ntimesteps*a_nindices
147  allocate(fvalues(nvalues))
148  CALL xf_read_scalar_values_at_indices_float(xscalaraid, a_nindices, a_indices, 1, &
149  ntimesteps, fvalues, error)
150  if (error < 0) then
151  return
152  endif
153 
154  ! output the data
155  WRITE(*,*) ''
156  WRITE(*,*) 'Reading scalar A indices'
157  id = 1;
158  do i = 1, ntimesteps
159  WRITE(*,*) 'Timestep: ', i
160  do j = 1, a_nindices
161  WRITE(*,*) 'index:', a_indices(j), ' value: ', fvalues(id)
162  id = id + 1
163  enddo
164  enddo
165  WRITE(*,*) ''
166 
167  deallocate(fvalues)
168 
169  return
170 
171 END SUBROUTINE
172 ! TD_READ_SCALARA_INDICES
173 
174 ! --------------------------------------------------------------------------
175 ! FUNCTION tdReadDatasets
176 ! PURPOSE Read a dataset group from an XMDF file and output information to
177 ! to a text file
178 ! NOTES
179 ! --------------------------------------------------------------------------
180 RECURSIVE SUBROUTINE td_read_datasets (a_xGroupId, a_FileUnit, error)
181 INTEGER, INTENT(IN) :: a_xgroupid
182 INTEGER, INTENT(IN) :: a_fileunit
183 INTEGER, INTENT(OUT) :: error
184 INTEGER npaths, nmaxpathlength, j
185 CHARACTER, ALLOCATABLE, DIMENSION(:) :: paths
186 CHARACTER(LEN=500) individualpath
187 INTEGER nstatus, i
188 INTEGER xscalarid, xvectorid, xmultiid
189 INTEGER nmultidatasets
190 
191 xscalarid = none
192 xvectorid = none
193 
194 nmultidatasets = 0
195 npaths = 0
196 nmaxpathlength = 0
197 
198  ! Look for scalar datasets
199 call xf_get_scalar_datasets_info(a_xgroupid, npaths, nmaxpathlength, nstatus)
200 if (nstatus >= 0 .AND. npaths > 0) then
201  allocate(paths(npaths*nmaxpathlength))
202  call xf_get_scalar_dataset_paths(a_xgroupid, npaths, nmaxpathlength, paths, &
203  error)
204 endif
205 if (nstatus < 0) then
206  error = -1
207  return
208 endif
209 
210  ! Output number and paths to scalar datasets
211 WRITE(a_fileunit,*) 'Number of Scalars ', npaths
212 do i=1, npaths
213  individualpath = ''
214  do j=1, nmaxpathlength-1
215  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
216  enddo
217  WRITE(a_fileunit,*) 'Reading scalar: ', individualpath(1:nmaxpathlength-1)
218  call xf_open_group(a_xgroupid, individualpath(1:nmaxpathlength-1), &
219  xscalarid, nstatus)
220  if (nstatus < 0) then
221  error = -1
222  return
223  endif
224 
225  call tdi_read_scalar(xscalarid, a_fileunit, nstatus)
226  call xf_close_group(xscalarid, error)
227  if (nstatus < 0) then
228  WRITE(*,*) 'Error reading scalar dataset.'
229  error = -1
230  return
231  endif
232 enddo
233 
234 if (allocated(paths)) deallocate(paths)
235  ! Look for vector datasets
236 call xf_get_vector_datasets_info(a_xgroupid, npaths, nmaxpathlength, nstatus)
237 if (nstatus >= 0 .AND. npaths > 0) then
238  allocate(paths(npaths*nmaxpathlength))
239  call xf_get_vector_dataset_paths(a_xgroupid, npaths, nmaxpathlength, paths, error)
240 endif
241 if (nstatus < 0) then
242  error = -1
243  return
244 endif
245 
246  ! Output number and paths to scalar datasets
247 WRITE(a_fileunit,*) 'Number of Vectors ', npaths
248 do i=2, npaths
249  do j=1, nmaxpathlength-1
250  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
251  enddo
252  WRITE(a_fileunit,*) 'Reading Vector: ', &
253  individualpath(1:nmaxpathlength-1)
254  call xf_open_group(a_xgroupid, individualpath(1:nmaxpathlength-1), &
255  xvectorid, nstatus)
256  if (nstatus < 0) then
257  error = -1
258  return
259  endif
260  call tdi_read_vector(xvectorid, a_fileunit, nstatus)
261  call xf_close_group(xvectorid, error)
262  if (nstatus < 0) then
263  WRITE(*,*) 'Error reading vector dataset.'
264  error = -1
265  return
266  endif
267 enddo
268 
269 if (allocated(paths)) deallocate(paths)
270 
271 ! find multidataset folders
272 call xf_get_grp_pths_sz_mlt_dsets(a_xgroupid, nmultidatasets, &
273  nmaxpathlength, nstatus)
274 if (nstatus >= 0 .AND. nmultidatasets > 0) then
275  allocate(paths(nmultidatasets*nmaxpathlength))
276  call xf_get_all_grp_paths_mlt_dsets(a_xgroupid, nmultidatasets, &
277  nmaxpathlength, paths, error)
278  if (nstatus < 0) then
279  error = -1
280  return
281  endif
282 
283  ! Output number and paths to multidatasets
284  WRITE(a_fileunit,*) 'Number of Multidatasets ', nmultidatasets
285  do i=2, nmultidatasets
286  individualpath = ''
287  do j=1, nmaxpathlength-1
288  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
289  enddo
290  WRITE(a_fileunit,*) 'Reading multidataset: ', &
291  individualpath(1:nmaxpathlength-1)
292  call xf_open_group(a_xgroupid, individualpath(1:nmaxpathlength-1), &
293  xmultiid, nstatus)
294  if (nstatus < 0) then
295  error = -1
296  return
297  endif
298 
299  call td_read_datasets(xmultiid, a_fileunit, nstatus)
300  call xf_close_group(xmultiid, error)
301  if (nstatus < 0) then
302  WRITE(*,*) 'Error reading multidatasets.'
303  error = -1
304  return
305  endif
306  enddo
307 endif
308 if (allocated(paths)) deallocate(paths)
309 
310 error = 1
311 return
312 
313 END SUBROUTINE
314 !tdReadDatasets
315 ! --------------------------------------------------------------------------
316 ! FUNCTION tdReadActivityScalarAIndex
317 ! PURPOSE Read all timestep values for a particular index
318 ! NOTES
319 ! --------------------------------------------------------------------------
320 SUBROUTINE td_read_activity_scalar_a_index(a_Filename, a_Index, error)
321 CHARACTER(LEN=*), INTENT(IN) :: a_filename
322 INTEGER, INTENT(IN) :: a_index
323 INTEGER, INTENT(OUT) :: error
324 INTEGER status
325 INTEGER xfileid, xdsetsid, xscalaraid
326 INTEGER ntimesteps, i
327 INTEGER, ALLOCATABLE :: bactive(:)
328 
329 xfileid = none
330 xdsetsid = none
331 xscalaraid = none
332 
333  ! open the file
334 call xf_open_file(a_filename, .true., xfileid, status)
335 if (status < 0) then
336  error = -1
337  return
338 endif
339 
340  ! open the dataset group
341 call xf_open_group(xfileid, datasets_location, xdsetsid, status)
342 if (status >= 0) then
343  call xf_open_group(xdsetsid, scalar_a_location, xscalaraid, status)
344 endif
345 if (status < 0) then
346  error = status
347  return
348 endif
349 
350  ! Find out the number of timesteps in the file
351 CALL xf_get_dataset_num_times(xscalaraid, ntimesteps, status)
352 if (status < 0) then
353  error = status
354  return
355 endif
356 
357 if (ntimesteps < 1) then
358  error = -1
359  return
360 endif
361 
362  ! Read the values for the index
363 allocate(bactive(ntimesteps))
364 call xf_read_active_vals_at_index(xscalaraid, a_index, 1, ntimesteps, &
365  bactive, status)
366  ! output the data
367 WRITE(*,*) ''
368 WRITE(*,*) 'Reading activity for scalar A slice at index: ', a_index
369 do i=1, ntimesteps
370  WRITE(*,*) bactive(i), ' '
371 enddo
372 
373 deallocate(bactive)
374 
375 error = status
376 return
377 
378 END SUBROUTINE
379 ! tdReadActivityScalarAIndex
380 
381 ! --------------------------------------------------------------------------
382 ! FUNCTION tdReadScalarAIndex
383 ! PURPOSE Read all timestep values for a particular index
384 ! NOTES
385 ! --------------------------------------------------------------------------
386 SUBROUTINE td_read_scalar_a_index (a_Filename, a_Index, error)
387 CHARACTER(LEN=*), INTENT(IN) :: a_filename
388 INTEGER, INTENT(IN) :: a_index
389 INTEGER, INTENT(OUT) :: error
390 INTEGER status
391 INTEGER xfileid, xdsetsid, xscalaraid
392 INTEGER ntimesteps, i
393 REAL, ALLOCATABLE :: fvalues(:)
394 
395 xfileid = none
396 xdsetsid = none
397 xscalaraid = none
398 
399  ! open the file
400 call xf_open_file(a_filename, .true., xfileid, status)
401 if (status < 0) then
402  error = -1
403  return
404 endif
405 
406  ! open the dataset group
407 call xf_open_group(xfileid, datasets_location, xdsetsid, status)
408 if (status >= 0) then
409  call xf_open_group(xdsetsid, scalar_a_location, xscalaraid, status)
410 endif
411 if (status < 0) then
412  error = status
413  return
414 endif
415 
416  ! Find out the number of timesteps in the file
417 call xf_get_dataset_num_times(xscalaraid, ntimesteps, status)
418 if (status < 0) then
419  error = status
420  return
421 endif
422 
423 if (ntimesteps < 1) then
424  error = -1
425  return
426 endif
427 
428  ! Read the values for the index
429 allocate (fvalues(ntimesteps))
430 call xf_read_scalar_values_at_index(xscalaraid, a_index, 1, ntimesteps, &
431  fvalues, status)
432 
433  ! output the data
434 WRITE(*,*) ''
435 WRITE(*,*) 'Reading scalar A slice at index: ', a_index
436 do i=1, ntimesteps
437  WRITE(*,*) fvalues(i), ' '
438 enddo
439 
440 deallocate(fvalues)
441 
442 error = status
443 return
444 
445 END SUBROUTINE
446 ! tdReadScalarAtIndex
447 
448 ! --------------------------------------------------------------------------
449 ! FUNCTION tdWriteScalarA
450 ! PURPOSE Write scalar Dataset to an HDF5 File
451 ! NOTES This tests dynamic data sets, and activity
452 ! This dataset is dynamic concentrations (mg/L) with output times
453 ! in minutes.
454 ! Dataset is for a mesh and so nActive is the number of elements
455 ! which is not the same as the nValues which would be number of nodes
456 ! reads/writes a reference time in julian days
457 ! --------------------------------------------------------------------------
458 SUBROUTINE td_write_scalar_a (a_Filename, a_Compression, error)
459  CHARACTER(LEN=*), INTENT(IN) :: a_filename
460  INTEGER, INTENT(IN) :: a_compression
461  INTEGER, INTENT(OUT) :: error
462  INTEGER xfileid, xdsetsid, xscalaraid, xcoordid
463  INTEGER nvalues, ntimes, nactive
464  REAL(DOUBLE) dtime, djulianreftime
465  INTEGER itimestep, iactive, ihpgnzone
466  REAL fvalues(10) ! nValues
467  INTEGER*1 bactivity(10) ! activity
468  INTEGER i, status
469 
470  ! initialize the data
471  nvalues = 10
472  ntimes = 3
473  nactive = 8
474  dtime = 0.0
475 
476  ! 5th item in data set is always inactive, others active
477  do iactive = 1, nactive
478  bactivity(iactive) = 1
479  enddo
480  bactivity(6) = 0
481 
482 
483  ! create the file
484  call xf_create_file(a_filename, .true., xfileid, error)
485  if (error .LT. 0) then
486  ! close the file
487  call xf_close_file(xfileid, error)
488  return
489  endif
490 
491  ! create the group where we will put all the datasets
492  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
493  if (status < 0) then
494  call xf_close_file(xfileid, error)
495  error = -1
496  return
497  endif
498 
499  ! Create the scalar A dataset group
500  call xf_create_scalar_dataset(xdsetsid, scalar_a_location, 'mg/L', &
501  ts_hours, a_compression, xscalaraid, status)
502  if (status .LT. 0) then
503  ! close the dataset
504  call xf_close_group(xscalaraid, error)
505  call xf_close_group(xdsetsid, error)
506  call xf_close_file(xfileid, error)
507  error = status
508  return
509  endif
510 
511  ! Add in a reftime. This is a julian day for:
512  ! noon July 1, 2003
513  djulianreftime = 2452822.0;
514  call xf_write_reftime(xscalaraid, djulianreftime, status)
515  if (status < 0) then
516  call xf_close_group(xscalaraid, error)
517  call xf_close_group(xdsetsid, error)
518  call xf_close_file(xfileid, error)
519  endif
520 
521  ! Loop through timesteps adding them to the file
522  do itimestep = 1, ntimes
523  ! We will have an 0.5 hour timestep
524  dtime = itimestep * 0.5
525 
526  fvalues(1) = dtime
527  do i = 2, nvalues
528  fvalues(i) = fvalues(i-1)*2.5
529  end do
530 
531  ! write the dataset array values
532  call xf_write_scalar_timestep(xscalaraid, dtime, nvalues, fvalues, error)
533  if (error .GE. 0) then
534  ! write activity array
535  call xf_write_activity_timestep(xscalaraid, nactive, bactivity, error)
536  end if
537  enddo
538 
539  ! Write Coordinate file - for ScalarA, we will set the coordinate system
540  ! to be Geographic HPGN, with HPGN settings written to the file.
541  call xf_create_coordinate_group(xfileid, xcoordid, status)
542  if (status < 0) then
543  call xf_close_group(xscalaraid, error)
544  call xf_close_group(xdsetsid, error)
545  call xf_close_file(xfileid, error)
546  error = -1
547  return
548  endif
549 
550  ! set HPGN Zone for test
551  ihpgnzone = 29 ! Utah
552  ! Write Coordinate Information to file
553  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic_hpgn, error)
554  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
555  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
556  call xf_set_vert_units(xcoordid, coord_units_meters, error)
557 
558  ! write additional information
559  call xf_set_hpgn_area(xcoordid, ihpgnzone, error)
560 
561  call xf_close_group(xcoordid, error)
562  xcoordid = 0;
563 
564  ! close the dataset
565  call xf_close_group(xscalaraid, error)
566  call xf_close_group(xdsetsid, error)
567  call xf_close_file(xfileid, error)
568 
569  return
570 END SUBROUTINE
571 ! tdWriteScalarA
572 
573 ! --------------------------------------------------------------------------
574 ! FUNCTION tdWriteScalarAPIECES
575 ! PURPOSE Write scalar Dataset to an HDF5 File
576 ! NOTES This tests dynamic data sets, and activity
577 ! This dataset is dynamic concentrations (mg/L) with output times
578 ! in minutes.
579 ! Dataset is for a mesh and so nActive is the number of elements
580 ! which is not the same as the nValues which would be number of nodes
581 ! reads/writes a reference time in julian days
582 ! --------------------------------------------------------------------------
583 SUBROUTINE td_write_scalar_a_pieces (a_Filename, a_Compression, error)
584  CHARACTER(LEN=*), INTENT(IN) :: a_filename
585  INTEGER, INTENT(IN) :: a_compression
586  INTEGER, INTENT(OUT) :: error
587  INTEGER xfileid, xdsetsid, xscalaraid, xcoordid
588  INTEGER nvalues, ntimes, nactive
589  REAL(DOUBLE) dtime, djulianreftime
590  INTEGER itimestep, iactive, ihpgnzone
591  REAL fvalues(10) ! nValues
592  INTEGER*1 bactivity(10) ! activity
593  INTEGER i, status
594  REAL minvalue, maxvalue
595  INTEGER timestepid
596  INTEGER activets
597 
598  ! initialize the data
599  nvalues = 10
600  ntimes = 3
601  nactive = 8
602  dtime = 0.0
603 
604  ! 5th item in data set is always inactive, others active
605  do iactive = 1, nactive
606  bactivity(iactive) = 1
607  enddo
608  bactivity(6) = 0
609 
610 
611  ! create the file
612  call xf_create_file(a_filename, .true., xfileid, error)
613  if (error .LT. 0) then
614  ! close the file
615  call xf_close_file(xfileid, error)
616  return
617  endif
618 
619  ! create the group where we will put all the datasets
620  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
621  if (status < 0) then
622  call xf_close_file(xfileid, error)
623  error = -1
624  return
625  endif
626 
627  ! Create the scalar A dataset group
628  call xf_create_scalar_dataset(xdsetsid, scalar_a_location, 'mg/L', &
629  ts_hours, a_compression, xscalaraid, status)
630  if (status .LT. 0) then
631  ! close the dataset
632  call xf_close_group(xscalaraid, error)
633  call xf_close_group(xdsetsid, error)
634  call xf_close_file(xfileid, error)
635  error = status
636  return
637  endif
638 
639  ! Add in a reftime. This is a julian day for:
640  ! noon July 1, 2003
641  djulianreftime = 2452822.0;
642  call xf_write_reftime(xscalaraid, djulianreftime, status)
643  if (status < 0) then
644  call xf_close_group(xscalaraid, error)
645  call xf_close_group(xdsetsid, error)
646  call xf_close_file(xfileid, error)
647  endif
648 
649  ! Loop through timesteps adding them to the file
650  do itimestep = 1, ntimes
651  ! We will have an 0.5 hour timestep
652  dtime = itimestep * 0.5
653 
654  fvalues(1) = dtime
655  minvalue = fvalues(1)
656  maxvalue = fvalues(1)
657  do i = 2, nvalues
658  fvalues(i) = fvalues(i-1)*2.5
659  minvalue = min(minvalue, fvalues(i))
660  maxvalue = max(maxvalue, fvalues(i))
661  end do
662 
663  ! write the dataset array values
664  call xf_initialize_scalar_timestep(xscalaraid, dtime, nvalues, minvalue, &
665  maxvalue, timestepid, error)
666 
667  ! write data in pairs
668  do i = 1, nvalues, +2
669  call xf_write_scalar_timestep_portion(xscalaraid, timestepid, 2, i, &
670  fvalues(i), error)
671  enddo
672 
673  if (error .GE. 0) then
674  ! write activity array
675  call xf_initialize_activity_timestep(xscalaraid, nactive, activets, error)
676 
677  do i = 1, nactive, +2
678  call xf_write_activity_timestep_portion(xscalaraid, activets, 2, &
679  i, bactivity(i), error)
680  enddo
681 
682  end if
683  enddo
684 
685  ! Write Coordinate file - for ScalarA, we will set the coordinate system
686  ! to be Geographic HPGN, with HPGN settings written to the file.
687  call xf_create_coordinate_group(xfileid, xcoordid, status)
688  if (status < 0) then
689  call xf_close_group(xscalaraid, error)
690  call xf_close_group(xdsetsid, error)
691  call xf_close_file(xfileid, error)
692  error = -1
693  return
694  endif
695 
696  ! set HPGN Zone for test
697  ihpgnzone = 29 ! Utah
698  ! Write Coordinate Information to file
699  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic_hpgn, error)
700  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
701  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
702  call xf_set_vert_units(xcoordid, coord_units_meters, error)
703 
704  ! write additional information
705  call xf_set_hpgn_area(xcoordid, ihpgnzone, error)
706 
707  call xf_close_group(xcoordid, error)
708  xcoordid = 0;
709 
710  ! close the dataset
711  call xf_close_group(xscalaraid, error)
712  call xf_close_group(xdsetsid, error)
713  call xf_close_file(xfileid, error)
714 
715  return
716 END SUBROUTINE
717 ! tdWriteScalarAPieces
718 
719 ! --------------------------------------------------------------------------
720 ! FUNCTION tdWriteScalarAPIECESAltMinMax
721 ! PURPOSE Write scalar Dataset to an HDF5 File
722 ! NOTES This tests dynamic data sets, and activity
723 ! This dataset is dynamic concentrations (mg/L) with output times
724 ! in minutes.
725 ! Dataset is for a mesh and so nActive is the number of elements
726 ! which is not the same as the nValues which would be number of nodes
727 ! reads/writes a reference time in julian days
728 ! --------------------------------------------------------------------------
729 SUBROUTINE td_write_scalar_a_pieces_alt_min_max (a_Filename, a_Compression, error)
730  CHARACTER(LEN=*), INTENT(IN) :: a_filename
731  INTEGER, INTENT(IN) :: a_compression
732  INTEGER, INTENT(OUT) :: error
733  INTEGER xfileid, xdsetsid, xscalaraid, xcoordid
734  INTEGER nvalues, ntimes, nactive
735  REAL(DOUBLE) dtime, djulianreftime
736  INTEGER itimestep, iactive, ihpgnzone
737  REAL fvalues(10) ! nValues
738  INTEGER*1 bactivity(10) ! activity
739  INTEGER i, status
740  REAL minvalue, maxvalue
741  INTEGER timestepid
742  INTEGER activets
743 
744  ! initialize the data
745  nvalues = 10
746  ntimes = 3
747  nactive = 8
748  dtime = 0.0
749 
750  ! 5th item in data set is always inactive, others active
751  do iactive = 1, nactive
752  bactivity(iactive) = 1
753  enddo
754  bactivity(6) = 0
755 
756 
757  ! create the file
758  call xf_create_file(a_filename, .true., xfileid, error)
759  if (error .LT. 0) then
760  ! close the file
761  call xf_close_file(xfileid, error)
762  return
763  endif
764 
765  ! create the group where we will put all the datasets
766  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
767  if (status < 0) then
768  call xf_close_file(xfileid, error)
769  error = -1
770  return
771  endif
772 
773  ! Create the scalar A dataset group
774  call xf_create_scalar_dataset(xdsetsid, scalar_a_location, 'mg/L', &
775  ts_hours, a_compression, xscalaraid, status)
776  if (status .LT. 0) then
777  ! close the dataset
778  call xf_close_group(xscalaraid, error)
779  call xf_close_group(xdsetsid, error)
780  call xf_close_file(xfileid, error)
781  error = status
782  return
783  endif
784 
785  ! Add in a reftime. This is a julian day for:
786  ! noon July 1, 2003
787  djulianreftime = 2452822.0;
788  call xf_write_reftime(xscalaraid, djulianreftime, status)
789  if (status < 0) then
790  call xf_close_group(xscalaraid, error)
791  call xf_close_group(xdsetsid, error)
792  call xf_close_file(xfileid, error)
793  endif
794 
795  ! Loop through timesteps adding them to the file
796  do itimestep = 1, ntimes
797  ! We will have an 0.5 hour timestep
798  dtime = itimestep * 0.5
799 
800  fvalues(1) = dtime
801  minvalue = fvalues(1)
802  maxvalue = fvalues(1)
803  do i = 2, nvalues
804  fvalues(i) = fvalues(i-1)*2.5
805  minvalue = min(minvalue, fvalues(i))
806  maxvalue = max(maxvalue, fvalues(i))
807  end do
808 
809  ! write the dataset array values
810  call xf_initialize_scalar_timestep(xscalaraid, dtime, nvalues, minvalue, &
811  maxvalue, timestepid, error)
812 
813  ! write data in pairs
814  do i = 1, nvalues, +2
815  call xf_write_scalar_timestep_portion(xscalaraid, timestepid, 2, i, &
816  fvalues(i), error)
817  enddo
818 
819  minvalue = 0.1111*timestepid
820  maxvalue = 1111*timestepid
821  call xf_set_dataset_timestep_min_max(xscalaraid, timestepid, minvalue, &
822  maxvalue, error)
823 
824  if (error .GE. 0) then
825  ! write activity array
826  call xf_initialize_activity_timestep(xscalaraid, nactive, activets, error)
827 
828  do i = 1, nactive, +2
829  call xf_write_activity_timestep_portion(xscalaraid, activets, 2, &
830  i, bactivity(i), error)
831  enddo
832 
833  end if
834  enddo
835 
836  ! Write Coordinate file - for ScalarA, we will set the coordinate system
837  ! to be Geographic HPGN, with HPGN settings written to the file.
838  call xf_create_coordinate_group(xfileid, xcoordid, status)
839  if (status < 0) then
840  call xf_close_group(xscalaraid, error)
841  call xf_close_group(xdsetsid, error)
842  call xf_close_file(xfileid, error)
843  error = -1
844  return
845  endif
846 
847  ! set HPGN Zone for test
848  ihpgnzone = 29 ! Utah
849  ! Write Coordinate Information to file
850  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic_hpgn, error)
851  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
852  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
853  call xf_set_vert_units(xcoordid, coord_units_meters, error)
854 
855  ! write additional information
856  call xf_set_hpgn_area(xcoordid, ihpgnzone, error)
857 
858  call xf_close_group(xcoordid, error)
859  xcoordid = 0;
860 
861  ! close the dataset
862  call xf_close_group(xscalaraid, error)
863  call xf_close_group(xdsetsid, error)
864  call xf_close_file(xfileid, error)
865 
866  return
867 END SUBROUTINE
868 ! TD_WRITE_SCALAR_A_PIECES_ALT_MIN_MAX
869 
870 
871 ! --------------------------------------------------------------------------
872 ! FUNCTION TD_WRITE_SCALAR_B
873 ! PURPOSE Write scalar Dataset to an HDF5 File
874 ! NOTES This tests dynamic data sets, and activity
875 ! This dataset is dynamic concentrations (mg/L) with output times
876 ! in minutes.
877 ! Dataset is for a mesh and so nActive is the number of elements
878 ! which is not the same as the nValues which would be number of nodes
879 ! reads/writes a reference time in julian days
880 ! --------------------------------------------------------------------------
881 SUBROUTINE td_write_scalar_b (a_Filename, a_Compression, a_Overwrite, error)
882  CHARACTER(LEN=*), INTENT(IN) :: a_filename
883  INTEGER, INTENT(IN) :: a_compression
884  LOGICAL, INTENT(IN) :: a_overwrite
885  INTEGER, INTENT(OUT) :: error
886  INTEGER xfileid, xdsetsid, xscalarbid, xcoordid
887  INTEGER nvalues, ntimes, nactive
888  REAL(DOUBLE) dtime, djulianreftime
889  INTEGER itimestep, iactive
890  REAL fvalues(10) ! nValues
891  INTEGER*1 bactivity(10) ! activity
892  INTEGER i, status
893 
894  ! initialize the data
895  nvalues = 10
896  ntimes = 3
897  nactive = 8
898  dtime = 0.0
899  i = 0
900 
901  ! 5th item in data set is always inactive, others active
902  do iactive = 1, nactive
903  bactivity(iactive) = 1
904  enddo
905  bactivity(6) = 0
906 
907  if (a_overwrite) then
908  ! open the already-existing file
909  call xf_open_file(a_filename, .false., xfileid, status)
910  if (status < 0) then
911  error = -1
912  return
913  endif
914  ! open the group where we have all the datasets
915  call xf_open_group(xfileid, datasets_location, xdsetsid, status)
916  if (status < 0) then
917  call xf_close_file(xfileid, error)
918  error = -1
919  return
920  endif
921  else
922  ! create the file
923  call xf_create_file(a_filename, .true., xfileid, error)
924  if (error .LT. 0) then
925  ! close the file
926  call xf_close_file(xfileid, error)
927  return
928  endif
929 
930  ! create the group where we will put all the datasets
931  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
932  if (status < 0) then
933  call xf_close_file(xfileid, error)
934  error = -1
935  return
936  endif
937  endif
938 
939  ! Create/Overwrite the scalar B dataset group
940  call xf_create_scalar_dataset(xdsetsid, scalar_b_location, 'mg/L', &
941  ts_hours, a_compression, xscalarbid, status)
942  if (status < 0) then
943  ! close the dataset
944  call xf_close_group(xscalarbid, error)
945  call xf_close_group(xdsetsid, error)
946  call xf_close_file(xfileid, error)
947  error = status
948  return
949  endif
950 
951  ! Add in a reftime. This is a julian day for:
952  ! noon July 1, 2003
953  djulianreftime = 2452822.0;
954  call xf_write_reftime(xscalarbid, djulianreftime, status)
955  if (status < 0) then
956  call xf_close_group(xscalarbid, error)
957  call xf_close_group(xdsetsid, error)
958  call xf_close_file(xfileid, error)
959  endif
960 
961  if (.NOT. a_overwrite) then
962  ! Loop through timesteps adding them to the file
963  do itimestep = 1, ntimes
964  ! We will have an 0.5 hour timestep
965  dtime = itimestep * 0.5
966 
967  fvalues(1) = dtime
968  do i = 2, nvalues
969  fvalues(i) = fvalues(i-1)*2.5
970  end do
971 
972  ! write the dataset array values
973  call xf_write_scalar_timestep(xscalarbid, dtime, nvalues, fvalues, error)
974  if (error .GE. 0) then
975  ! write activity array
976  call xf_write_activity_timestep(xscalarbid, nactive, bactivity, error)
977  end if
978  if (error < 0) then
979  call xf_close_group(xscalarbid, error)
980  call xf_close_group(xdsetsid, error)
981  call xf_close_file(xfileid, error)
982  endif
983  enddo
984  else
985  ! Loop through timesteps adding them to the file
986  do itimestep = 1, ntimes
987  ! We will have an 1.5 hour timestep
988  dtime = itimestep * 1.5
989 
990  fvalues(1) = dtime
991  do i = 2, nvalues
992  fvalues(i) = fvalues(i-1)*1.5
993  end do
994 
995  ! write the dataset array values
996  call xf_write_scalar_timestep(xscalarbid, dtime, nvalues, fvalues, error)
997  if (error .GE. 0) then
998  ! write activity array
999  call xf_write_activity_timestep(xscalarbid, nactive, bactivity, error)
1000  end if
1001  if (error < 0) then
1002  call xf_close_group(xscalarbid, error)
1003  call xf_close_group(xdsetsid, error)
1004  call xf_close_file(xfileid, error)
1005  endif
1006  enddo
1007  endif
1008 
1009  if (.NOT. a_overwrite) then
1010  ! Write Coordinate file - for ScalarB, we will set the coordinate system
1011  ! to be UTM, with UTM Zone settings written to the file.
1012  call xf_create_coordinate_group(xfileid, xcoordid, status)
1013  if (status < 0) then
1014  call xf_close_group(xscalarbid, error)
1015  call xf_close_group(xdsetsid, error)
1016  call xf_close_file(xfileid, error)
1017  error = -1
1018  return
1019  endif
1020 
1021  ! Write Coord Info to file
1022  call xf_set_horiz_datum(xcoordid, horiz_datum_utm, error)
1023  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
1024 
1025  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
1026  call xf_set_vert_units(xcoordid, coord_units_meters, error)
1027 
1028  ! write additional information - we'll use the max value for this test
1029  call xf_set_utm_zone(xcoordid, utm_zone_max, error)
1030 
1031  call xf_close_group(xcoordid, error)
1032  xcoordid = 0
1033  endif
1034 
1035  ! close the dataset
1036  call xf_close_group(xscalarbid, error)
1037  call xf_close_group(xdsetsid, error)
1038  call xf_close_file(xfileid, error)
1039 
1040  error = 1
1041  return
1042 END SUBROUTINE
1043 ! tdWriteScalarB
1044 !------------------------------------------------------------------------------
1045 ! FUNCTION TD_WRITE_COORDS_TO_MULTI
1046 ! PURPOSE Write coordinate system to a multidataset file
1047 ! NOTES
1048 !------------------------------------------------------------------------------
1049 SUBROUTINE td_write_coords_to_multi (a_xFileId, error)
1050 INTEGER, INTENT(IN) :: a_xfileid
1051 INTEGER, INTENT(OUT) :: error
1052 INTEGER xcoordid
1053 INTEGER status
1054 
1055  ! Write Coordinate file - for Multidatasets, we will set the coordinate system
1056  ! to be UTM, with UTM Zone settings written to the file.
1057  call xf_create_coordinate_group(a_xfileid, xcoordid, status)
1058  if (status < 0) then
1059  error = status
1060  return
1061  endif
1062 
1063  ! Write Coord Info to file
1064  call xf_set_horiz_datum(xcoordid, horiz_datum_utm, error)
1065  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
1066 
1067  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
1068  call xf_set_vert_units(xcoordid, coord_units_meters, error)
1069 
1070  ! write additional information - we'll use the max value for this test
1071  call xf_set_utm_zone(xcoordid, utm_zone_max, error)
1072 
1073  call xf_close_group(xcoordid, error)
1074  xcoordid = 0
1075 
1076  return
1077 END SUBROUTINE
1078 
1079 ! --------------------------------------------------------------------------
1080 ! FUNCTION tdWriteScalarAToMulti
1081 ! PURPOSE Write scalar Dataset to an HDF5 File
1082 ! NOTES This tests dynamic data sets, and activity
1083 ! This dataset is dynamic concentrations (mg/L) with output times
1084 ! in minutes.
1085 ! Dataset is for a mesh and so nActive is the number of elements
1086 ! which is not the same as the nValues which would be number of nodes
1087 ! reads/writes a reference time in julian days
1088 ! --------------------------------------------------------------------------
1089 SUBROUTINE td_write_scalar_a_to_multi (a_GroupID, status)
1090  ! CHARACTER(LEN=*), INTENT(IN) :: a_Filename
1091  ! INTEGER, INTENT(IN) :: a_Compression
1092  ! INTEGER, INTENT(OUT) :: error
1093  INTEGER xfileid, xdsetsid, xscalaraid
1094  INTEGER a_groupid
1095  INTEGER nvalues, ntimes, nactive
1096  REAL(DOUBLE) dtime, djulianreftime
1097  INTEGER itimestep, iactive
1098  REAL fvalues(10) ! nValues
1099  INTEGER*1 bactivity(10) ! activity
1100  INTEGER i, status
1101 
1102  ! initialize the data
1103  nvalues = 10
1104  ntimes = 3
1105  nactive = 8
1106  dtime = 0.0
1107 
1108  ! 5th item in data set is always inactive, others active
1109  do iactive = 1, nactive
1110  bactivity(iactive) = 1
1111  enddo
1112  bactivity(6) = 0
1113 
1114  ! Create the scalar A dataset group
1115  call xf_create_scalar_dataset(a_groupid, scalar_a_location, 'mg/L', &
1116  ts_hours, none, xscalaraid, status)
1117  if (status .LT. 0) then
1118  ! close the dataset
1119  call xf_close_group(xscalaraid, status)
1120  call xf_close_group(xdsetsid, status)
1121  call xf_close_file(xfileid, status)
1122  return
1123  endif
1124 
1125  ! Add in a reftime. This is a julian day for:
1126  ! noon July 1, 2003
1127  djulianreftime = 2452822.0;
1128  call xf_write_reftime(xscalaraid, djulianreftime, status)
1129  if (status < 0) then
1130  call xf_close_group(xscalaraid, status)
1131  call xf_close_group(xdsetsid, status)
1132  call xf_close_file(xfileid, status)
1133  endif
1134 
1135  ! Loop through timesteps adding them to the file
1136  do itimestep = 1, ntimes
1137  ! We will have an 0.5 hour timestep
1138  dtime = itimestep * 0.5
1139 
1140  fvalues(1) = dtime
1141  do i = 2, nvalues
1142  fvalues(i) = fvalues(i-1)*2.5
1143  end do
1144 
1145  ! write the dataset array values
1146  call xf_write_scalar_timestep(xscalaraid, dtime, nvalues, fvalues, status)
1147  if (status .GE. 0) then
1148  ! write activity array
1149  call xf_write_activity_timestep(xscalaraid, nactive, bactivity, status)
1150  end if
1151  enddo
1152 
1153  ! close the dataset
1154  call xf_close_group(xscalaraid, status)
1155  !call XF_CLOSE_GROUP(a_GroupID, status)
1156  !call XF_CLOSE_FILE(a_FileID, status)
1157 
1158  return
1159 END SUBROUTINE
1160 ! tdWriteScalarAToMulti
1161 ! --------------------------------------------------------------------------
1162 ! FUNCTION tdReadVector2DAIndex
1163 ! PURPOSE Read all timestep values for a particular index
1164 ! NOTES
1165 ! --------------------------------------------------------------------------
1166 SUBROUTINE td_read_vector2d_a_index (a_Filename, a_Index, error)
1167 CHARACTER(LEN=*), INTENT(IN) :: a_filename
1168 INTEGER, INTENT(IN) :: a_index
1169 INTEGER, INTENT(OUT) :: error
1170 INTEGER status
1171 INTEGER xfileid, xdsetsid, xvector2da
1172 INTEGER ntimesteps, i
1173 REAL, ALLOCATABLE :: fvalues(:)
1174 
1175 xfileid = none
1176 xdsetsid = none
1177 xvector2da = none
1178 
1179  ! open the file
1180 call xf_open_file(a_filename, .true., xfileid, status)
1181 if (status < 0) then
1182  error = -1
1183  return
1184 endif
1185 
1186  ! open the dataset group
1187 call xf_open_group(xfileid, datasets_location, xdsetsid, status)
1188 if (status >= 0) then
1189  call xf_open_group(xdsetsid, vector2d_a_location, xvector2da, status)
1190 endif
1191 if (status < 0) then
1192  error = status
1193  return
1194 endif
1195 
1196  ! Find out the number of timesteps in the file
1197 call xf_get_dataset_num_times(xvector2da, ntimesteps, status)
1198 if (status < 0) then
1199  error = status
1200  return
1201 endif
1202 
1203 if (ntimesteps < 1) then
1204  error = -1
1205  return
1206 endif
1207 
1208  ! Read the values for the index
1209 allocate(fvalues(ntimesteps*2))
1210 call xf_read_vector_values_at_index(xvector2da, a_index, 1, ntimesteps, 2, &
1211  fvalues, status)
1212 
1213  ! output the data
1214 WRITE(*,*) ''
1215 WRITE(*,*) 'Reading vector 2D A slice at index: ', a_index
1216 do i=1, ntimesteps
1217  WRITE(*,*) fvalues(i*2-1), ' ', fvalues(i*2)
1218 enddo
1219 WRITE(*,*) ''
1220 
1221 deallocate(fvalues)
1222 
1223 error = status
1224 return
1225 
1226 END SUBROUTINE
1227 !tdReadVector2DAIndex
1228 
1229 ! --------------------------------------------------------------------------
1230 ! FUNCTION tdWriteVector2D_A
1231 ! PURPOSE Write scalar Dataset to an HDF5 File
1232 ! NOTES This tests dynamic data sets, and activity
1233 ! This dataset is dynamic concentrations (mg/L) with output times
1234 ! in minutes.
1235 ! Dataset is for a mesh and so nActive is the number of elements
1236 ! which is not the same as the nValues which would be number of nodes
1237 ! reads/writes a reference time in julian days
1238 ! --------------------------------------------------------------------------
1239 SUBROUTINE td_write_vector2d_a (a_Filename, a_Compression, error)
1240  CHARACTER(LEN=*), INTENT(IN) :: a_filename
1241  INTEGER, INTENT(IN) :: a_compression
1242  INTEGER, INTENT(OUT) :: error
1243  INTEGER xfileid, xdsetsid, xvector2d_a, xcoordid
1244  INTEGER nvalues, ntimes, ncomponents, nactive
1245  REAL(DOUBLE) dtime
1246  INTEGER itimestep, iactive
1247  REAL, DIMENSION(2, 100) :: fvalues ! nComponents, nValues
1248  INTEGER*1 bactivity(100) ! activity
1249  INTEGER i, j, status
1250  INTEGER ihpgnzone
1251 
1252  ! initialize the data
1253  ncomponents = 2
1254  nvalues = 100
1255  ntimes = 6
1256  nactive = 75
1257  dtime = 0.0
1258 
1259  ! 5th item in data set is always inactive, others active
1260  bactivity(1) = 0
1261  do iactive = 2, nactive
1262  if (mod(iactive-1, 3) == 0) then
1263  bactivity(iactive) = 0
1264  else
1265  bactivity(iactive) = 1
1266  endif
1267  enddo
1268 
1269  ! create the file
1270  call xf_create_file(a_filename, .true., xfileid, error)
1271  if (error .LT. 0) then
1272  ! close the dataset
1273  call xf_close_file(xfileid, error)
1274  return
1275  endif
1276 
1277  ! create the group where we will put all the datasets
1278  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
1279  if (status < 0) then
1280  call xf_close_file(xfileid, error)
1281  error = -1
1282  return
1283  endif
1284 
1285  ! Create the vector dataset group
1286  call xf_create_vector_dataset(xdsetsid, vector2d_a_location, 'ft/s', &
1287  ts_seconds, a_compression, xvector2d_a, status)
1288  if (status .LT. 0) then
1289  ! close the dataset
1290  call xf_close_group(xvector2d_a, error)
1291  call xf_close_group(xdsetsid, error)
1292  call xf_close_file(xfileid, error)
1293  error = status
1294  return
1295  endif
1296 
1297  ! Loop through timesteps adding them to the file
1298  do itimestep = 1, ntimes
1299  ! We will have an 0.5 hour timestep
1300  dtime = itimestep * 0.5
1301 
1302  do i = 1, nvalues
1303  do j = 1, ncomponents
1304  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
1305  end do
1306  end do
1307 
1308  ! write the dataset array values
1309  call xf_write_vector_timestep(xvector2d_a, dtime, nvalues, ncomponents, &
1310  fvalues, error)
1311  if (error .GE. 0) then
1312  ! write activity array
1313  call xf_write_activity_timestep(xvector2d_a, nactive, bactivity, error)
1314  end if
1315  enddo
1316 
1317  ! Write Coordinate file - for Vector2D_A, we will set the coordinate system
1318  ! to be Geographic HPGN, with HPGN settings written to the file.
1319  call xf_create_coordinate_group(xfileid, xcoordid, status)
1320  if (status < 0) then
1321  call xf_close_group(xvector2d_a, error)
1322  call xf_close_group(xdsetsid, error)
1323  call xf_close_file(xfileid, error)
1324  error = -1
1325  return
1326  endif
1327 
1328  ! set HPGN info for test
1329  ihpgnzone = 29 ! Utah
1330 
1331  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic_hpgn, error)
1332  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
1333  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
1334  call xf_set_vert_units(xcoordid, coord_units_meters, error)
1335 
1336  ! write additional information
1337  call xf_set_hpgn_area(xcoordid, ihpgnzone, error)
1338 
1339  call xf_close_group(xcoordid, error)
1340  xcoordid = 0
1341 
1342  ! close the dataset
1343  call xf_close_group(xvector2d_a, error)
1344  call xf_close_group(xdsetsid, error)
1345  call xf_close_file(xfileid, error)
1346 
1347  return
1348 END SUBROUTINE
1349 ! tdWriteVector2D_A
1350 
1351 ! --------------------------------------------------------------------------
1352 ! FUNCTION tdWriteVector2D_A_PIECES
1353 ! PURPOSE Write scalar Dataset to an HDF5 File
1354 ! NOTES This tests dynamic data sets, and activity
1355 ! This dataset is dynamic concentrations (mg/L) with output times
1356 ! in minutes.
1357 ! Dataset is for a mesh and so nActive is the number of elements
1358 ! which is not the same as the nValues which would be number of nodes
1359 ! reads/writes a reference time in julian days
1360 ! --------------------------------------------------------------------------
1361 SUBROUTINE td_write_vector2d_a_pieces (a_Filename, a_Compression, error)
1362  CHARACTER(LEN=*), INTENT(IN) :: a_filename
1363  INTEGER, INTENT(IN) :: a_compression
1364  INTEGER, INTENT(OUT) :: error
1365  INTEGER xfileid, xdsetsid, xvector2d_a, xcoordid
1366  INTEGER nvalues, ntimes, ncomponents, nactive
1367  REAL(DOUBLE) dtime
1368  INTEGER itimestep, iactive
1369  REAL*4, DIMENSION(2, 100) :: fvalues ! nComponents, nValues
1370  INTEGER*1 bactivity(100) ! activity
1371  INTEGER i, j, status
1372  INTEGER ihpgnzone
1373  INTEGER nvaluestowrite, ncomponentstowrite, startcomponent
1374  REAL*4 minvalue, maxvalue
1375  INTEGER timeid
1376  REAL*8 mag
1377 
1378  ! initialize the data
1379  ncomponents = 2
1380  nvalues = 100
1381  ntimes = 6
1382  nactive = 75
1383  dtime = 0.0
1384 
1385  ! 5th item in data set is always inactive, others active
1386  bactivity(1) = 0
1387  do iactive = 2, nactive
1388  if (mod(iactive-1, 3) == 0) then
1389  bactivity(iactive) = 0
1390  else
1391  bactivity(iactive) = 1
1392  endif
1393  enddo
1394 
1395  ! create the file
1396  call xf_create_file(a_filename, .true., xfileid, error)
1397  if (error .LT. 0) then
1398  ! close the dataset
1399  call xf_close_file(xfileid, error)
1400  return
1401  endif
1402 
1403  ! create the group where we will put all the datasets
1404  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
1405  if (status < 0) then
1406  call xf_close_file(xfileid, error)
1407  error = -1
1408  return
1409  endif
1410 
1411  ! Create the vector dataset group
1412  call xf_create_vector_dataset(xdsetsid, vector2d_a_location, 'ft/s', &
1413  ts_seconds, a_compression, xvector2d_a, status)
1414  if (status .LT. 0) then
1415  ! close the dataset
1416  call xf_close_group(xvector2d_a, error)
1417  call xf_close_group(xdsetsid, error)
1418  call xf_close_file(xfileid, error)
1419  error = status
1420  return
1421  endif
1422 
1423  ! Loop through timesteps adding them to the file
1424  do itimestep = 1, ntimes
1425  ! We will have an 0.5 hour timestep
1426  dtime = itimestep * 0.5
1427 
1428  do i = 1, nvalues
1429  do j = 1, ncomponents
1430  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
1431  end do
1432  end do
1433 
1434  minvalue = 99999.0
1435  maxvalue = 0.0
1436  do i = 1, nvalues
1437  mag = 0.0
1438  do j = 1, ncomponents
1439  mag = mag + fvalues(j, i)**2
1440  end do
1441  mag = mag**0.5
1442 
1443  minvalue = min(minvalue, mag)
1444  maxvalue = max(maxvalue, mag)
1445  end do
1446 
1447  ! write the dataset array values
1448  call xf_initialize_vector_timestep(xvector2d_a, dtime, nvalues, ncomponents, &
1449  minvalue, maxvalue, timeid, error)
1450 
1451  nvaluestowrite = 2
1452  ncomponentstowrite = 2
1453  startcomponent = 1
1454 
1455  do i = 1, nvalues, +2
1456  call xf_write_vector_timestep_portion(xvector2d_a, timeid, nvaluestowrite, &
1457  ncomponentstowrite, i, startcomponent, fvalues(1, i), error)
1458  enddo
1459 
1460  if (error .GE. 0) then
1461  ! write activity array
1462  call xf_write_activity_timestep(xvector2d_a, nactive, bactivity, error)
1463  end if
1464  enddo
1465 
1466  ! Write Coordinate file - for Vector2D_A, we will set the coordinate system
1467  ! to be Geographic HPGN, with HPGN settings written to the file.
1468  call xf_create_coordinate_group(xfileid, xcoordid, status)
1469  if (status < 0) then
1470  call xf_close_group(xvector2d_a, error)
1471  call xf_close_group(xdsetsid, error)
1472  call xf_close_file(xfileid, error)
1473  error = -1
1474  return
1475  endif
1476 
1477  ! set HPGN info for test
1478  ihpgnzone = 29 ! Utah
1479 
1480  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic_hpgn, error)
1481  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
1482  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
1483  call xf_set_vert_units(xcoordid, coord_units_meters, error)
1484 
1485  ! write additional information
1486  call xf_set_hpgn_area(xcoordid, ihpgnzone, error)
1487 
1488  call xf_close_group(xcoordid, error)
1489  xcoordid = 0
1490 
1491  ! close the dataset
1492  call xf_close_group(xvector2d_a, error)
1493  call xf_close_group(xdsetsid, error)
1494  call xf_close_file(xfileid, error)
1495 
1496  return
1497 END SUBROUTINE
1498 ! tdWriteVector2D_A_PIECES
1499 
1500 ! --------------------------------------------------------------------------
1501 ! FUNCTION TD_WRITE_VECTOR2D_B
1502 ! PURPOSE Write scalar Dataset to an HDF5 File
1503 ! NOTES This tests dynamic data sets, and activity
1504 ! This dataset is dynamic concentrations (mg/L) with output times
1505 ! in minutes.
1506 ! Dataset is for a mesh and so nActive is the number of elements
1507 ! which is not the same as the nValues which would be number of nodes
1508 ! reads/writes a reference time in julian days
1509 ! --------------------------------------------------------------------------
1510 SUBROUTINE td_write_vector2d_b (a_Filename, a_Compression, a_Overwrite, error)
1511  CHARACTER(LEN=*), INTENT(IN) :: a_filename
1512  INTEGER, INTENT(IN) :: a_compression
1513  LOGICAL, INTENT(IN) :: a_overwrite
1514  INTEGER, INTENT(OUT) :: error
1515  INTEGER xfileid, xdsetsid, xvector2d_b, xcoordid
1516  INTEGER nvalues, ntimes, ncomponents, nactive
1517  REAL(DOUBLE) dtime
1518  INTEGER itimestep, iactive
1519  REAL, DIMENSION(2, 100) :: fvalues
1520  INTEGER*1 bactivity(100)
1521  INTEGER i, j, status
1522 
1523  ! initialize the data
1524  ncomponents = 2
1525  nvalues = 100
1526  ntimes = 6
1527  nactive = 75
1528  dtime = 0.0
1529 
1530  ! 5th item in data set is always inactive, others active
1531  bactivity(1) = 0
1532  do iactive = 2, nactive
1533  if (mod(iactive-1, 3) == 0) then
1534  bactivity(iactive) = 0
1535  else
1536  bactivity(iactive) = 1
1537  endif
1538  enddo
1539 
1540  if (a_overwrite) then
1541  ! open the already-existing file
1542  call xf_open_file(a_filename, .false., xfileid, status)
1543  if (status < 0) then
1544  error = -1
1545  return
1546  endif
1547  ! open the group where we have all the datasets
1548  call xf_open_group(xfileid, datasets_location, xdsetsid, status)
1549  if (status < 0) then
1550  call xf_close_file(xfileid, error)
1551  error = -1
1552  return
1553  endif
1554  else
1555  ! create the file
1556  call xf_create_file(a_filename, .true., xfileid, error)
1557  if (error .LT. 0) then
1558  ! close the dataset
1559  call xf_close_file(xfileid, error)
1560  return
1561  endif
1562 
1563  ! create the group where we will put all the datasets
1564  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
1565  if (status < 0) then
1566  call xf_close_file(xfileid, error)
1567  error = -1
1568  return
1569  endif
1570  endif
1571 
1572  ! Create/Overwrite the vector dataset group
1573  call xf_create_vector_dataset(xdsetsid, vector2d_b_location, 'ft/s', &
1574  ts_seconds, a_compression, xvector2d_b, status)
1575  if (status .LT. 0) then
1576  ! close the dataset
1577  call xf_close_group(xvector2d_b, error)
1578  call xf_close_group(xdsetsid, error)
1579  call xf_close_file(xfileid, error)
1580  error = status
1581  return
1582  endif
1583 
1584  if (.NOT. a_overwrite) then
1585  ! Loop through timesteps adding them to the file
1586  do itimestep = 1, ntimes
1587  ! We will have an 0.5 hour timestep
1588  dtime = itimestep * 0.5
1589  do i = 1, nvalues
1590  do j = 1, ncomponents
1591  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
1592  end do
1593  end do
1594  ! write the dataset array values
1595  call xf_write_vector_timestep(xvector2d_b, dtime, nvalues, ncomponents, &
1596  fvalues, error)
1597  if (error .GE. 0) then
1598  ! write activity array
1599  call xf_write_activity_timestep(xvector2d_b, nactive, bactivity, error)
1600  end if
1601  if (error < 0) then
1602  call xf_close_group(xvector2d_b, error)
1603  call xf_close_group(xdsetsid, error)
1604  call xf_close_file(xfileid, error)
1605  endif
1606  enddo
1607  else
1608  ! Loop through timesteps adding them to the file
1609  do itimestep = 1, ntimes
1610  ! We will have an 1.5 hour timestep
1611  dtime = itimestep * 1.5
1612  do i = 1, nvalues
1613  do j = 1, ncomponents
1614  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
1615  end do
1616  end do
1617  ! write the dataset array values
1618  call xf_write_vector_timestep(xvector2d_b, dtime, nvalues, ncomponents, &
1619  fvalues, error)
1620  if (error .GE. 0) then
1621  ! write activity array
1622  call xf_write_activity_timestep(xvector2d_b, nactive, bactivity, error)
1623  end if
1624  if (error < 0) then
1625  call xf_close_group(xvector2d_b, error)
1626  call xf_close_group(xdsetsid, error)
1627  call xf_close_file(xfileid, error)
1628  endif
1629  enddo
1630  endif
1631 
1632  if (.NOT. a_overwrite) then
1633  ! Write Coordinate file - for ScalarB, we will set the coordinate system
1634  ! to be UTM, with UTM Zone settings written to the file.
1635  call xf_create_coordinate_group(xfileid, xcoordid, status)
1636  if (status < 0) then
1637  call xf_close_group(xvector2d_b, error)
1638  call xf_close_group(xdsetsid, error)
1639  call xf_close_file(xfileid, error)
1640  error = -1
1641  return
1642  endif
1643 
1644  ! write the coordinate data to the file
1645  call xf_set_horiz_datum(xcoordid, horiz_datum_utm, error)
1646  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
1647  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
1648  call xf_set_vert_units(xcoordid, coord_units_meters, error)
1649 
1650  ! write additional information - we'll use the max UTM zone for the test
1651  call xf_set_utm_zone(xcoordid, utm_zone_max, error)
1652 
1653  call xf_close_group(xcoordid, error)
1654  xcoordid = 0
1655  endif
1656 
1657  ! close the dataset
1658  call xf_close_group(xvector2d_b, error)
1659  call xf_close_group(xdsetsid, error)
1660  call xf_close_file(xfileid, error)
1661 
1662  return
1663 END SUBROUTINE
1664 ! tdWriteVector2D_B
1665 
1666 ! --------------------------------------------------------------------------
1667 ! FUNCTION tdWriteVector2D_AToMulti
1668 ! PURPOSE Write scalar Dataset to an HDF5 File
1669 ! NOTES This tests dynamic data sets, and activity
1670 ! This dataset is dynamic concentrations (mg/L) with output times
1671 ! in minutes.
1672 ! Dataset is for a mesh and so nActive is the number of elements
1673 ! which is not the same as the nValues which would be number of nodes
1674 ! reads/writes a reference time in julian days
1675 ! --------------------------------------------------------------------------
1676 SUBROUTINE td_write_vector2d_a_to_multi (a_FileID, a_GroupID, status)
1677  INTEGER xvector2d_a
1678  INTEGER a_fileid, a_groupid
1679  INTEGER nvalues, ntimes, ncomponents, nactive
1680  REAL(DOUBLE) dtime
1681  INTEGER itimestep, iactive
1682  REAL, DIMENSION(2, 100) :: fvalues ! nComponents, nValues
1683  INTEGER*1 bactivity(100) ! activity
1684  INTEGER i, j, status
1685 
1686  ! initialize the data
1687  ncomponents = 2
1688  nvalues = 100
1689  ntimes = 6
1690  nactive = 75
1691  dtime = 0.0
1692 
1693  ! 5th item in data set is always inactive, others active
1694  bactivity(1) = 0
1695  do iactive = 2, nactive
1696  if (mod(iactive-1, 3) == 0) then
1697  bactivity(iactive) = 0
1698  else
1699  bactivity(iactive) = 1
1700  endif
1701  enddo
1702 
1703  ! Create the vector dataset group
1704  call xf_create_vector_dataset(a_groupid, vector2d_a_location, 'ft/s', &
1705  ts_seconds, none, xvector2d_a, status)
1706  if (status .LT. 0) then
1707  ! close the dataset
1708  call xf_close_group(xvector2d_a, status)
1709  call xf_close_group(a_groupid, status)
1710  call xf_close_file(a_fileid, status)
1711  return
1712  endif
1713 
1714  ! Loop through timesteps adding them to the file
1715  do itimestep = 1, ntimes
1716  ! We will have an 0.5 hour timestep
1717  dtime = itimestep * 0.5
1718 
1719  do i = 1, nvalues
1720  do j = 1, ncomponents
1721  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
1722  end do
1723  end do
1724 
1725  ! write the dataset array values
1726  call xf_write_vector_timestep(xvector2d_a, dtime, nvalues, ncomponents, &
1727  fvalues, status)
1728  if (status .GE. 0) then
1729  ! write activity array
1730  call xf_write_activity_timestep(xvector2d_a, nactive, bactivity, status)
1731  end if
1732  enddo
1733 
1734  ! close the dataset
1735  call xf_close_group(xvector2d_a, status)
1736  return
1737 END SUBROUTINE
1738 ! tdWriteVector2D_AToMulti
1739  integer function strlen(st)
1740  integer i
1741  character st*(*)
1742  i = len(st)
1743  do while (st(i:i) .eq. ' ' .or. ichar(st(i:i)) .eq. 0)
1744  i = i - 1
1745  enddo
1746  strlen = i
1747  return
1748  end function
1749 ! --------------------------------------------------------------------------
1750 ! FUNCTION tdiReadScalar
1751 ! PURPOSE Read a scalar from an XMDF file and output information to
1752 ! to a text file
1753 ! NOTES
1754 ! --------------------------------------------------------------------------
1755 SUBROUTINE tdi_read_scalar (a_xScalarId, FileUnit, error)
1756  INTEGER, INTENT(IN) :: a_xscalarid
1757  INTEGER, INTENT(IN) :: fileunit
1758  INTEGER, INTENT(OUT) :: error
1759  INTEGER ntimes, nvalues, nactive
1760  LOGICAL*2 busereftime
1761  INTEGER itime
1762  CHARACTER(LEN=100) timeunits, units
1763  REAL(DOUBLE), ALLOCATABLE :: times(:)
1764  REAL, ALLOCATABLE :: values(:), minimums(:), maximums(:)
1765  INTEGER, ALLOCATABLE :: active(:)
1766  REAL(DOUBLE) reftime
1767 ntimes = none
1768 nvalues = none
1769 nactive = none
1770 
1771  ! read the time units
1772  call xf_get_dataset_time_units(a_xscalarid, timeunits, error)
1773  if (error < 0) return
1774 
1775  WRITE(fileunit,*) 'Time units: ', timeunits(1:strlen(timeunits))
1776 
1777  ! see if we are using a reftime
1778  call xf_use_reftime(a_xscalarid, busereftime, error)
1779  if (error < 0 ) then
1780  return
1781  endif
1782  if (busereftime) then
1783  call xf_read_reftime(a_xscalarid, reftime, error)
1784  if (error < 0 ) then
1785  return
1786  endif
1787  error = 1;
1788  !WRITE(FileUnit,*) 'Reftime: ', Reftime
1789  endif
1790 
1791  ! read the units
1792  call xf_get_dataset_units(a_xscalarid, units, error)
1793  if (error < 0) return
1794  WRITE(fileunit,*) 'Units: ', units(1:strlen(units))
1795 
1796  ! read in the number of values and number of active values
1797  call xf_get_dataset_numvals(a_xscalarid, nvalues, error)
1798  if (error .GE. 0) then
1799  call xf_get_dataset_numactive(a_xscalarid, nactive, error)
1800  endif
1801  if (error .LT. 0) return
1802 
1803  if (nvalues <= 0) then
1804  WRITE(fileunit, *) 'No data to read in.'
1805  error = -1
1806  return
1807  endif
1808 
1809  ! read in the number of times
1810  call xf_get_dataset_num_times(a_xscalarid, ntimes, error)
1811  if (error < 0) then
1812  return
1813  endif
1814 
1815  ! Read in the individual time values
1816  allocate(times(ntimes))
1817 
1818  call xf_get_dataset_times(a_xscalarid, ntimes, times, error)
1819  if (error < 0) return
1820 
1821  ! Read in the minimum and maximum values
1822  allocate(minimums(ntimes))
1823  allocate(maximums(ntimes))
1824 
1825  call xf_get_dataset_mins(a_xscalarid, ntimes, minimums, error)
1826  if (error >= 0) then
1827  call xf_get_dataset_maxs(a_xscalarid, ntimes, maximums, error)
1828  endif
1829  if (error < 0) then
1830  deallocate(times)
1831  deallocate(minimums)
1832  deallocate(maximums)
1833  return
1834  endif
1835 
1836  allocate(values(nvalues))
1837  if (nactive .GT. 0) then
1838  allocate(active(nactive))
1839  endif
1840 
1841  WRITE(fileunit,*) 'Number Timesteps: ', ntimes
1842  WRITE(fileunit,*) 'Number Values: ', nvalues
1843  WRITE(fileunit,*) 'Number Active: ', nactive
1844  WRITE(fileunit,*) ''
1845 
1846  ! loop through the timesteps, read the values and active values and write
1847  ! them to the text file
1848  do itime = 1, ntimes
1849  call xf_read_scalar_values_timestep(a_xscalarid, itime, nvalues, values, error)
1850  if (error >= 0 .AND. nactive > 0) then
1851  call xf_read_activity_timestep(a_xscalarid, itime, nactive, active, error)
1852  endif
1853 
1854  ! Write the time, min, max, values and active values to the text output
1855  ! file.
1856  WRITE(fileunit,*) 'Timestep at ', times(itime)
1857  WRITE(fileunit,*) 'Min: ', minimums(itime)
1858  WRITE(fileunit,*) 'Max: ', maximums(itime)
1859 
1860  WRITE(fileunit,*) 'Values:'
1861  WRITE(fileunit,*) values(1:nvalues)
1862  WRITE(fileunit,*) ''
1863 
1864  WRITE(fileunit,*) 'Activity:'
1865  WRITE(fileunit,*) active(1:nactive)
1866  WRITE(fileunit,*) ''
1867  end do
1868 
1869  if (allocated(times)) then
1870  deallocate(times)
1871  endif
1872 
1873  if (allocated(minimums)) then
1874  deallocate(minimums)
1875  endif
1876 
1877  if (allocated(maximums)) then
1878  deallocate(maximums)
1879  endif
1880 
1881  if (allocated(values)) then
1882  deallocate(values)
1883  endif
1884 
1885  if (allocated(active)) then
1886  deallocate(active)
1887  endif
1888 
1889  return
1890 END SUBROUTINE
1891 ! tdiReadScalar
1892 
1893 ! --------------------------------------------------------------------------
1894 ! FUNCTION TDI_READ_VECTOR
1895 ! PURPOSE Read a vector from an XMDF file and output information to
1896 ! to a text file
1897 ! NOTES
1898 ! --------------------------------------------------------------------------
1899 SUBROUTINE tdi_read_vector (a_xVectorId, FileUnit, error)
1900  INTEGER, INTENT(IN) :: a_xvectorid
1901  INTEGER, INTENT(IN) :: fileunit
1902  INTEGER, INTENT(OUT) :: error
1903  INTEGER ntimes, nvalues, nactive, ncomponents
1904  INTEGER itime, i
1905  LOGICAL*2 busereftime
1906  CHARACTER(LEN=100) timeunits
1907  REAL(DOUBLE), ALLOCATABLE :: times(:)
1908  REAL, ALLOCATABLE, DIMENSION (:, :) :: values
1909  REAL, ALLOCATABLE :: minimums(:), maximums(:)
1910  INTEGER, ALLOCATABLE :: active(:)
1911  REAL(DOUBLE) reftime
1912 
1913 ntimes = none
1914 nvalues = none
1915 nactive = none
1916 ncomponents = none
1917 
1918  ! read the time units
1919  call xf_get_dataset_time_units(a_xvectorid, timeunits, error)
1920  if (error < 0) return
1921 
1922  WRITE(fileunit,*) 'Time units: ', timeunits(1:len_trim(timeunits))
1923 
1924  ! see if we are using a reftime
1925  call xf_use_reftime(a_xvectorid, busereftime, error)
1926  if (error < 0) then
1927  return
1928  endif
1929  if (busereftime) then
1930  call xf_read_reftime(a_xvectorid, reftime, error)
1931  if (error < 0) then
1932  return
1933  endif
1934  !WRITE(FileUnit,*) 'Reftime: ', Reftime
1935  endif
1936 
1937  ! read in the number of values and number of active values
1938  call xf_get_dataset_numvals(a_xvectorid, nvalues, error)
1939  if (error .GE. 0) then
1940  call xf_get_dataset_numcomponents(a_xvectorid, ncomponents, error)
1941  if (error .GE. 0) then
1942  call xf_get_dataset_numactive(a_xvectorid, nactive, error)
1943  endif
1944  endif
1945  if (error .LT. 0) return
1946 
1947  if (nvalues <= 0) then
1948  WRITE(fileunit, *) 'No data to read in.'
1949  error = -1
1950  return
1951  endif
1952 
1953  ! read in the number of times
1954  call xf_get_dataset_num_times(a_xvectorid, ntimes, error)
1955  if (error < 0) then
1956  return
1957  endif
1958 
1959  ! Read in the individual time values
1960  allocate(times(ntimes))
1961 
1962  call xf_get_dataset_times(a_xvectorid, ntimes, times, error)
1963  if (error < 0) return
1964 
1965  ! Read in the minimum and maximum values
1966  allocate(minimums(ntimes))
1967  allocate(maximums(ntimes))
1968 
1969  call xf_get_dataset_mins(a_xvectorid, ntimes, minimums, error)
1970  if (error >= 0) then
1971  call xf_get_dataset_maxs(a_xvectorid, ntimes, maximums, error)
1972  endif
1973  if (error < 0) then
1974  deallocate(times)
1975  deallocate(minimums)
1976  deallocate(maximums)
1977  return
1978  endif
1979 
1980  allocate(values(ncomponents, nvalues))
1981  if (nactive .GT. 0) then
1982  allocate(active(nactive))
1983  endif
1984 
1985  WRITE(fileunit,*) 'Number Timesteps: ', ntimes
1986  WRITE(fileunit,*) 'Number Values: ', nvalues
1987  WRITE(fileunit,*) 'Number Components: ', ncomponents
1988  WRITE(fileunit,*) 'Number Active: ', nactive
1989 
1990  ! loop through the timesteps, read the values and active values and write
1991  ! them to the text file
1992  do itime = 1, ntimes
1993  call xf_read_vector_values_timestep(a_xvectorid, itime, nvalues, &
1994  ncomponents, values, error)
1995  if (error >= 0 .AND. nactive > 0) then
1996  call xf_read_activity_timestep(a_xvectorid, itime, nactive, active, error)
1997  endif
1998 
1999  ! Write the time, min, max, values and active values to the text output
2000  ! file.
2001  WRITE(fileunit,*) ''
2002  WRITE(fileunit,*) 'Timestep at ', times(itime)
2003  WRITE(fileunit,*) 'Min: ', minimums(itime)
2004  WRITE(fileunit,*) 'Max: ', maximums(itime)
2005 
2006  WRITE(fileunit,*) 'Values:'
2007  do i=1, nvalues
2008  WRITE(fileunit,*) values(1:ncomponents,i:i)
2009  enddo
2010  WRITE(fileunit,*) ''
2011 
2012  WRITE(fileunit,*) 'Activity:'
2013  WRITE(fileunit,*) active(1:nactive)
2014  WRITE(fileunit,*) ''
2015  WRITE(fileunit,*) ''
2016 
2017  end do
2018 
2019  if (allocated(times)) then
2020  deallocate(times)
2021  endif
2022 
2023  if (allocated(minimums)) then
2024  deallocate(minimums)
2025  endif
2026 
2027  if (allocated(maximums)) then
2028  deallocate(maximums)
2029  endif
2030 
2031  if (allocated(values)) then
2032  deallocate(values)
2033  endif
2034 
2035  if (allocated(active)) then
2036  deallocate(active)
2037  endif
2038 
2039  return
2040 END SUBROUTINE
2041 ! tdiReadVector
2042 
2043 END MODULE testdatasets

TestDatasets.f90 tests datasets

1 MODULE testgeompaths
2 
3 USE testdatasets
4 USE xmdf
5 USE errordefinitions
6 USE xmdfdefs
7 
8 CONTAINS
9 
10 SUBROUTINE tm_read_test_paths(Filename, OutFilename, error)
11  CHARACTER(LEN=*), INTENT(IN) :: filename, outfilename
12  INTEGER, INTENT(OUT) :: error
13  INTEGER outunit
14  INTEGER :: xfileid, xgroupid
15  INTEGER ngroups, nmaxpathlength, ndims, npaths, ntimes
16  INTEGER i, j, nstatus
17  CHARACTER,ALLOCATABLE :: paths(:)
18  REAL(DOUBLE), ALLOCATABLE :: times(:), locs(:)
19  REAL(DOUBLE) nullval
20  CHARACTER(LEN=BIG_STRING_SIZE) :: individualpath
21  INTEGER startloc
22  INTEGER, DIMENSION(2) :: pathindices
23  INTEGER itime, ipath
24 
25  ! open the XMDF file
26  CALL xf_open_file(filename, .false., xfileid, nstatus)
27  if (nstatus < 0) then
28  WRITE(*) 'Unable to open XMDF file TM_READ_TEST_PATHS.'
29  error = nstatus;
30  return
31  endif
32 
33  ! open the Output file
34  outunit = 79
35  OPEN(unit=outunit, file=outfilename, status='REPLACE', action='WRITE', &
36  iostat = error)
37  if (outunit == 0) then
38  call xf_close_file(xfileid, error)
39  error = -1
40  return
41  endif
42 
43  ! find the geomotric path groups
44  ! Get the number and paths of datasets in the file.
45  CALL xf_grp_pths_sz_for_geom_pths(xfileid, ngroups, &
46  nmaxpathlength, nstatus)
47  if (nstatus >= 0 .AND. ngroups > 0) then
48  allocate (paths(nmaxpathlength*ngroups))
49  CALL xf_grp_pths_for_geom_pths(xfileid, ngroups, &
50  nmaxpathlength, paths, nstatus)
51  endif
52  if (nstatus < 0) then
53  CALL xf_close_file(xfileid, nstatus)
54  error = -1
55  return
56  endif
57  ! Report the number and paths to individual geom paths groups in the file.
58  WRITE(outunit,*) 'Number of geometric paths in file: ', ngroups
59 
60  do i = 1, ngroups
61  startloc = (i-1)*nmaxpathlength + 1
62  individualpath = ''
63  do j = 1, nmaxpathlength - 1
64  individualpath(j:j) = paths(startloc+j-1)
65  enddo
66  WRITE(outunit,*) 'Reading particles in group: ', &
67  individualpath(1:len_trim(individualpath))
68 
69  CALL xf_open_group(xfileid, individualpath, xgroupid, nstatus)
70  if (nstatus >= 0) then
71  ! read the dimensionality of the paths
72  CALL xf_get_path_dimensionality(xgroupid, ndims, nstatus)
73  if (nstatus >= 0) then
74  WRITE(outunit,*) 'Group dimensionality:', ndims
75  ! read the number of paths
76  CALL xf_get_number_of_paths(xgroupid, npaths, nstatus)
77  if (nstatus >= 0) then
78  WRITE(outunit,*) 'Number of paths in group:', npaths
79  ! read the number of timesteps
80  CALL xf_get_number_of_times(xgroupid, ntimes, nstatus)
81  if (nstatus >= 0) then
82  WRITE(outunit,*) 'Number of timesteps in group:', ntimes
83  ! allocate times array
84  allocate(times(ntimes))
85  CALL xf_get_path_times_array(xgroupid, ntimes, times, nstatus)
86  if (nstatus >= 0) then
87  ! get space for the data location values
88  allocate(locs(npaths*ndims))
89  ! read null value
90  CALL xf_get_path_null_val(xgroupid, nullval, nstatus)
91  if (nstatus >= 0) then
92  WRITE(outunit,*) 'Null value:', nullval
93  do itime=1, ntimes
94  WRITE(outunit,*) 'Timestep: ', times(itime)
95  ! read the data for this timestep
96  CALL xf_read_path_locations_at_time(xgroupid, itime, &
97  1, npaths, locs, nstatus)
98  if (nstatus >= 0) then
99  ! write the data for this timestep
100  WRITE(outunit,*) ' X Y'
101  if (ndims == 3) then
102  WRITE(outunit,*) ' Z'
103  endif
104  WRITE(outunit,*) ''
105  do j=1, npaths
106  if (locs(j*ndims) == nullval) then
107  WRITE(outunit,*) 'Particle not active yet'
108  else
109  WRITE(outunit,*) locs((j-1)*ndims+1), ' ', locs((j-1)*ndims+2)
110  if (ndims == 3) then
111  WRITE(outunit,*) ' ', locs((j-1)*ndims+3)
112  endif
113  WRITE(outunit,*) ''
114  endif
115  enddo ! Loop through the paths
116  endif ! if timestep read
117  enddo ! Loop through timesteps
118  endif ! if null value read
119  if (allocated(locs)) deallocate (locs)
120  endif ! if we get get the times array
121  ! get space for the data location values - 1 particle
122  ! all times
123  allocate(locs(ntimes*ndims))
124  do ipath=1, npaths
125  ! read the values for this particle for all timesteps
126  CALL xf_read_path_locs_for_part(xgroupid, ipath, &
127  1, ntimes, locs, nstatus)
128  if (nstatus >= 0) then
129  ! write the data for this path
130  WRITE(outunit,*) 'Time X Y'
131  if (ndims == 3) then
132  WRITE(outunit,*) ' Z'
133  endif
134  WRITE(outunit, *) ''
135  do j=1, ntimes
136  if (locs((j-1)*ndims+1) .NE. nullval) then
137  WRITE(outunit,*) times(j), ' ', locs((j-1)*ndims + 1), ' ', &
138  locs((j-1)*ndims+2)
139  if (ndims == 3) then
140  WRITE(outunit,*) ' ', locs((j-1)*ndims+3)
141  endif
142  WRITE(outunit,*) ''
143  endif
144  enddo ! Loop through the times
145  endif ! if read path locations for particle
146  enddo ! Loop through the paths
147  if (allocated(locs)) deallocate(locs)
148  endif ! if allocation of locs suceeded
149  ! get space for the data location values - 2 particle
150  ! all times
151  allocate(locs(ntimes*ndims*2))
152  pathindices(1) = 1
153  pathindices(2) = npaths
154  CALL xf_read_path_locs_for_parts(xgroupid, 2, &
155  pathindices, 1, ntimes, locs, nstatus)
156  if (nstatus >= 0) then
157  ! write the data for these 2 paths
158  if (ndims == 3) then
159  WRITE(outunit,*) 'Timestep X1 Y1 Z1 Xn Yn Zn'
160  else
161  WRITE(outunit,*) 'Timestep X1 Y1 Xn Yn'
162  endif
163  do j=1, ntimes
164  if (ndims == 3) then
165  WRITE(outunit,*) times(j), ' ', locs((j-1)*2*ndims+1), ' ', &
166  locs((j-1)*2*ndims+2),' ', locs((j-1)*2*ndims+3),' ', &
167  locs((j-1)*2*ndims+4), locs((j-1)*2*ndims+5), ' ',locs((j-1)*2*ndims+6)
168  else
169  WRITE(outunit,*) times(j),' ', locs((j-1)*2*ndims + 1), ' ', &
170  locs((j-1)*2*ndims+2), ' ',locs((j-1)*2*ndims+3)
171  endif
172  enddo
173  else
174  WRITE(*,*) 'Error reading path locations for multiple particles'
175  endif
176  endif
177  if (allocated(locs)) deallocate(locs)
178  if (allocated(times)) deallocate(times)
179  CALL xf_close_group(xgroupid, nstatus)
180  if (nstatus < 0) then
181  WRITE(*)'Error reading geometric paths..'
182  endif
183  endif
184  endif
185  enddo ! loop through groups
186  ! free the paths
187  if (allocated(paths)) deallocate(paths)
188  ! close the files
189  close(outunit)
190  CALL xf_close_file(xfileid, nstatus)
191  error = nstatus;
192  return
193 
194 END SUBROUTINE
195 ! tmReadTestPaths
196 
197 SUBROUTINE tm_write_test_paths(Filename, Compression, error)
198  CHARACTER(LEN=*), INTENT(IN) :: filename
199  INTEGER, INTENT(IN) :: compression
200  INTEGER, INTENT(OUT) :: error
201  INTEGER npaths
202  REAL(DOUBLE), DIMENSION(6) :: plocs
203  REAL, DIMENSION(2) :: pspeeds
204  REAL(DOUBLE) nullval
205 
206  REAL nullvalsingleprec
207  INTEGER xfileid, xpathgroupid, xspeedid, xpropgroupid
208  INTEGER status
209  REAL(DOUBLE) timeval
210 
211  npaths = 0
212  nullval = -99999.9d0
213  nullvalsingleprec = -99999.9
214 
215  ! create the file
216  call xf_create_file(filename, .true., xfileid, status)
217  if (status < 0) then
218  error = -1
219  return
220  endif
221 
222  ! create the group to store the particle paths
223  call xf_create_geometric_path_group(xfileid, "particles", &
224  'abcdefglasdfjoaieur', compression, xpathgroupid, nullval, status)
225  if (status < 0) then
226  error = -1
227  return
228  endif
229 
230  ! create the data set to store the speed
231  call xf_create_scalar_dset_extndbl(xpathgroupid, 'Vmag', 'm/s', &
232  ts_seconds, nullvalsingleprec, compression, xspeedid, status)
233  if (status < 0) then
234  error = -1
235  return
236  endif
237 
238  call xf_create_property_group(xspeedid, xpropgroupid, status)
239  if (status < 0) then
240  call xf_close_group(xspeedid, status)
241  error = -1
242  return
243  endif
244  call xf_write_property_float(xpropgroupid, prop_null_value, 1, &
245  nullvalsingleprec, -1, status)
246  call xf_close_group(xpropgroupid, status)
247 
248  ! Setup the arrays for the path group at timestep 0
249  ! particle location at the first timestep
250  npaths = 1
251  plocs(1) = 1.0
252  plocs(2) = 2.0
253  plocs(3) = 3.0
254 
255  ! store the particles for the first timestep
256  timeval = 0.0
257  call xf_write_particle_timestep(xpathgroupid, 3, timeval, npaths, plocs, status)
258  if (status < 0) then
259  error = -1
260  return
261  endif
262  ! set up and store the speed at timestep 0
263  pspeeds(1) = 1.1
264  timeval = 0.0
265  call xf_write_scalar_timestep(xspeedid, timeval, 1, pspeeds, status)
266  if (status < 0) then
267  error = -1
268  return
269  endif
270 
271  ! Setup the arrays for the path group at timestep 1
272  ! particle location at the first timestep
273  plocs(1) = 4.0
274  plocs(2) = 5.0
275  plocs(3) = 6.0
276 
277  ! store the particles for the second timestep
278  timeval = 1.0
279  call xf_write_particle_timestep(xpathgroupid, 3, timeval, npaths, plocs, status)
280  if (status < 0) then
281  error = -1
282  return
283  endif
284 
285  ! set up and store the speed at timestep 2
286  pspeeds(1) = 1.2
287  call xf_write_scalar_timestep(xspeedid, timeval, 1, pspeeds, status)
288  if (status < 0) then
289  error = -1
290  return
291  endif
292 
293  ! Setup the arrays for the path group at timestep 3-add a particle
294  ! particle location at the first timestep
295  npaths = 2
296  plocs(1) = 7.0
297  plocs(2) = 8.0
298  plocs(3) = 9.0
299  plocs(4) = -1.0
300  plocs(5) = -2.0
301  plocs(6) = -3.0
302 
303  ! store the particles for the timestep 3
304  timeval = 2.0
305  call xf_write_particle_timestep(xpathgroupid, 3, timeval, npaths, plocs, status)
306  if (status < 0) then
307  error = -1
308  return
309  endif
310 
311  ! extend the data set for speed
312  call xf_extend_scalar_dataset(xspeedid, 2, status)
313  if (status < 0) then
314  error = -1
315  return
316  endif
317 
318  ! set up and store the speed at timestep 4
319  pspeeds(1) = 1.3
320  pspeeds(2) = 2.1
321  call xf_write_scalar_timestep(xspeedid, timeval, 2, pspeeds, status)
322  if (status < 0) then
323  error = -1
324  return
325  endif
326 
327  ! Setup the arrays for the path group at timestep 3-inactive particle(static)
328  ! particle location at the first timestep
329  plocs(1) = 7.0
330  plocs(2) = 8.0
331  plocs(3) = 9.0
332  plocs(4) = -4.0
333  plocs(5) = -5.0
334  plocs(6) = -6.0
335 
336  ! store the particles for timestep 4
337  timeval = 3.0
338  call xf_write_particle_timestep(xpathgroupid, 3, timeval, npaths, plocs, status)
339  if (status < 0) then
340  error = -1
341  return
342  endif
343 
344  ! set up and store the speed at timestep 4
345  pspeeds(1) = nullval
346  pspeeds(2) = 2.2
347  call xf_write_scalar_timestep(xspeedid, timeval, 2, pspeeds, status)
348  if (status < 0) then
349  error = -1
350  return
351  endif
352 
353  ! close the resources
354  call xf_close_group(xspeedid, status)
355  call xf_close_group(xpathgroupid, status)
356  call xf_close_file(xfileid, status)
357 
358  error = 1
359  return
360 END SUBROUTINE
361 
362 END MODULE testgeompaths

TestGeomPaths.f90 tests geometry paths

1 MODULE testgrid
2 
3 USE errordefinitions
4 USE xmdfdefs
5 USE xmdf
6 
7 CHARACTER(LEN=*),PARAMETER :: grid_cart2d_group_name = 'Grid Cart2D Group'
8 CHARACTER(LEN=*),PARAMETER :: grid_curv2d_group_name = 'Grid Curv2D Group'
9 CHARACTER(LEN=*),PARAMETER :: grid_cart3d_group_name = 'Grid Cart3D Group'
10 
11 CONTAINS
12 
13 !****************************
14 ! ---------------------------------------------------------------------------
15 ! FUNCTION TG_READ_GRID
16 ! PURPOSE Read a grid and write data to a text file
17 ! NOTES
18 ! ---------------------------------------------------------------------------
19 SUBROUTINE tg_read_grid(a_Id, a_Outfile, error)
20 INTEGER, INTENT(IN) :: a_id
21 INTEGER, INTENT(IN) :: a_outfile
22 INTEGER, INTENT(OUT) :: error
23 INTEGER ngridtype, nextrudetype, ndims, ncellsi, ncellsj
24 INTEGER ncellsk, nlayers, norientation, nvalsi, nvalsj
25 INTEGER nvalsk, nextrudevals, ncomporigin, nudir
26 CHARACTER(265) strgridtype, strextrudetype
27 LOGICAL*2 bdefined
28 REAL(DOUBLE) dorigin(3), dbearing, ddip, droll
29 REAL(DOUBLE), ALLOCATABLE :: dextrudevals(:), dcoordi(:), dcoordj(:)
30 REAL(DOUBLE), ALLOCATABLE :: dcoordk(:)
31 INTEGER i
32 
33 ngridtype = 0
34 nextrudetype = 0
35 ndims = 0
36 ncellsi = 0
37 ncellsj = 0
38 ncellsk = 0
39 nlayers = 0
40 norientation = 0
41 nvalsi = 0
42 nvalsj = 0
43 nvalsk = 0
44 nextrudevals = 0
45 ncomporigin = 1
46 nudir = 1
47 bdefined = .false.
48 dbearing = 0.0
49 ddip = 0.0
50 droll =0.0
51 i = 0
52 error = 1
53 
54  ! Grid type
55 call xf_get_grid_type(a_id, ngridtype, error)
56 if (error < 0) then
57  return
58 endif
59 SELECT CASE (ngridtype)
60  CASE (grid_type_cartesian)
61  strgridtype = 'Cartesian'
62  CASE (grid_type_curvilinear)
63  strgridtype = 'Curvilinear'
64  CASE (grid_type_cartesian_extruded)
65  strgridtype = 'Cartesian extruded'
66  CASE (grid_type_curvilinear_extruded)
67  strgridtype = 'Curvilinear extruded'
68  CASE default
69  WRITE(*,*) 'Invalid grid type'
70  error = -1
71 END SELECT
72 WRITE(a_outfile,*) 'The grid type is: ', strgridtype(1:len_trim(strgridtype))
73 
74  ! Number of dimensions
75 call xf_get_number_of_dimensions(a_id, ndims, error)
76 if (error .LT. 0) then
77  return
78 endif
79 if (ndims == 2) then
80  WRITE(a_outfile,*) 'The grid is two-dimensional'
81 elseif (ndims == 3) then
82  WRITE(a_outfile,*) 'The grid is three-dimensional'
83 else
84  WRITE(*,*) 'The grid dimensions are invalid'
85  error = -1
86 endif
87 
88  ! Extrusion type if applicable
89 if (ngridtype .EQ. grid_type_cartesian_extruded .OR. &
90  ngridtype .EQ. grid_type_curvilinear_extruded) then
91  call xf_get_extrusion_type(a_id, nextrudetype, error)
92  if (error < 0) then
93  return
94  endif
95  SELECT CASE (nextrudetype)
96  case (extrude_sigma)
97  strextrudetype = 'Sigma stretch'
98  case (extrude_cartesian)
99  strextrudetype = 'Cartesian'
100  case (extrude_curv_at_corners)
101  strextrudetype = 'Curvilinear at Corners'
102  case (extrude_curv_at_cells)
103  strextrudetype = 'Curvilinear at Cells'
104  END SELECT
105  WRITE(a_outfile,*) 'The grid is extruding using: ', &
106  strextrudetype(1:len_trim(strextrudetype))
107 endif
108 
109  ! Origin
110 call xf_origin_defined(a_id, bdefined, error)
111 if (error < 0) then
112  return
113 endif
114 if (bdefined) then
115  call xf_get_origin(a_id, dorigin(1), dorigin(2), dorigin(3), error)
116  if (error < 0) then
117  return
118  endif
119  WRITE(a_outfile,*) 'The grid origin is ', dorigin(1), ' ',&
120  dorigin(2), ' ', dorigin(3)
121 endif
122 
123  ! Orientation
124 call xf_get_orientation(a_id, norientation, error)
125 if (error < 0) then
126  return
127 endif
128 if (norientation == orientation_right_hand) then
129  WRITE(a_outfile,*) 'The grid has a right hand orientation'
130 elseif (norientation == orientation_left_hand) then
131  WRITE(a_outfile,*) 'The grid has a left hand orientation'
132 else
133  WRITE(*,*) 'Invalid grid orientation';
134  error = -1
135  return
136 endif
137 
138  ! Bearing
139 call xf_bearing_defined(a_id, bdefined, error)
140 if (error < 0) then
141  return
142 endif
143 if (bdefined) then
144  call xf_get_bearing(a_id, dbearing, error)
145  if (error < 0) then
146  return
147  endif
148  WRITE(a_outfile,*) 'The grid bearing is ', dbearing
149 endif
150 
151  ! Dip
152 call xf_dip_defined(a_id, bdefined, error)
153 if (error < 0) then
154  return
155 endif
156 if (bdefined) then
157  call xf_get_dip(a_id, ddip, error)
158  if (error < 0) then
159  return
160  endif
161  WRITE(a_outfile,*) 'The grid Dip is ', ddip
162 endif
163 
164 if (ndims == 3) then
165  ! Roll
166  call xf_roll_defined(a_id, bdefined, error)
167  if (error < 0) then
168  return
169  endif
170  if (bdefined) then
171  call xf_get_roll(a_id, droll, error)
172  if (error < 0) then
173  return
174  endif
175  WRITE(a_outfile,*) 'The grid Roll is ', droll
176  endif
177 endif
178 
179  ! Computational origin
180 call xf_computational_origin_defined(a_id, bdefined, error)
181 if (error < 0) then
182  return
183 endif
184 if (bdefined) then
185  call xf_get_computational_origin(a_id, ncomporigin, error)
186  if (error < 0) then
187  return
188  endif
189  WRITE(a_outfile,*) 'The grid Computational Origin is ', ncomporigin
190 else
191  WRITE(a_outfile,*) 'The grid Computational Origin is not defined';
192 endif
193 
194 
195  ! U Direction
196 call xf_get_u_direction_defined(a_id, bdefined, error)
197 if (error < 0) then
198  return
199 endif
200 if (bdefined) then
201  call xf_get_u_direction(a_id, nudir, error)
202  if (error < 0) then
203  return
204  endif
205  WRITE(a_outfile,*) 'The grid U Direction is ', nudir
206 else
207  WRITE(a_outfile,*) 'The grid U Direction is not defined'
208 endif
209 
210  ! number of cells in each direction
211 call xf_get_number_cells_in_i(a_id, ncellsi, error)
212 if (error >= 0) then
213  call xf_get_number_cells_in_j(a_id, ncellsj, error)
214  if ((error >= 0) .AND. (ndims == 3)) then
215  call xf_get_number_cells_in_k(a_id, ncellsk, error)
216  endif
217 endif
218 if (error < 0) then
219  return
220 endif
221 WRITE(a_outfile,*) 'Number of cells in I ', ncellsi
222 WRITE(a_outfile,*) 'Number of cells in J ', ncellsj
223 if (ndims == 3) then
224  WRITE(a_outfile,*) 'Number of cells in K ', ncellsk
225 endif
226 
227  ! Grid coordinates
228 if (ngridtype == grid_type_cartesian .OR. &
229  ngridtype == grid_type_cartesian_extruded) then
230  nvalsi = ncellsi
231  nvalsj = ncellsj
232  if (ndims == 3) then
233  nvalsk = ncellsk
234  endif
235 elseif (ngridtype == grid_type_curvilinear .OR. &
236  ngridtype == grid_type_curvilinear_extruded) then
237  if (ndims == 3) then
238  ! three dimensions
239  nvalsk = (ncellsi + 1) * (ncellsj + 1) * (ncellsk + 1)
240  nvalsj = nvalsk
241  nvalsi = nvalsj
242  else
243  ! two dimensions
244  nvalsj = (ncellsi + 1) * (ncellsj + 1)
245  nvalsi = nvalsj
246  endif
247 else
248  WRITE(*,*) 'Invalid grid type'
249  error = -1
250  return
251 endif
252 
253 ALLOCATE(dcoordi(nvalsi))
254 ALLOCATE(dcoordj(nvalsj))
255 if (ndims == 3) then
256  ALLOCATE(dcoordk(nvalsk))
257 endif
258 
259 call xf_get_grid_coords_i(a_id, nvalsi, dcoordi, error)
260 if (error >= 0) then
261  call xf_get_grid_coords_j(a_id, nvalsj, dcoordj, error)
262  if ((error >= 0) .AND. (ndims == 3)) then
263  call xf_get_grid_coords_k(a_id, nvalsk, dcoordk, error)
264  endif
265 endif
266 if (error < 0) then
267  WRITE(*,*) 'Error reading coordinates'
268  error = -1
269  return
270 endif
271 
272 WRITE(a_outfile,*) 'The Coordinates in direction I:'
273 do i = 1, nvalsi
274  if (mod(i,5) == 0) then
275  WRITE(a_outfile,*) ''
276  endif
277  WRITE(a_outfile,*) dcoordi(i)
278 enddo
279 WRITE(a_outfile,*) ''
280 
281 WRITE(a_outfile,*) 'The Coordinates in direction J:'
282 do i = 1, nvalsj
283  if (mod(i,5) == 0) then
284  WRITE(a_outfile,*) ''
285  endif
286  WRITE(a_outfile,*) dcoordj(i)
287 enddo
288 WRITE(a_outfile,*) ''
289 
290 if (ndims == 3) then
291  WRITE(a_outfile,*) 'The Coordinates in direction K:'
292  do i = 1, nvalsk
293  if (mod(i,5) == 0) then
294  WRITE(a_outfile,*) ''
295  endif
296  WRITE(a_outfile,*) dcoordk(i)
297  enddo
298 endif
299 WRITE(a_outfile,*) ''
300 
301 if (ALLOCATED(dcoordi)) DEALLOCATE(dcoordi)
302 if (ALLOCATED(dcoordj)) DEALLOCATE(dcoordj)
303 if (ALLOCATED(dcoordk)) DEALLOCATE(dcoordk)
304 
305 ! // Extrude data
306 if (ngridtype .EQ. grid_type_cartesian_extruded .OR. &
307  ngridtype .EQ. grid_type_curvilinear_extruded) then
308  call xf_get_extrude_num_layers(a_id, nlayers, error)
309  if (error < 0) then
310  return
311  endif
312 
313  SELECT CASE(nextrudetype)
314  case (extrude_sigma)
315  nextrudevals = nlayers
316  case (extrude_curv_at_corners)
317  nextrudevals = (ncellsi + 1) * (ncellsj + 1) * nlayers
318  case (extrude_curv_at_cells)
319  nextrudevals = ncellsi * ncellsj * nlayers
320  END SELECT
321 
322  ALLOCATE(dextrudevals(nextrudevals))
323 
324  call xf_get_extrude_values(a_id, nextrudevals, dextrudevals, error)
325  if (error < 0) then
326  return
327  endif
328 
329  WRITE(*,*) 'The extrude values are:'
330  do i = 1, nextrudevals
331  if (mod(i,5) == 0) then
332  WRITE(a_outfile,*) ''
333  endif
334  WRITE(a_outfile,*) dextrudevals(i)
335  enddo
336  if (ALLOCATED(dextrudevals)) DEALLOCATE(dextrudevals)
337 endif
338 
339 return
340 
341 END SUBROUTINE tg_read_grid
342 
343 !****************************
344 !----------------------------------------------------------------------------
345 ! FUNCTION TG_WRITE_TEST_GRID_CART_2D
346 ! PURPOSE Write a file that contains data for a 2D Cartesian Grid
347 ! NOTES A picture of the grid is in the file (TestGridCart2D.gif)
348 ! returns TRUE on success and FALSE on failure
349 !----------------------------------------------------------------------------
350 SUBROUTINE tg_write_test_grid_cart_2d(Filename, error)
351 CHARACTER(LEN=*), INTENT(IN) :: filename
352 INTEGER, INTENT(OUT) :: error
353 INTEGER ndimensions
354 INTEGER ncellsi, ncellsj
355 INTEGER ngridtype
356 INTEGER ncomporigin, nudir
357 REAL(DOUBLE) doriginx, doriginy, doriginz
358 INTEGER norientation
359 REAL(DOUBLE) dbearing
360 REAL(DOUBLE) planesi(5), planesj(5)
361 INTEGER i, j, ispczone
362 INTEGER xfileid, xgridid, xcoordid
363 INTEGER status
364 INTEGER tmpout1, tmpout2
365 
366  ndimensions = 2
367  ncellsi = 5
368  ncellsj = 5
369  ngridtype = grid_type_cartesian
370  ncomporigin = 4
371  nudir = -2
372  doriginx = 10.0
373  doriginy = 10.0
374  doriginz = 0.0
375  norientation = orientation_right_hand
376  dbearing = 45.0
377  xfileid = none
378  xgridid = none
379 
380  ! Fill in the grid plane data with a constant size of 30
381  do i = 1, ncellsi
382  planesi(i) = i*30.0
383  enddo
384  do j = 1, ncellsj
385  planesj(j) = j*30.0
386  enddo
387 
388  ! create the file
389  call xf_create_file(filename, .true., xfileid, status)
390  if (status < 0) then
391  error = -1
392  return
393  endif
394 
395  ! create the group to store the grid
396  call xf_create_group_for_grid(xfileid, grid_cart2d_group_name, xgridid, status)
397  if (status < 0) then
398  call xf_close_file(xfileid, error)
399  error = -1
400  return
401  endif
402 
403  ! Write the grid information to the file
404  call xf_set_grid_type(xgridid, ngridtype, tmpout1)
405  call xf_set_number_of_dimensions(xgridid, ndimensions, tmpout2)
406 
407  if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
408  call xf_close_group(xgridid, error)
409  call xf_close_file(xfileid, error)
410  error = -1
411  return
412  endif
413 
414  ! set origin and orientation
415  call xf_set_origin(xgridid, doriginx, doriginy, doriginz, tmpout1)
416  call xf_set_orientation(xgridid, norientation, tmpout2)
417 
418  if ((tmpout1 < 0) .OR. (tmpout2 .LT. 0)) then
419  call xf_close_group(xgridid, error)
420  call xf_close_file(xfileid, error)
421  error = -1
422  return
423  endif
424 
425  ! Set bearing
426  call xf_set_bearing(xgridid, dbearing, tmpout1)
427  if (tmpout1 < 0) then
428  call xf_close_group(xgridid, error)
429  call xf_close_file(xfileid, error)
430  error = -1
431  return
432  endif
433 
434  ! Set computational origin
435  call xf_set_computational_origin(xgridid, ncomporigin, tmpout1)
436  if (tmpout1 < 0) then
437  call xf_close_group(xgridid, error)
438  call xf_close_file(xfileid, error)
439  error = -1
440  return
441  endif
442 
443  ! Set u direction
444  call xf_set_u_direction(xgridid, nudir, tmpout1)
445  if (tmpout1 < 0) then
446  call xf_close_group(xgridid, error)
447  call xf_close_file(xfileid, error)
448  error = -1
449  return
450  endif
451 
452  ! Write the grid geometry to the file
453  ! Set the number of cells in each direction
454  call xf_set_number_cells_in_i(xgridid, ncellsi, tmpout1)
455  call xf_set_number_cells_in_j(xgridid, ncellsj, tmpout2)
456  if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
457  call xf_close_group(xgridid, error)
458  call xf_close_file(xfileid, error)
459  error = -1
460  return
461  endif
462 
463  ! Set the grid plane locations
464  call xf_set_grid_coords_i(xgridid, ncellsi, planesi, tmpout1)
465  call xf_set_grid_coords_j(xgridid, ncellsj, planesj, tmpout2)
466  if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
467  call xf_close_group(xgridid, error)
468  call xf_close_file(xfileid, error)
469  error = -1
470  return
471  endif
472 
473  ! Write Coordinate file - for GridCart2D, we will set the coordinate system
474  ! to be State Plane NAD27.
475  call xf_create_coordinate_group(xfileid, xcoordid, status)
476  if (status < 0) then
477  call xf_close_group(xgridid, error)
478  call xf_close_file(xfileid, error)
479  error = status
480  return
481  endif
482 
483  ispczone = 3601 ! Oregon North
484 
485  call xf_set_horiz_datum(xcoordid, horiz_datum_state_plane_nad27, error)
486  call xf_set_horiz_units(xcoordid, coord_units_us_feet, error)
487 
488  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
489  call xf_set_vert_units(xcoordid, coord_units_us_feet, error)
490 
491  ! write additional information
492  call xf_set_spc_zone(xcoordid, ispczone, error)
493 
494  call xf_close_group(xcoordid, error)
495  xcoordid = 0
496 
497  ! release memory
498 call xf_close_group(xgridid, error)
499 call xf_close_file(xfileid, error)
500 return
501 
502 END SUBROUTINE tg_write_test_grid_cart_2d
503 
504 !****************************
505 !------------------------------------------------------------------------------
506 ! FUNCTION tgWriteTestGridCurv2D
507 ! PURPOSE Write a file that contains data for a 2D Curvilinear Grid
508 ! NOTES A picture of the grid is TestGridCurv2D.gif
509 ! returns TRUE on success and FALSE on failure
510 !------------------------------------------------------------------------------
511 SUBROUTINE tg_write_test_grid_curv_2d(Filename, Compression, error)
512 CHARACTER(LEN=*), INTENT(IN) :: filename
513 INTEGER, INTENT(IN) :: compression
514 INTEGER, INTENT(OUT) :: error
515 INTEGER ndimensions
516 INTEGER ncomporigin, nudir
517 INTEGER ncellsi, ncellsj
518 INTEGER ncells
519 INTEGER ncorners
520 INTEGER ngridtype
521 REAL(DOUBLE) xvals(16), yvals(16) !There are 16 corners
522 INTEGER i
523 INTEGER xfileid, xgridid, xpropid, xdatasetsid, xscalarid
524 INTEGER xcoordid
525 REAL(DOUBLE) dnullvalue(1)
526 INTEGER norientation
527 REAL fdsetcellvals(6) !For cell-centered dataset
528 REAL fdsetcornervals(12) !For corner-centered dataset
529 REAL(DOUBLE) tempdouble, dcpplat, dcpplon
530 INTEGER*1 bdsetcellactive(6)
531 INTEGER*1 bdsetcorneractive(12)
532 INTEGER status
533 INTEGER tmpout1, tmpout2
534 
535  ndimensions = 2
536  ncomporigin = 1
537  nudir = 1;
538  ncellsi = 2
539  ncellsj = 3
540  ncells = ncellsi*ncellsj
541  ncorners = (ncellsi + 1)*(ncellsj + 1)
542  ngridtype = grid_type_curvilinear
543  xfileid = none
544  xgridid = none
545  xpropid = none
546  xdatasetsid = none
547  xscalarid = none
548  dnullvalue(1) = -999.0
549  norientation = orientation_right_hand
550 
551  ! There is no cell in the top right corner so we have a NullValue for
552  ! the top right corner
553 
554  ! xValues row by row
555  xvals(1) = 0.0
556  xvals(2) = 7.5
557  xvals(3) = 15.0
558  xvals(4) = 2.5
559  xvals(5) = 10.0
560  xvals(6) = 17.5
561  xvals(7) = 3.5
562  xvals(8) = 11.0
563  xvals(9) = 18.5
564  xvals(10) = 0.0
565  xvals(11) = 7.5
566  xvals(12) = dnullvalue(1)
567 
568  ! yValues row by row
569  yvals(1) = 0.0
570  yvals(2) = 0.0
571  yvals(3) = 0.0
572  yvals(4) = 10.0
573  yvals(5) = 10.0
574  yvals(6) = 10.0
575  yvals(7) = 20.0
576  yvals(8) = 20.0
577  yvals(9) = 20.0
578  yvals(10) = 30.0
579  yvals(11) = 30.0
580  yvals(12) = dnullvalue(1)
581 
582  ! cell centered velocity dataset values
583  fdsetcellvals(1) = 2.1
584  fdsetcellvals(2) = 2.0
585  fdsetcellvals(3) = 1.9
586  fdsetcellvals(4) = 2.3
587  fdsetcellvals(5) = 2.5
588  fdsetcellvals(6) = dnullvalue(1)
589 
590  ! all are active except the last value
591  do i = 1, ncells
592  bdsetcellactive(i) = 1
593  enddo
594  bdsetcellactive(ncells) = 0
595 
596  ! corner centered elevation dataset values
597  fdsetcornervals(1) = 1.0
598  fdsetcornervals(2) = 0.8
599  fdsetcornervals(3) = 1.2
600  fdsetcornervals(4) = 1.4
601  fdsetcornervals(5) = 1.8
602  fdsetcornervals(6) = 2.2
603  fdsetcornervals(7) = 1.8
604  fdsetcornervals(8) = 1.4
605  fdsetcornervals(9) = 2.0
606  fdsetcornervals(10) = 1.0
607  fdsetcornervals(11) = 1.8
608  fdsetcornervals(12) = 2.2
609 
610  ! all are active except the last value
611  do i = 1, ncorners
612  bdsetcorneractive(i) = 1
613  enddo
614  bdsetcorneractive(ncorners) = 0
615 
616  ! create the file
617  call xf_create_file(filename, .true., xfileid, status)
618  if (status < 0) then
619  error = -1
620  return
621  endif
622 
623  ! create the group to store the grid
624  call xf_create_group_for_grid(xfileid, grid_curv2d_group_name, xgridid, status)
625  if (status < 0) then
626  call xf_close_file(xfileid, error)
627  error = -1
628  return
629  endif
630 
631  ! Write the grid information to the file
632  call xf_set_grid_type(xgridid, ngridtype, tmpout1)
633  call xf_set_number_of_dimensions(xgridid, ndimensions, tmpout2)
634  if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
635  call xf_close_group(xgridid, error)
636  call xf_close_file(xfileid, error)
637  error = -1
638  return
639  endif
640 
641  ! set orientation
642  call xf_set_orientation(xgridid, norientation, tmpout1)
643  if (tmpout1 < 0 ) then
644  call xf_close_group(xgridid, error)
645  call xf_close_file(xfileid, error)
646  error = -1
647  return
648  endif
649 
650  ! Set computational origin
651  call xf_set_computational_origin(xgridid, ncomporigin, tmpout1)
652  if (tmpout1 < 0 ) then
653  call xf_close_group(xgridid, error)
654  call xf_close_file(xfileid, error)
655  error = -1
656  return
657  endif
658 
659  ! Set u direction
660  call xf_set_u_direction(xgridid, nudir, tmpout1)
661  if (tmpout1 < 0 ) then
662  call xf_close_group(xgridid, error)
663  call xf_close_file(xfileid, error)
664  error = -1
665  return
666  endif
667 
668  ! Write the grid geometry to the file
669  ! Set the number of cells in each direction
670  call xf_set_number_cells_in_i(xgridid, ncellsi, tmpout1)
671  call xf_set_number_cells_in_j(xgridid, ncellsj, tmpout2)
672  if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
673  call xf_close_group(xgridid, error)
674  call xf_close_file(xfileid, error)
675  error = -1
676  return
677  endif
678  ! Set a NullValue. This is used to identify locations in the grid that are
679  ! not being used. In our case no geometry is defined for the top right
680  ! corner.
681  call xf_create_grid_property_group(xgridid, xpropid, tmpout1)
682  if (xpropid < 0) then
683  call xf_close_group(xgridid, error)
684  call xf_close_file(xfileid, error)
685  error = -1
686  return
687  endif
688 
689  call xf_write_property_double(xpropid, prop_null_value, 1, dnullvalue, none, &
690  tmpout1)
691  if (tmpout1 < 0) then
692  call xf_close_group(xpropid, error)
693  call xf_close_group(xgridid, error)
694  call xf_close_file(xfileid, error)
695  error = -1
696  return
697  endif
698  call xf_close_group(xpropid, error)
699 
700  ! Set the grid plane locations
701  call xf_set_grid_coords_i(xgridid, ncorners, xvals, tmpout1)
702  call xf_set_grid_coords_j(xgridid, ncorners, yvals, tmpout2)
703  if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
704  call xf_close_group(xgridid, error)
705  call xf_close_file(xfileid, error)
706  error = -1
707  return
708  endif
709 
710  ! Create the datasets group
711  call xf_create_generic_group(xgridid, 'Datasets', xdatasetsid, tmpout1)
712  if (tmpout1 < 0) then
713  call xf_close_group(xgridid, error)
714  call xf_close_file(xfileid, error)
715  error = -1
716  return
717  endif
718 
719  ! Create the cell-centered dataset
720  call xf_create_scalar_dataset(xdatasetsid, 'Velocity Mag', 'ft/s', ts_minutes, &
721  compression, xscalarid, tmpout1)
722  if (tmpout1 < 0) then
723  call xf_close_group(xdatasetsid, error)
724  call xf_close_group(xgridid, error)
725  call xf_close_file(xfileid, error)
726  error = -1
727  return
728  endif
729 
730  ! specify that the dataset is cell-centered
731  call xf_scalar_data_location(xscalarid, grid_loc_center, tmpout1)
732  if (tmpout1 < 0) then
733  call xf_close_group(xscalarid, error)
734  call xf_close_group(xdatasetsid, error)
735  call xf_close_group(xgridid, error)
736  call xf_close_file(xfileid, error)
737  error = -1
738  return
739  endif
740 
741  ! Write the data
742  ! set a temporary variable to pass into the function
743  tempdouble = 0.0
744  call xf_write_scalar_timestep(xscalarid, tempdouble, ncells, fdsetcellvals, tmpout1)
745  call xf_write_activity_timestep(xscalarid, ncells, bdsetcellactive, tmpout2)
746  if ((tmpout1 < 0) .OR. (tmpout1 < 0)) then
747  call xf_close_group(xscalarid, error)
748  call xf_close_group(xdatasetsid, error)
749  call xf_close_group(xgridid, error)
750  call xf_close_file(xfileid, error)
751  error = -1
752  return
753  endif
754 
755  ! close the cell-centered dataset
756  call xf_close_group(xscalarid, error)
757 
758  ! Create the corner-centered dataset
759  call xf_create_scalar_dataset(xdatasetsid, 'elevation', 'ft', ts_minutes, &
760  compression, xscalarid, tmpout1)
761  if (tmpout1 < 0) then
762  call xf_close_group(xdatasetsid, error)
763  call xf_close_group(xgridid, error)
764  call xf_close_file(xfileid, error)
765  error = -1
766  return
767  endif
768 
769  ! specify that the dataset is corner-centered
770  call xf_scalar_data_location(xscalarid, grid_loc_corner, tmpout1)
771  if (tmpout1 < 0) then
772  call xf_close_group(xscalarid, error)
773  call xf_close_group(xdatasetsid, error)
774  call xf_close_group(xgridid, error)
775  call xf_close_file(xfileid, error)
776  error = -1
777  return
778  endif
779 
780  ! Write the data
781  ! set a temporary variable to pass into the function
782  tempdouble = 0.0
783  call xf_write_scalar_timestep(xscalarid, tempdouble, ncorners, fdsetcornervals, tmpout1)
784  call xf_write_activity_timestep(xscalarid, ncorners, bdsetcorneractive, tmpout2)
785  if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
786  call xf_close_group(xscalarid, error)
787  call xf_close_group(xdatasetsid, error)
788  call xf_close_group(xgridid, error)
789  call xf_close_file(xfileid, error)
790  error = -1
791  return
792  endif
793 
794  ! Write Coordinate file - for GridCurv2D, we will set the coordinate system
795  ! to be CPP, with CPP Latitude and CPP Longitude settings written
796  ! to the file.
797  call xf_create_coordinate_group(xfileid, xcoordid, status)
798  if (status < 0) then
799  call xf_close_group(xscalarid, error)
800  call xf_close_group(xdatasetsid, error)
801  call xf_close_group(xgridid, error)
802  call xf_close_file(xfileid, error)
803  error = status
804  return
805  endif
806 
807  dcpplat = 56.0; ! Made-up value
808  dcpplon = 23.0; ! Made-up value
809 
810  call xf_set_horiz_datum(xcoordid, horiz_datum_cpp, error)
811  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
812 
813  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
814  call xf_set_vert_units(xcoordid, coord_units_meters, error)
815 
816  ! write additional information
817  call xf_set_cpp_lat(xcoordid, dcpplat, error)
818  call xf_set_cpp_lon(xcoordid, dcpplon, error)
819 
820  call xf_close_group(xcoordid, error)
821  xcoordid = 0
822 
823  ! release memory
824 call xf_close_group(xscalarid, error)
825 call xf_close_group(xdatasetsid, error)
826 call xf_close_group(xgridid, error)
827 call xf_close_file(xfileid, error)
828 return
829 
830 END SUBROUTINE tg_write_test_grid_curv_2d
831 
832 !****************************
833 !------------------------------------------------------------------------------
834 ! FUNCTION tgWriteTestGridCart3D
835 ! PURPOSE Write a file that contains data for a 2D Cartesian Grid
836 ! NOTES A picture of the grid is in the file (TestGridCart2D.gif)
837 ! returns TRUE on success and FALSE on failure
838 !------------------------------------------------------------------------------
839 SUBROUTINE tg_write_test_grid_cart_3d(Filename, Compression, error)
840 CHARACTER(LEN=*), INTENT(IN) :: filename
841 INTEGER, INTENT(IN) :: compression
842 INTEGER, INTENT(OUT) :: error
843 INTEGER ndimensions
844 INTEGER ncomporigin, nudir
845 INTEGER ncellsi, ncellsj, ncellsk
846 INTEGER ngridtype
847 REAL(DOUBLE) doriginx, doriginy, doriginz
848 INTEGER norientation
849 REAL(DOUBLE) dbearing, ddip, droll
850 REAL(DOUBLE) planesi(5), planesj(5), planesk(3)
851 INTEGER i, j, status, ispczone
852 INTEGER xfileid, xgridid, xpropid, xcoordid
853 INTEGER ncells
854 INTEGER active(75)
855 INTEGER tmpout1, tmpout2, tmpout3
856 
857 ndimensions = 3
858 ncomporigin = 8
859 nudir = -2
860 ncellsi = 5
861 ncellsj = 5
862 ncellsk = 3
863 ngridtype = grid_type_cartesian
864 doriginx = 10.0
865 doriginy = 10.0
866 doriginz = 0.0
867 norientation = orientation_right_hand
868 dbearing = 45.0
869 ddip = 0.0
870 droll = 0.0
871 xfileid = none
872 xgridid = none
873 xpropid = none
874 ncells = ncellsi * ncellsj * ncellsk
875 
876  ! Fill in the grid plane data with a constant size of 30
877 do i = 1, ncellsi
878  planesi(i) = i*30.0
879 enddo
880 do j = 1, ncellsj
881  planesj(j) = j*30.0
882 enddo
883 do j = 1, ncellsk
884  planesk(j) = j*30.0
885 enddo
886 
887  ! fill in the activity array
888  ! default array to active
889 do i = 1, ncells
890  active(i) = 1
891 enddo
892 
893  ! two cells are inactive (identified by array index)
894  ! i = 0, j = 0, k = 0 and i = 4, j = 4, k = 0
895 active(1) = 0
896 active(4*ncellsj*ncellsk+4*ncellsk+1) = 0
897 
898  ! create the file
899 call xf_create_file(filename, .true., xfileid, tmpout1)
900 if (tmpout1 < 0) then
901  error = -1
902  return
903 endif
904 
905  ! create the group to store the grid
906 call xf_create_group_for_grid(xfileid, grid_cart3d_group_name, xgridid, tmpout1)
907 if (tmpout1 < 0) then
908  call xf_close_file(xfileid, error)
909  error = -1
910  return
911 endif
912 
913  ! Write the grid information to the file
914 call xf_set_grid_type(xgridid, ngridtype, tmpout1)
915 call xf_set_number_of_dimensions(xgridid, ndimensions, tmpout2)
916 if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
917  call xf_close_group(xgridid, error)
918  call xf_close_file(xfileid, error)
919  error = -1
920  return
921 endif
922 
923  ! set origin and orientation
924 call xf_set_origin(xgridid, doriginx, doriginy, doriginz, tmpout1)
925 call xf_set_orientation(xgridid, norientation, tmpout2)
926 if ((tmpout1 < 0) .OR. (tmpout2 < 0)) then
927  call xf_close_group(xgridid, error)
928  call xf_close_file(xfileid, error)
929  error = -1
930  return
931 endif
932 
933  ! Set bearing, dip and roll
934 call xf_set_bearing(xgridid, dbearing, tmpout1)
935 call xf_set_dip(xgridid, ddip, tmpout2)
936 call xf_set_roll(xgridid, droll, tmpout3)
937 if ((tmpout1 < 0) .OR. (tmpout2 < 0) .OR. (tmpout3 < 0)) then
938  call xf_close_group(xgridid, error)
939  call xf_close_file(xfileid, error)
940  error = -1
941  return
942 endif
943 
944  ! Set computational origin
945 call xf_set_computational_origin(xgridid, ncomporigin, tmpout1)
946 if (tmpout1 < 0) then
947  call xf_close_group(xgridid, error)
948  call xf_close_file(xfileid, error)
949  error = -1
950  return
951 endif
952 
953  ! Set u direction
954 call xf_set_u_direction(xgridid, nudir, tmpout1)
955 if (tmpout1 < 0) then
956  call xf_close_group(xgridid, error)
957  call xf_close_file(xfileid, error)
958  error = -1
959  return
960 endif
961 
962  ! Write the grid geometry to the file
963  ! Set the number of cells in each direction
964 call xf_set_number_cells_in_i(xgridid, ncellsi, tmpout1)
965 call xf_set_number_cells_in_j(xgridid, ncellsj, tmpout2)
966 call xf_set_number_cells_in_k(xgridid, ncellsk, tmpout3)
967 if ((tmpout1 < 0) .OR. (tmpout2 < 0) .OR. (tmpout3 < 0)) then
968  call xf_close_group(xgridid, error)
969  call xf_close_file(xfileid, error)
970  error = -1
971  return
972 endif
973 
974  ! Set the grid plane locations
975 call xf_set_grid_coords_i(xgridid, ncellsi, planesi, tmpout1)
976 call xf_set_grid_coords_j(xgridid, ncellsj, planesj, tmpout2)
977 call xf_set_grid_coords_k(xgridid, ncellsk, planesk, tmpout3)
978 
979 if ((tmpout1 < 0) .OR. (tmpout2 < 0) .OR. (tmpout3 < 0)) then
980  call xf_close_group(xgridid, error)
981  call xf_close_file(xfileid, error)
982  error = -1
983  return
984 endif
985 
986  ! Write the activity array
987 call xf_create_grid_cell_prop_grp(xgridid, xpropid, tmpout1)
988 if (xpropid < 0) then
989  call xf_close_group(xgridid, error)
990  call xf_close_file(xfileid, error)
991  error = -1
992  return
993 endif
994 
995 call xf_write_property_int(xpropid, prop_activity, ncells, active, &
996  compression, tmpout1)
997 
998 if (tmpout1 < 0) then
999  call xf_close_group(xpropid, error)
1000  call xf_close_group(xgridid, error)
1001  call xf_close_file(xfileid, error)
1002  error = -1
1003  return
1004 endif
1005 
1006 call xf_close_group(xpropid, error)
1007 
1008  ! Write Coordinate file - for GridCart3D, we will set the coordinate system
1009  ! to be State Plane NAD27.
1010  call xf_create_coordinate_group(xfileid, xcoordid, status)
1011  if (status < 0) then
1012  call xf_close_group(xgridid, error)
1013  call xf_close_file(xfileid, error)
1014  error = status
1015  endif
1016 
1017  ispczone = 3601; ! Oregon North
1018 
1019  call xf_set_horiz_datum(xcoordid, horiz_datum_state_plane_nad27, error)
1020  call xf_set_horiz_units(xcoordid, coord_units_us_feet, error)
1021 
1022  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
1023  call xf_set_vert_units(xcoordid, coord_units_us_feet, error)
1024 
1025  ! write additional information
1026  call xf_set_spc_zone(xcoordid, ispczone, error)
1027 
1028  call xf_close_group(xcoordid, error)
1029  xcoordid = 0
1030 
1031  ! release memory
1032 call xf_close_group(xgridid, error)
1033 call xf_close_file(xfileid, error)
1034 return
1035 
1036 END SUBROUTINE tg_write_test_grid_cart_3d
1037 
1038 END MODULE TestGrid

TestGrid.f90 tests grids

1 MODULE testmesh
2 
3 USE testdatasets
4 USE xmdf
5 USE errordefinitions
6 USE xmdfdefs
7 
8 CHARACTER(LEN=*), PARAMETER :: mesh_a_group_name = 'MeshA Group'
9 CHARACTER(LEN=*), PARAMETER :: mesh_b_group_name = 'MeshB Group'
10 
11 
12 CONTAINS
13 
14 !****************************
15 ! ---------------------------------------------------------------------------
16 ! FUNCTION TM_READ_MESH
17 ! PURPOSE Read a mesh file and write a text output file
18 ! NOTES
19 ! ---------------------------------------------------------------------------
20 SUBROUTINE tm_read_mesh (xGroupId, a_OutFile, error)
21 INTEGER, INTENT(IN) :: xgroupid
22 INTEGER, INTENT(IN) :: a_outfile
23 INTEGER, INTENT(OUT) :: error
24 INTEGER nelems, nnodes, nnodesperelem, nelemtype, nnodeid
25 LOGICAL*2 elementsonetype
26 INTEGER status, i, j
27 INTEGER strtype, uinttype, inttype, dbltype, floattype
28 INTEGER xpropgroupid
29 INTEGER, ALLOCATABLE :: elemtypes(:), nodesinelem(:)
30 REAL(DOUBLE), ALLOCATABLE :: xnodelocs(:), ynodelocs(:), znodelocs(:)
31 
32  ! Get the number of elements, nodes, and Maximum number of nodes per element
33  call xf_get_number_of_elements(xgroupid, nelems, status)
34  if (status >= 0) then
35  call xf_get_number_of_nodes(xgroupid, nnodes, status)
36  if (status >= 0) then
37  call xf_get_max_nodes_in_elem(xgroupid, nnodesperelem, status)
38  endif
39  endif
40 
41  if (status < 0) then
42  error = -1
43  return
44  endif
45 
46  ! Do Element information first
47  WRITE(a_outfile,*) 'Number of Elements: ', nelems
48 
49  ! Element types
50  call xf_are_all_elems_same_type(xgroupid, elementsonetype, status)
51  if (status < 0) then
52  error = -1
53  return
54  endif
55 
56  if (elementsonetype) then
57  call xf_read_elem_types_single_value(xgroupid, nelemtype, status)
58  WRITE(a_outfile,*) 'All elements are type ', nelemtype
59  else
60  allocate (elemtypes(nelems))
61  call xf_read_elem_types(xgroupid, nelems, elemtypes, status)
62  if (status < 0) then
63  error = -1
64  return
65  endif
66  WRITE(a_outfile,*) 'Element Types:'
67  do i=1, nelems
68  WRITE(a_outfile,*) 'Elem ', i, ', ', 'Type ', elemtypes(i)
69  enddo
70  deallocate (elemtypes)
71  endif
72 
73  ! Nodes in each element
74  allocate (nodesinelem(nelems*nnodesperelem))
75  call xf_read_elem_node_ids(xgroupid, nelems, nnodesperelem, nodesinelem, error)
76  if (error < 0) then
77  WRITE (*,*) 'Error reading mesh'
78  return
79  endif
80 
81  do i=1, nelems
82  WRITE(a_outfile,*) 'Elem: ', i, ' - '
83  do j=1, nnodesperelem
84  nnodeid = nodesinelem((i-1)*nnodesperelem + j)
85  if (nnodeid > 0) then ! -1 is for unused array locations
86  WRITE(a_outfile,*) nnodeid, ' '
87  endif
88  enddo
89  WRITE (a_outfile,*) ''
90  enddo
91 
92  if (allocated(nodesinelem)) deallocate (nodesinelem)
93 
94  ! NodeLocations
95  allocate (xnodelocs(nnodes))
96  allocate (ynodelocs(nnodes))
97  allocate (znodelocs(nnodes))
98 
99  call xf_read_x_node_locations(xgroupid, nnodes, xnodelocs, status)
100  if (status >= 0) then
101  call xf_read_y_node_locations(xgroupid, nnodes, ynodelocs, status)
102  if (status >= 0) then
103  call xf_read_z_node_locations(xgroupid, nnodes, znodelocs, status)
104  endif
105  endif
106 
107  WRITE(a_outfile,*) 'Node Locations:'
108  do i=1, nnodes
109  WRITE(a_outfile,*) 'Node: ', i, ' Location: ', xnodelocs(i), ' ', &
110  ynodelocs(i), ' ', znodelocs(i)
111  enddo
112 
113  deallocate (xnodelocs)
114  deallocate (ynodelocs)
115  deallocate (znodelocs)
116 
117  ! Open Properties Group
118  call xf_open_group(xgroupid, 'PROPERTIES', xpropgroupid, status)
119  if (status < 0) then
120  WRITE(a_outfile,*) ''
121  WRITE(a_outfile,*) 'Properties Group not found'
122  WRITE(a_outfile,*) ''
123  error = -1
124  return
125  endif
126 
127  call xf_get_property_type(xpropgroupid, 'String', strtype, status)
128  call xf_get_property_type(xpropgroupid, 'UInt', uinttype, status)
129  call xf_get_property_type(xpropgroupid, 'Int', inttype, status)
130  call xf_get_property_type(xpropgroupid, 'Double', dbltype, status)
131  call xf_get_property_type(xpropgroupid, 'Float', floattype, status)
132 
133  ! Property Types:
134  WRITE(a_outfile,*) ''
135  if (strtype == xf_type_string) then
136  WRITE(a_outfile,*) 'String Property Type Read Correctly'
137  else
138  WRITE(a_outfile,*) 'Error in Getting String Property Type'
139  endif
140  if (uinttype == xf_type_uint) then
141  WRITE(a_outfile,*) 'Unsigned Integer Property Type Read Correctly'
142  else
143  WRITE(a_outfile,*) 'Error in Getting Unsigned Integer Property Type'
144  endif
145  if (inttype == xf_type_int) then
146  WRITE(a_outfile,*) 'Integer Property Type Read Correctly'
147  else
148  WRITE(a_outfile,*) 'Error in Getting Integer Property Type'
149  endif
150  if (dbltype == xf_type_double) then
151  WRITE(a_outfile,*) 'Double Property Type Read Correctly'
152  else
153  WRITE(a_outfile,*) 'Error in Getting Double Property Type'
154  endif
155  if (floattype == xf_type_float) then
156  WRITE(a_outfile,*) 'Float Property Type Read Correctly'
157  else
158  WRITE(a_outfile,*) 'Error in Getting Float Property Type'
159  endif
160  WRITE(a_outfile,*) ''
161 
162  error = 0
163  return
164 
165 END SUBROUTINE
166 
167 !****************************
168 !------------------------------------------------------------------------------
169 ! FUNCTION TM_WRITE_TEST_MESH_A
170 ! PURPOSE Write a file that contains data for an all tri mesh
171 ! NOTES A picture of the mesh is in the file (TestMeshA.gif)
172 ! error equals TRUE on success and FALSE on failure
173 !------------------------------------------------------------------------------
174 SUBROUTINE tm_write_test_mesh_a (Filename, Compression, error)
175 CHARACTER(LEN=*), INTENT(IN) :: filename
176 INTEGER, INTENT(IN) :: compression
177 INTEGER, INTENT(OUT) :: error
178 INTEGER nelements, nnodes
179 INTEGER xfileid, xmeshid, xpropgroupid, xcoordid
180 REAL(DOUBLE), DIMENSION(5) :: dnodelocsx
181 REAL(DOUBLE), DIMENSION(5) :: dnodelocsy
182 REAL(DOUBLE), DIMENSION(5) :: dnodelocsz
183 INTEGER, DIMENSION(3,3) :: ielementnodes
184 INTEGER status, elemtype, propint(1), iellipse
185 CHARACTER(LEN=BIG_STRING_SIZE) propstring
186 INTEGER propuint(1)
187 REAL(DOUBLE) propdouble(1), dmajorr, dminorr
188 REAL propfloat
189 
190  ! set Element type to EL_TYPE_TRI_LINEAR
191  elemtype = el_type_tri_linear
192 
193  nelements = 3
194  nnodes = 5
195  xfileid = none
196  xmeshid = none
197 
198  ! Setup the arrays for the mesh data
199  ! nodes
200  dnodelocsx(1) = 0.0
201  dnodelocsx(2) = 5.0
202  dnodelocsx(3) = 0.0
203  dnodelocsx(4) = 5.0
204  dnodelocsx(5) = 7.5
205 
206  dnodelocsy(1) = 5.0
207  dnodelocsy(2) = 5.0
208  dnodelocsy(3) = 0.0
209  dnodelocsy(4) = 0.0
210  dnodelocsy(5) = 2.5
211 
212  dnodelocsz(1) = 0.0
213  dnodelocsz(2) = 0.0
214  dnodelocsz(3) = 0.0
215  dnodelocsz(4) = 0.0
216  dnodelocsz(5) = 0.0
217 
218  ! nodes for each element must be counter-clockwize
219  ielementnodes(1,1) = 1
220  ielementnodes(2,1) = 3
221  ielementnodes(3,1) = 2
222 
223  ielementnodes(1,2) = 2
224  ielementnodes(2,2) = 3
225  ielementnodes(3,2) = 4
226 
227  ielementnodes(1,3) = 5
228  ielementnodes(2,3) = 2
229  ielementnodes(3,3) = 4
230 
231  ! create the file
232  call xf_create_file(filename, .true., xfileid, status)
233  if (status < 0) then
234  ! close the resources
235  call xf_close_file(xfileid, error)
236  error = -1
237  return
238  endif
239 
240  ! create the group to store the mesh
241  call xf_create_group_for_mesh(xfileid, mesh_a_group_name, xmeshid, status)
242  if (status < 0) then
243  ! close the resources
244  call xf_close_group(xmeshid, error)
245  call xf_close_file(xfileid, error)
246  error = -1
247  return
248  endif
249 
250  ! Element types - all are linear triangles
251  call xf_set_all_elems_same_type(xmeshid, elemtype, status)
252  if (status < 0) then
253  ! close the resources
254  call xf_close_group(xmeshid, error)
255  call xf_close_file(xfileid, error)
256  error = -1
257  return
258  endif
259 
260  ! node information
261  call xf_set_number_of_nodes(xmeshid, nnodes, status)
262  if (status < 0) then
263  ! close the resources
264  call xf_close_group(xmeshid, error)
265  call xf_close_file(xfileid, error)
266  error = -1
267  return
268  endif
269 
270  call xf_write_x_node_locations(xmeshid, nnodes, dnodelocsx, compression, status)
271  if (status < 0) then
272  ! close the resources
273  call xf_close_group(xmeshid, error)
274  call xf_close_file(xfileid, error)
275  error = -1
276  return
277  endif
278 
279  call xf_write_y_node_locations(xmeshid, nnodes, dnodelocsy, status)
280  if (status < 0) then
281  ! close the resources
282  call xf_close_group(xmeshid, error)
283  call xf_close_file(xfileid, error)
284  error = -1
285  return
286  endif
287 
288  call xf_write_z_node_locations(xmeshid, nnodes, dnodelocsz, status)
289  if (status < 0) then
290  ! close the resources
291  call xf_close_group(xmeshid, error)
292  call xf_close_file(xfileid, error)
293  error = -1
294  return
295  endif
296 
297  ! element information
298  call xf_set_number_of_elements(xmeshid, nelements, status)
299  if (status < 0) then
300  ! close the resources
301  call xf_close_group(xmeshid, error)
302  call xf_close_file(xfileid, error)
303  error = -1
304  return
305  endif
306 
307  ! Write the node array ids for the nodes in each element
308  call xf_write_elem_node_ids(xmeshid, nelements, 3, ielementnodes, &
309  compression, status)
310  if (status < 0) then
311  ! close the resources
312  call xf_close_group(xmeshid, error)
313  call xf_close_file(xfileid, error)
314  error = -1
315  return
316  endif
317 
318  ! Write the Property File
319  call xf_create_mesh_property_group(xmeshid, xpropgroupid, status)
320  if (status < 0) then
321  call xf_close_group(xmeshid, error)
322  call xf_close_file(xfileid, error)
323  error = -1
324  return
325  endif
326 
327  propstring = 'Property String'
328  propuint = 5
329  propint = -5
330  propdouble = 5.6789012345d0
331  propfloat = 5.6789
332 
333  call xf_write_property_string(xpropgroupid, 'String', propstring, status)
334  call xf_write_property_uint(xpropgroupid, 'UInt', 1, propuint, none, status)
335  call xf_write_property_int(xpropgroupid, 'Int', 1, propint, none, status)
336  call xf_write_property_double(xpropgroupid, 'Double', 1, &
337  propdouble, none, status)
338  call xf_write_property_float(xpropgroupid, 'Float', 1, &
339  propfloat, none, status)
340 
341  ! Write Coordinate file - for MeshA, we will set the coordinate system to be
342  ! Geogrpahic, with Latitude, Longitude, and user-defined ellipsoid settings
343  ! written to the file.
344  call xf_create_coordinate_group(xfileid, xcoordid, status)
345  if (status < 0) then
346  call xf_close_group(xpropgroupid, error)
347  call xf_close_group(xmeshid, error)
348  call xf_close_file(xfileid, error)
349  error = status
350  return
351  endif
352 
353  ! set coordinate values
354  iellipse = 32 ! User defined
355  dmajorr = 45.0 ! Made up
356  dminorr = 32.0 ! Made up
357 
358  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic, error)
359  call xf_set_horiz_units(xcoordid, coord_units_us_feet, error)
360 
361  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
362  call xf_set_vert_units(xcoordid, coord_units_us_feet, error)
363 
364  ! write additional information
365  call xf_set_ellipse(xcoordid, iellipse, error)
366  call xf_set_lat(xcoordid, latitude_north, error)
367  call xf_set_lon(xcoordid, longitude_east, error)
368  call xf_set_major_r(xcoordid, dmajorr, error)
369  call xf_set_minor_r(xcoordid, dminorr, error)
370 
371  call xf_close_group(xcoordid, error)
372  xcoordid = 0
373 
374  ! close the resources
375  call xf_close_group(xpropgroupid, error)
376  call xf_close_group(xmeshid, error)
377  call xf_close_file(xfileid, error)
378 
379  return
380 
381 END SUBROUTINE
382 
383 !****************************
384 !------------------------------------------------------------------------------
385 ! FUNCTION TM_WRITE_TEST_MESH_B
386 ! PURPOSE Write a file that contains data for an mixed quad/tri linear mesh
387 ! NOTES A picture of the mesh is in the file (TestMeshB.gif)
388 ! error equals TRUE on success and FALSE on failure
389 !------------------------------------------------------------------------------
390 SUBROUTINE tm_write_test_mesh_b (Filename, Compression, error)
391 CHARACTER(LEN=*), INTENT(IN) :: filename
392 INTEGER, INTENT(IN) :: compression
393 INTEGER, INTENT(OUT) :: error
394 INTEGER nelements, nnodes, nmaxnodeperelem
395 INTEGER xfileid, xmeshid, xpropgroupid, xcoordid
396 REAL(DOUBLE), DIMENSION(5) :: dnodelocsx
397 REAL(DOUBLE), DIMENSION(5) :: dnodelocsy
398 REAL(DOUBLE), DIMENSION(5) :: dnodelocsz
399 INTEGER, DIMENSION(4,2) :: ielementnodes
400 INTEGER, DIMENSION(2) :: ielementtypes
401 INTEGER status, propint(1), iellipse
402 CHARACTER(LEN=BIG_STRING_SIZE) propstring
403 INTEGER propuint(1)
404 REAL(DOUBLE) propdouble(1)
405 REAL propfloat
406 
407  nelements = 2
408  nnodes = 5
409  nmaxnodeperelem = 4
410  xfileid = none
411  xmeshid = none
412 
413  ! Setup the arrays for the mesh data
414  ! nodes
415  dnodelocsx(1) = 0.0
416  dnodelocsx(2) = 5.0
417  dnodelocsx(3) = 0.0
418  dnodelocsx(4) = 5.0
419  dnodelocsx(5) = 7.5
420 
421  dnodelocsy(1) = 5.0
422  dnodelocsy(2) = 5.0
423  dnodelocsy(3) = 0.0
424  dnodelocsy(4) = 0.0
425  dnodelocsy(5) = 2.5
426 
427  dnodelocsz(1) = 0.0
428  dnodelocsz(2) = 0.0
429  dnodelocsz(3) = 0.0
430  dnodelocsz(4) = 0.0
431  dnodelocsz(5) = 0.0
432 
433  ! nodes for each element must be counter-clockwize
434  ielementnodes(1,1) = 1
435  ielementnodes(2,1) = 3
436  ielementnodes(3,1) = 4
437  ielementnodes(4,1) = 2
438 
439  ielementnodes(1,2) = 2
440  ielementnodes(2,2) = 4
441  ielementnodes(3,2) = 5
442  ielementnodes(4,2) = none;
443 
444  ielementtypes(1) = el_type_quadrilateral_linear;
445  ielementtypes(2) = el_type_tri_linear;
446 
447  ! create the file
448  call xf_create_file(filename, .true., xfileid, status)
449  if (status < 0) then
450  ! close the resources
451  call xf_close_file(xfileid, error)
452  error = -1
453  return
454  endif
455 
456  ! create the group to store the mesh
457  call xf_create_group_for_mesh(xfileid, mesh_b_group_name, xmeshid, status)
458  if (status < 0) then
459  ! close the resources
460  call xf_close_group(xmeshid, error)
461  call xf_close_file(xfileid, error)
462  error = -1
463  return
464  endif
465 
466  ! node information
467  call xf_set_number_of_nodes(xmeshid, nnodes, status)
468  if (status < 0) then
469  ! close the resources
470  call xf_close_group(xmeshid, error)
471  call xf_close_file(xfileid, error)
472  error = -1
473  return
474  endif
475 
476  call xf_write_x_node_locations(xmeshid, nnodes, dnodelocsx, &
477  compression, status)
478  if (status < 0) then
479  ! close the resources
480  call xf_close_group(xmeshid, error)
481  call xf_close_file(xfileid, error)
482  error = -1
483  return
484  endif
485 
486  call xf_write_y_node_locations(xmeshid, nnodes, dnodelocsy, status)
487  if (status < 0) then
488  ! close the resources
489  call xf_close_group(xmeshid, error)
490  call xf_close_file(xfileid, error)
491  error = -1
492  return
493  endif
494 
495  call xf_write_z_node_locations(xmeshid, nnodes, dnodelocsz, status)
496  if (status < 0) then
497  ! close the resources
498  call xf_close_group(xmeshid, error)
499  call xf_close_file(xfileid, error)
500  error = -1
501  return
502  endif
503 
504  ! element information
505  call xf_set_number_of_elements(xmeshid, nelements, status)
506  if (status < 0) then
507  ! close the resources
508  call xf_close_group(xmeshid, error)
509  call xf_close_file(xfileid, error)
510  error = -1
511  return
512  endif
513 
514  ! Element types
515  call xf_write_elem_types(xmeshid, nelements, ielementtypes, &
516  compression, status)
517  if (status < 0) then
518  ! close the resources
519  call xf_close_group(xmeshid, error)
520  call xf_close_file(xfileid, error)
521  error = -1
522  return
523  endif
524 
525  ! Write the node array ids for the nodes in each element
526  call xf_write_elem_node_ids(xmeshid, nelements, nmaxnodeperelem, &
527  ielementnodes, compression, status)
528  if (status < 0) then
529  ! close the resources
530  call xf_close_group(xmeshid, error)
531  call xf_close_file(xfileid, error)
532  error = -1
533  return
534  endif
535 
536  ! Write the Property File
537  call xf_create_mesh_property_group(xmeshid, xpropgroupid, status)
538  if (status < 0) then
539  call xf_close_group(xmeshid, error)
540  call xf_close_group(xfileid, error)
541  error = -1
542  return
543  endif
544 
545  propstring = 'String Property'
546  propuint = 2
547  propint = -2
548  propdouble = 2.3456789012d0
549  propfloat = 2.3456
550 
551  call xf_write_property_string(xpropgroupid, 'String', propstring, status)
552  call xf_write_property_uint(xpropgroupid, 'UInt', 1, propuint, none, status)
553  call xf_write_property_int(xpropgroupid, 'Int', 1, propint, none, status)
554  call xf_write_property_double(xpropgroupid, 'Double', 1, &
555  propdouble, none, status)
556  call xf_write_property_float(xpropgroupid, 'Float', 1, &
557  propfloat, none, status)
558 
559  ! Write Coordinate file - for MeshB, we will set the coordinate system to be
560  ! Geogrpahic, with Latitude, Longitude, and standard ellipsoid settings
561  ! written to the file.
562  call xf_create_coordinate_group(xfileid, xcoordid, status)
563  if (status < 0) then
564  call xf_close_group(xpropgroupid, error)
565  call xf_close_group(xmeshid, error)
566  call xf_close_file(xfileid, error)
567  error = status
568  return
569  endif
570 
571  ! set coordinate values
572  iellipse = 21 ! International 1924
573 
574  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic, error)
575  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
576 
577  call xf_set_vert_datum(xcoordid, vert_datum_ngvd_88, error)
578  call xf_set_vert_units(xcoordid, coord_units_meters, error)
579 
580  ! write additional information
581  call xf_set_ellipse(xcoordid, iellipse, error)
582  call xf_set_lat(xcoordid, latitude_south, error)
583  call xf_set_lon(xcoordid, longitude_west, error)
584 
585  call xf_close_group(xcoordid, error)
586  xcoordid = 0
587 
588  ! close the resources
589  call xf_close_group(xpropgroupid, error)
590  call xf_close_group(xmeshid, error)
591  call xf_close_file(xfileid, error)
592 
593  return
594 
595 END SUBROUTINE
596 
597 END MODULE testmesh

TestMesh.f90 tests meshes

1 MODULE testtimestep
2 
3 USE xmdf
4 
5 CHARACTER(LEN=*), PARAMETER :: datasets_location = 'Datasets'
6 CHARACTER(LEN=*), PARAMETER :: scalar_a_location = 'Scalars/ScalarA'
7 CHARACTER(LEN=*), PARAMETER :: scalar_b_location = 'Scalars/ScalarB'
8 CHARACTER(LEN=*), PARAMETER :: vector2d_a_location = 'Vectors/Vector2D_A'
9 CHARACTER(LEN=*), PARAMETER :: vector2d_b_location = 'Vectors/Vector2D_B'
10 
11 CONTAINS
12 
13 ! --------------------------------------------------------------------------
14 ! FUNCTION ttiTestNumTimes
15 ! PURPOSE Change the NumTimes to truncate timesteps
16 ! NOTES
17 ! --------------------------------------------------------------------------
18 RECURSIVE SUBROUTINE tti_test_num_times (a_DatasetId, a_Itimestep, error)
19 INTEGER, INTENT(IN) :: a_datasetid
20 INTEGER, INTENT(IN) :: a_itimestep
21 INTEGER, INTENT(OUT) :: error
22 INTEGER numtimes
23 INTEGER itimestep
24 
25 error = 1
26 
27 ! Fortran loops are 1 to 3 but C is 0 to 2, etc
28 itimestep = a_itimestep - 1
29 
30  ! truncate just written timestep and test error conditions
31 if (1 == itimestep .OR. 3 == itimestep .OR. 5 == itimestep) then
32  ! Test setting NumTimes after end of dataset
33  call xf_set_dataset_num_times( a_datasetid, itimestep + 2, error )
34  if (error >= 0) then
35  WRITE(*,*) 'ERROR1: XF_SET_DATASET_NUM_TIMES must return ERROR.'
36  endif
37 
38  if (1 == itimestep) then
39  itimestep = 1;
40  endif
41  if (3 == itimestep) then
42  itimestep = 2;
43  endif
44  if (5 == itimestep) then
45  itimestep = 3;
46  endif
47 
48  ! Write actual NumTimes
49  call xf_set_dataset_num_times( a_datasetid, itimestep, error )
50  if (error < 0) then
51  WRITE(*,*) 'ERROR2: xfSetDatasetNumTimes must NOT return error.'
52  endif
53 
54  ! Test setting NumTimes after end step.
55  call xf_set_dataset_num_times( a_datasetid, itimestep + 1, error )
56  if (error >= 0) then
57  WRITE(*,*) 'ERROR3: xfSetDatasetNumTimes must return ERROR.'
58  endif
59 
60  ! Test reading NumTimes
61  call xf_get_dataset_num_times( a_datasetid, numtimes, error )
62  if (error < 0) then
63  WRITE(*,*) 'ERROR4: xfSetDatasetNumTimes must NOT return error.'
64  endif
65  if (numtimes .NE. itimestep) then
66  WRITE(*,*) 'ERROR5: xfGetDatasetNumTimes must return CORRECT NumTimes.'
67  endif
68 endif
69 
70 return
71 
72 END SUBROUTINE
73 !ttiTestNumTimes
74 ! --------------------------------------------------------------------------
75 ! FUNCTION ttReadDatasets
76 ! PURPOSE Read a dataset group from an XMDF file and output information to
77 ! to a text file
78 ! NOTES
79 ! --------------------------------------------------------------------------
80 RECURSIVE SUBROUTINE tt_read_datasets (a_xGroupId, a_FileUnit, error)
81 INTEGER, INTENT(IN) :: a_xgroupid
82 INTEGER, INTENT(IN) :: a_fileunit
83 INTEGER, INTENT(OUT) :: error
84 INTEGER npaths, nmaxpathlength, j
85 CHARACTER, ALLOCATABLE, DIMENSION(:) :: paths
86 CHARACTER(LEN=500) individualpath
87 INTEGER nstatus, i
88 INTEGER xscalarid, xvectorid, xmultiid
89 INTEGER nmultidatasets
90 
91 xscalarid = none
92 xvectorid = none
93 
94 nmultidatasets = 0
95 npaths = 0
96 nmaxpathlength = 0
97 
98  ! Look for scalar datasets
99 call xf_get_scalar_datasets_info(a_xgroupid, npaths, nmaxpathlength, nstatus)
100 if (nstatus >= 0 .AND. npaths > 0) then
101  allocate(paths(npaths*nmaxpathlength))
102  call xf_get_scalar_dataset_paths(a_xgroupid, npaths, nmaxpathlength, paths, &
103  error)
104 endif
105 if (nstatus < 0) then
106  error = -1
107  return
108 endif
109 
110  ! Output number and paths to scalar datasets
111 WRITE(a_fileunit,*) 'Number of Scalars ', npaths
112 do i=2, npaths
113  individualpath = ''
114  do j=1, nmaxpathlength-1
115  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
116  enddo
117  WRITE(a_fileunit,*) 'Reading scalar: ', individualpath(1:nmaxpathlength-1)
118  call xf_open_group(a_xgroupid, individualpath(1:nmaxpathlength-1), &
119  xscalarid, nstatus)
120  if (nstatus < 0) then
121  error = -1
122  return
123  endif
124 
125  call tti_read_scalar(xscalarid, a_fileunit, nstatus)
126  call xf_close_group(xscalarid, error)
127  if (nstatus < 0) then
128  WRITE(*,*) 'Error reading scalar dataset.'
129  error = -1
130  return
131  endif
132 enddo
133 
134 if (allocated(paths)) deallocate(paths)
135  ! Look for vector datasets
136 call xf_get_vector_datasets_info(a_xgroupid, npaths, nmaxpathlength, nstatus)
137 if (nstatus >= 0 .AND. npaths > 0) then
138  allocate(paths(npaths*nmaxpathlength))
139  call xf_get_vector_dataset_paths(a_xgroupid, npaths, nmaxpathlength, paths, error)
140 endif
141 if (nstatus < 0) then
142  error = -1
143  return
144 endif
145 
146  ! Output number and paths to scalar datasets
147 WRITE(a_fileunit,*) 'Number of Vectors ', npaths
148 do i=2, npaths
149  do j=1, nmaxpathlength-1
150  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
151  enddo
152  WRITE(a_fileunit,*) 'Reading Vector: ', &
153  individualpath(1:nmaxpathlength-1)
154  call xf_open_group(a_xgroupid, individualpath(1:nmaxpathlength-1), &
155  xvectorid, nstatus)
156  if (nstatus < 0) then
157  error = -1
158  return
159  endif
160  call tti_read_vector(xvectorid, a_fileunit, nstatus)
161  call xf_close_group(xvectorid, error)
162  if (nstatus < 0) then
163  WRITE(*,*) 'Error reading vector dataset.'
164  error = -1
165  return
166  endif
167 enddo
168 
169 if (allocated(paths)) deallocate(paths)
170 
171 ! find multidataset folders
172 call xf_get_grp_pths_sz_mlt_dsets(a_xgroupid, nmultidatasets, &
173  nmaxpathlength, nstatus)
174 if (nstatus >= 0 .AND. nmultidatasets > 0) then
175  allocate(paths(nmultidatasets*nmaxpathlength))
176  call xf_get_all_grp_paths_mlt_dsets(a_xgroupid, nmultidatasets, &
177  nmaxpathlength, paths, error)
178  if (nstatus < 0) then
179  error = -1
180  return
181  endif
182 
183  ! Output number and paths to multidatasets
184  WRITE(a_fileunit,*) 'Number of Multidatasets ', nmultidatasets
185  do i=2, nmultidatasets
186  individualpath = ''
187  do j=1, nmaxpathlength-1
188  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
189  enddo
190  WRITE(a_fileunit,*) 'Reading multidataset: ', &
191  individualpath(1:nmaxpathlength-1)
192  call xf_open_group(a_xgroupid, individualpath(1:nmaxpathlength-1), &
193  xmultiid, nstatus)
194  if (nstatus < 0) then
195  error = -1
196  return
197  endif
198 
199  call tt_read_datasets(xmultiid, a_fileunit, nstatus)
200  call xf_close_group(xmultiid, error)
201  if (nstatus < 0) then
202  WRITE(*,*) 'Error reading multidatasets.'
203  error = -1
204  return
205  endif
206  enddo
207 endif
208 if (allocated(paths)) deallocate(paths)
209 
210 error = 1
211 return
212 
213 END SUBROUTINE
214 !ttReadDatasets
215 ! --------------------------------------------------------------------------
216 ! FUNCTION ttReadActivityScalarAIndex
217 ! PURPOSE Read all timestep values for a particular index
218 ! NOTES
219 ! --------------------------------------------------------------------------
220 SUBROUTINE tt_read_activity_scalar_a_index(a_Filename, a_Index, error)
221 CHARACTER(LEN=*), INTENT(IN) :: a_filename
222 INTEGER, INTENT(IN) :: a_index
223 INTEGER, INTENT(OUT) :: error
224 INTEGER status
225 INTEGER xfileid, xdsetsid, xscalaraid
226 INTEGER ntimesteps, i
227 INTEGER, ALLOCATABLE :: bactive(:)
228 
229 xfileid = none
230 xdsetsid = none
231 xscalaraid = none
232 
233  ! open the file
234 call xf_open_file(a_filename, .true., xfileid, status)
235 if (status < 0) then
236  error = -1
237  return
238 endif
239 
240  ! open the dataset group
241 call xf_open_group(xfileid, datasets_location, xdsetsid, status)
242 if (status >= 0) then
243  call xf_open_group(xdsetsid, scalar_a_location, xscalaraid, status)
244 endif
245 if (status < 0) then
246  error = status
247  return
248 endif
249 
250  ! Find out the number of timesteps in the file
251 CALL xf_get_dataset_num_times(xscalaraid, ntimesteps, status)
252 if (status < 0) then
253  error = status
254  return
255 endif
256 
257 if (ntimesteps < 1) then
258  error = -1
259  return
260 endif
261 
262  ! Read the values for the index
263 allocate(bactive(ntimesteps))
264 call xf_read_active_vals_at_index(xscalaraid, a_index, 1, ntimesteps, &
265  bactive, status)
266  ! output the data
267 WRITE(*,*) ''
268 WRITE(*,*) 'Reading activity for scalar A slice at index: ', a_index
269 do i=1, ntimesteps
270  WRITE(*,*) bactive(i), ' '
271 enddo
272 
273 deallocate(bactive)
274 
275 error = status
276 return
277 
278 END SUBROUTINE
279 ! ttReadActivityScalarAIndex
280 
281 ! --------------------------------------------------------------------------
282 ! FUNCTION ttReadScalarAIndex
283 ! PURPOSE Read all timestep values for a particular index
284 ! NOTES
285 ! --------------------------------------------------------------------------
286 SUBROUTINE tt_read_scalar_a_index (a_Filename, a_Index, error)
287 CHARACTER(LEN=*), INTENT(IN) :: a_filename
288 INTEGER, INTENT(IN) :: a_index
289 INTEGER, INTENT(OUT) :: error
290 INTEGER status
291 INTEGER xfileid, xdsetsid, xscalaraid
292 INTEGER ntimesteps, i
293 REAL, ALLOCATABLE :: fvalues(:)
294 
295 xfileid = none
296 xdsetsid = none
297 xscalaraid = none
298 
299  ! open the file
300 call xf_open_file(a_filename, .true., xfileid, status)
301 if (status < 0) then
302  error = -1
303  return
304 endif
305 
306  ! open the dataset group
307 call xf_open_group(xfileid, datasets_location, xdsetsid, status)
308 if (status >= 0) then
309  call xf_open_group(xdsetsid, scalar_a_location, xscalaraid, status)
310 endif
311 if (status < 0) then
312  error = status
313  return
314 endif
315 
316  ! Find out the number of timesteps in the file
317 call xf_get_dataset_num_times(xscalaraid, ntimesteps, status)
318 if (status < 0) then
319  error = status
320  return
321 endif
322 
323 if (ntimesteps < 1) then
324  error = -1
325  return
326 endif
327 
328  ! Read the values for the index
329 allocate (fvalues(ntimesteps))
330 call xf_read_scalar_values_at_index(xscalaraid, a_index, 1, ntimesteps, &
331  fvalues, status)
332 
333  ! output the data
334 WRITE(*,*) ''
335 WRITE(*,*) 'Reading scalar A slice at index: ', a_index
336 do i=1, ntimesteps
337  WRITE(*,*) fvalues(i), ' '
338 enddo
339 
340 deallocate(fvalues)
341 
342 error = status
343 return
344 
345 END SUBROUTINE
346 ! ttReadScalarAtIndex
347 
348 ! --------------------------------------------------------------------------
349 ! FUNCTION ttWriteScalarA
350 ! PURPOSE Write scalar Dataset to an HDF5 File
351 ! NOTES This tests dynamic data sets, and activity
352 ! This dataset is dynamic concentrations (mg/L) with output times
353 ! in minutes.
354 ! Dataset is for a mesh and so nActive is the number of elements
355 ! which is not the same as the nValues which would be number of nodes
356 ! reads/writes a reference time in julian days
357 ! --------------------------------------------------------------------------
358 SUBROUTINE tt_write_scalar_a (a_Filename, a_Compression, error)
359  CHARACTER(LEN=*), INTENT(IN) :: a_filename
360  INTEGER, INTENT(IN) :: a_compression
361  INTEGER, INTENT(OUT) :: error
362  INTEGER xfileid, xdsetsid, xscalaraid, xcoordid
363  INTEGER nvalues, ntimes, nactive
364  REAL(DOUBLE) dtime, djulianreftime
365  INTEGER itimestep, iactive, ihpgnzone
366  REAL fvalues(10) ! nValues
367  INTEGER*1 bactivity(10) ! activity
368  INTEGER i, status
369 
370  ! initialize the data
371  nvalues = 10
372  ntimes = 3
373  nactive = 8
374  dtime = 0.0
375 
376  ! 5th item in data set is always inactive, others active
377  do iactive = 1, nactive
378  bactivity(iactive) = 1
379  enddo
380  bactivity(6) = 0
381 
382 
383  ! create the file
384  call xf_create_file(a_filename, .true., xfileid, error)
385  if (error .LT. 0) then
386  ! close the file
387  call xf_close_file(xfileid, error)
388  return
389  endif
390 
391  ! create the group where we will put all the datasets
392  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
393  if (status < 0) then
394  call xf_close_file(xfileid, error)
395  error = -1
396  return
397  endif
398 
399  ! Create the scalar A dataset group
400  call xf_create_scalar_dataset(xdsetsid, scalar_a_location, 'mg/L', &
401  ts_hours, a_compression, xscalaraid, status)
402  if (status .LT. 0) then
403  ! close the dataset
404  call xf_close_group(xscalaraid, error)
405  call xf_close_group(xdsetsid, error)
406  call xf_close_file(xfileid, error)
407  error = status
408  return
409  endif
410 
411  ! Add in a reftime. This is a julian day for:
412  ! noon July 1, 2003
413  djulianreftime = 2452822.0;
414  call xf_write_reftime(xscalaraid, djulianreftime, status)
415  if (status < 0) then
416  call xf_close_group(xscalaraid, error)
417  call xf_close_group(xdsetsid, error)
418  call xf_close_file(xfileid, error)
419  endif
420 
421  ! Loop through timesteps adding them to the file
422  do itimestep = 1, ntimes
423  ! We will have an 0.5 hour timestep
424  dtime = itimestep * 0.5
425 
426  fvalues(1) = dtime
427  do i = 2, nvalues
428  fvalues(i) = fvalues(i-1)*2.5
429  end do
430 
431  ! write the dataset array values
432  call xf_write_scalar_timestep(xscalaraid, dtime, nvalues, fvalues, error)
433  if (error .GE. 0) then
434  ! write activity array
435  call xf_write_activity_timestep(xscalaraid, nactive, bactivity, error)
436  end if
437 
438  call tti_test_num_times(xscalaraid, itimestep, error)
439  enddo
440 
441  ! Write Coordinate file - for ScalarA, we will set the coordinate system
442  ! to be Geographic HPGN, with HPGN settings written to the file.
443  call xf_create_coordinate_group(xfileid, xcoordid, status)
444  if (status < 0) then
445  call xf_close_group(xscalaraid, error)
446  call xf_close_group(xdsetsid, error)
447  call xf_close_file(xfileid, error)
448  error = -1
449  return
450  endif
451 
452  ! set HPGN Zone for test
453  ihpgnzone = 29 ! Utah
454  ! Write Coordinate Information to file
455  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic_hpgn, error)
456  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
457  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
458  call xf_set_vert_units(xcoordid, coord_units_meters, error)
459 
460  ! write additional information
461  call xf_set_hpgn_area(xcoordid, ihpgnzone, error)
462 
463  call xf_close_group(xcoordid, error)
464  xcoordid = 0;
465 
466  ! close the dataset
467  call xf_close_group(xscalaraid, error)
468  call xf_close_group(xdsetsid, error)
469  call xf_close_file(xfileid, error)
470 
471  return
472 END SUBROUTINE
473 ! ttWriteScalarA
474 
475 ! --------------------------------------------------------------------------
476 ! FUNCTION TT_WRITE_SCALAR_B
477 ! PURPOSE Write scalar Dataset to an HDF5 File
478 ! NOTES This tests dynamic data sets, and activity
479 ! This dataset is dynamic concentrations (mg/L) with output times
480 ! in minutes.
481 ! Dataset is for a mesh and so nActive is the number of elements
482 ! which is not the same as the nValues which would be number of nodes
483 ! reads/writes a reference time in julian days
484 ! --------------------------------------------------------------------------
485 SUBROUTINE tt_write_scalar_b (a_Filename, a_Compression, a_Overwrite, error)
486  CHARACTER(LEN=*), INTENT(IN) :: a_filename
487  INTEGER, INTENT(IN) :: a_compression
488  LOGICAL, INTENT(IN) :: a_overwrite
489  INTEGER, INTENT(OUT) :: error
490  INTEGER xfileid, xdsetsid, xscalarbid, xcoordid
491  INTEGER nvalues, ntimes, nactive
492  REAL(DOUBLE) dtime, djulianreftime
493  INTEGER itimestep, iactive
494  REAL fvalues(10) ! nValues
495  INTEGER*1 bactivity(10) ! activity
496  INTEGER i, status
497 
498  ! initialize the data
499  nvalues = 10
500  ntimes = 3
501  nactive = 8
502  dtime = 0.0
503  i = 0
504 
505  ! 5th item in data set is always inactive, others active
506  do iactive = 1, nactive
507  bactivity(iactive) = 1
508  enddo
509  bactivity(6) = 0
510 
511  if (a_overwrite) then
512  ! open the already-existing file
513  call xf_open_file(a_filename, .false., xfileid, status)
514  if (status < 0) then
515  error = -1
516  return
517  endif
518  ! open the group where we have all the datasets
519  call xf_open_group(xfileid, datasets_location, xdsetsid, status)
520  if (status < 0) then
521  call xf_close_file(xfileid, error)
522  error = -1
523  return
524  endif
525  else
526  ! create the file
527  call xf_create_file(a_filename, .true., xfileid, error)
528  if (error .LT. 0) then
529  ! close the file
530  call xf_close_file(xfileid, error)
531  return
532  endif
533 
534  ! create the group where we will put all the datasets
535  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
536  if (status < 0) then
537  call xf_close_file(xfileid, error)
538  error = -1
539  return
540  endif
541  endif
542 
543  ! Create/Overwrite the scalar B dataset group
544  call xf_create_scalar_dataset(xdsetsid, scalar_b_location, 'mg/L', &
545  ts_hours, a_compression, xscalarbid, status)
546  if (status < 0) then
547  ! close the dataset
548  call xf_close_group(xscalarbid, error)
549  call xf_close_group(xdsetsid, error)
550  call xf_close_file(xfileid, error)
551  error = status
552  return
553  endif
554 
555  ! Add in a reftime. This is a julian day for:
556  ! noon July 1, 2003
557  djulianreftime = 2452822.0;
558  call xf_write_reftime(xscalarbid, djulianreftime, status)
559  if (status < 0) then
560  call xf_close_group(xscalarbid, error)
561  call xf_close_group(xdsetsid, error)
562  call xf_close_file(xfileid, error)
563  endif
564 
565  if (.NOT. a_overwrite) then
566  ! Loop through timesteps adding them to the file
567  do itimestep = 1, ntimes
568  ! We will have an 0.5 hour timestep
569  dtime = itimestep * 0.5
570 
571  fvalues(1) = dtime
572  do i = 2, nvalues
573  fvalues(i) = fvalues(i-1)*2.5
574  end do
575 
576  ! write the dataset array values
577  call xf_write_scalar_timestep(xscalarbid, dtime, nvalues, fvalues, error)
578  if (error .GE. 0) then
579  ! write activity array
580  call xf_write_activity_timestep(xscalarbid, nactive, bactivity, error)
581  end if
582  if (error < 0) then
583  call xf_close_group(xscalarbid, error)
584  call xf_close_group(xdsetsid, error)
585  call xf_close_file(xfileid, error)
586  endif
587  enddo
588  else
589  ! Loop through timesteps adding them to the file
590  do itimestep = 1, ntimes
591  ! We will have an 1.5 hour timestep
592  dtime = itimestep * 1.5
593 
594  fvalues(1) = dtime
595  do i = 2, nvalues
596  fvalues(i) = fvalues(i-1)*1.5
597  end do
598 
599  ! write the dataset array values
600  call xf_write_scalar_timestep(xscalarbid, dtime, nvalues, fvalues, error)
601  if (error .GE. 0) then
602  ! write activity array
603  call xf_write_activity_timestep(xscalarbid, nactive, bactivity, error)
604  end if
605  if (error < 0) then
606  call xf_close_group(xscalarbid, error)
607  call xf_close_group(xdsetsid, error)
608  call xf_close_file(xfileid, error)
609  endif
610  enddo
611  endif
612 
613  if (.NOT. a_overwrite) then
614  ! Write Coordinate file - for ScalarB, we will set the coordinate system
615  ! to be UTM, with UTM Zone settings written to the file.
616  call xf_create_coordinate_group(xfileid, xcoordid, status)
617  if (status < 0) then
618  call xf_close_group(xscalarbid, error)
619  call xf_close_group(xdsetsid, error)
620  call xf_close_file(xfileid, error)
621  error = -1
622  return
623  endif
624 
625  ! Write Coord Info to file
626  call xf_set_horiz_datum(xcoordid, horiz_datum_utm, error)
627  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
628 
629  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
630  call xf_set_vert_units(xcoordid, coord_units_meters, error)
631 
632  ! write additional information - we'll use the max value for this test
633  call xf_set_utm_zone(xcoordid, utm_zone_max, error)
634 
635  call xf_close_group(xcoordid, error)
636  xcoordid = 0
637  endif
638 
639  ! close the dataset
640  call xf_close_group(xscalarbid, error)
641  call xf_close_group(xdsetsid, error)
642  call xf_close_file(xfileid, error)
643 
644  error = 1
645  return
646 END SUBROUTINE
647 ! ttWriteScalarB
648 !------------------------------------------------------------------------------
649 ! FUNCTION TT_WRITE_COORDS_TO_MULTI
650 ! PURPOSE Write coordinate system to a multidataset file
651 ! NOTES
652 !------------------------------------------------------------------------------
653 SUBROUTINE tt_write_coords_to_multi (a_xFileId, error)
654 INTEGER, INTENT(IN) :: a_xfileid
655 INTEGER, INTENT(OUT) :: error
656 INTEGER xcoordid
657 INTEGER status
658 
659  ! Write Coordinate file - for Multidatasets, we will set the coordinate system
660  ! to be UTM, with UTM Zone settings written to the file.
661  call xf_create_coordinate_group(a_xfileid, xcoordid, status)
662  if (status < 0) then
663  error = status
664  return
665  endif
666 
667  ! Write Coord Info to file
668  call xf_set_horiz_datum(xcoordid, horiz_datum_utm, error)
669  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
670 
671  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
672  call xf_set_vert_units(xcoordid, coord_units_meters, error)
673 
674  ! write additional information - we'll use the max value for this test
675  call xf_set_utm_zone(xcoordid, utm_zone_max, error)
676 
677  call xf_close_group(xcoordid, error)
678  xcoordid = 0
679 
680  return
681 END SUBROUTINE
682 
683 ! --------------------------------------------------------------------------
684 ! FUNCTION ttWriteScalarAToMulti
685 ! PURPOSE Write scalar Dataset to an HDF5 File
686 ! NOTES This tests dynamic data sets, and activity
687 ! This dataset is dynamic concentrations (mg/L) with output times
688 ! in minutes.
689 ! Dataset is for a mesh and so nActive is the number of elements
690 ! which is not the same as the nValues which would be number of nodes
691 ! reads/writes a reference time in julian days
692 ! --------------------------------------------------------------------------
693 SUBROUTINE tt_write_scalar_a_to_multi (a_GroupID, status)
694  ! CHARACTER(LEN=*), INTENT(IN) :: a_Filename
695  ! INTEGER, INTENT(IN) :: a_Compression
696  ! INTEGER, INTENT(OUT) :: error
697  INTEGER xfileid, xdsetsid, xscalaraid
698  INTEGER a_groupid
699  INTEGER nvalues, ntimes, nactive
700  REAL(DOUBLE) dtime, djulianreftime
701  INTEGER itimestep, iactive
702  REAL fvalues(10) ! nValues
703  INTEGER*1 bactivity(10) ! activity
704  INTEGER i, status
705 
706  ! initialize the data
707  nvalues = 10
708  ntimes = 3
709  nactive = 8
710  dtime = 0.0
711 
712  ! 5th item in data set is always inactive, others active
713  do iactive = 1, nactive
714  bactivity(iactive) = 1
715  enddo
716  bactivity(6) = 0
717 
718  ! Create the scalar A dataset group
719  call xf_create_scalar_dataset(a_groupid, scalar_a_location, 'mg/L', &
720  ts_hours, none, xscalaraid, status)
721  if (status .LT. 0) then
722  ! close the dataset
723  call xf_close_group(xscalaraid, status)
724  call xf_close_group(xdsetsid, status)
725  call xf_close_file(xfileid, status)
726  return
727  endif
728 
729  ! Add in a reftime. This is a julian day for:
730  ! noon July 1, 2003
731  djulianreftime = 2452822.0;
732  call xf_write_reftime(xscalaraid, djulianreftime, status)
733  if (status < 0) then
734  call xf_close_group(xscalaraid, status)
735  call xf_close_group(xdsetsid, status)
736  call xf_close_file(xfileid, status)
737  endif
738 
739  ! Loop through timesteps adding them to the file
740  do itimestep = 1, ntimes
741  ! We will have an 0.5 hour timestep
742  dtime = itimestep * 0.5
743 
744  fvalues(1) = dtime
745  do i = 2, nvalues
746  fvalues(i) = fvalues(i-1)*2.5
747  end do
748 
749  ! write the dataset array values
750  call xf_write_scalar_timestep(xscalaraid, dtime, nvalues, fvalues, status)
751  if (status .GE. 0) then
752  ! write activity array
753  call xf_write_activity_timestep(xscalaraid, nactive, bactivity, status)
754  end if
755 
756  call tti_test_num_times(xscalaraid, itimestep, status)
757  enddo
758 
759  ! close the dataset
760  call xf_close_group(xscalaraid, status)
761  !call XF_CLOSE_GROUP(a_GroupID, status)
762  !call XF_CLOSE_FILE(a_FileID, status)
763 
764  return
765 END SUBROUTINE
766 ! ttWriteScalarAToMulti
767 ! --------------------------------------------------------------------------
768 ! FUNCTION ttReadVector2DAIndex
769 ! PURPOSE Read all timestep values for a particular index
770 ! NOTES
771 ! --------------------------------------------------------------------------
772 SUBROUTINE tt_read_vector2d_a_index (a_Filename, a_Index, error)
773 CHARACTER(LEN=*), INTENT(IN) :: a_filename
774 INTEGER, INTENT(IN) :: a_index
775 INTEGER, INTENT(OUT) :: error
776 INTEGER status
777 INTEGER xfileid, xdsetsid, xvector2da
778 INTEGER ntimesteps, i
779 REAL, ALLOCATABLE :: fvalues(:)
780 
781 xfileid = none
782 xdsetsid = none
783 xvector2da = none
784 
785  ! open the file
786 call xf_open_file(a_filename, .true., xfileid, status)
787 if (status < 0) then
788  error = -1
789  return
790 endif
791 
792  ! open the dataset group
793 call xf_open_group(xfileid, datasets_location, xdsetsid, status)
794 if (status >= 0) then
795  call xf_open_group(xdsetsid, vector2d_a_location, xvector2da, status)
796 endif
797 if (status < 0) then
798  error = status
799  return
800 endif
801 
802  ! Find out the number of timesteps in the file
803 call xf_get_dataset_num_times(xvector2da, ntimesteps, status)
804 if (status < 0) then
805  error = status
806  return
807 endif
808 
809 if (ntimesteps < 1) then
810  error = -1
811  return
812 endif
813 
814  ! Read the values for the index
815 allocate(fvalues(ntimesteps*2))
816 call xf_read_vector_values_at_index(xvector2da, a_index, 1, ntimesteps, 2, &
817  fvalues, status)
818 
819  ! output the data
820 WRITE(*,*) ''
821 WRITE(*,*) 'Reading vector 2D A slice at index: ', a_index
822 do i=1, ntimesteps
823  WRITE(*,*) fvalues(i*2-1), ' ', fvalues(i*2)
824 enddo
825 WRITE(*,*) ''
826 
827 deallocate(fvalues)
828 
829 error = status
830 return
831 
832 END SUBROUTINE
833 !ttReadVector2DAIndex
834 
835 ! --------------------------------------------------------------------------
836 ! FUNCTION ttWriteVector2D_A
837 ! PURPOSE Write scalar Dataset to an HDF5 File
838 ! NOTES This tests dynamic data sets, and activity
839 ! This dataset is dynamic concentrations (mg/L) with output times
840 ! in minutes.
841 ! Dataset is for a mesh and so nActive is the number of elements
842 ! which is not the same as the nValues which would be number of nodes
843 ! reads/writes a reference time in julian days
844 ! --------------------------------------------------------------------------
845 SUBROUTINE tt_write_vector2d_a (a_Filename, a_Compression, error)
846  CHARACTER(LEN=*), INTENT(IN) :: a_filename
847  INTEGER, INTENT(IN) :: a_compression
848  INTEGER, INTENT(OUT) :: error
849  INTEGER xfileid, xdsetsid, xvector2d_a, xcoordid
850  INTEGER nvalues, ntimes, ncomponents, nactive
851  REAL(DOUBLE) dtime
852  INTEGER itimestep, iactive
853  REAL, DIMENSION(2, 100) :: fvalues ! nComponents, nValues
854  INTEGER*1 bactivity(100) ! activity
855  INTEGER i, j, status
856  INTEGER ihpgnzone
857 
858  ! initialize the data
859  ncomponents = 2
860  nvalues = 100
861  ntimes = 6
862  nactive = 75
863  dtime = 0.0
864 
865  ! 5th item in data set is always inactive, others active
866  bactivity(1) = 0
867  do iactive = 2, nactive
868  if (mod(iactive-1, 3) == 0) then
869  bactivity(iactive) = 0
870  else
871  bactivity(iactive) = 1
872  endif
873  enddo
874 
875  ! create the file
876  call xf_create_file(a_filename, .true., xfileid, error)
877  if (error .LT. 0) then
878  ! close the dataset
879  call xf_close_file(xfileid, error)
880  return
881  endif
882 
883  ! create the group where we will put all the datasets
884  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
885  if (status < 0) then
886  call xf_close_file(xfileid, error)
887  error = -1
888  return
889  endif
890 
891  ! Create the vector dataset group
892  call xf_create_vector_dataset(xdsetsid, vector2d_a_location, 'ft/s', &
893  ts_seconds, a_compression, xvector2d_a, status)
894  if (status .LT. 0) then
895  ! close the dataset
896  call xf_close_group(xvector2d_a, error)
897  call xf_close_group(xdsetsid, error)
898  call xf_close_file(xfileid, error)
899  error = status
900  return
901  endif
902 
903  ! Loop through timesteps adding them to the file
904  do itimestep = 1, ntimes
905  ! We will have an 0.5 hour timestep
906  dtime = itimestep * 0.5
907 
908  do i = 1, nvalues
909  do j = 1, ncomponents
910  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
911  end do
912  end do
913 
914  ! write the dataset array values
915  call xf_write_vector_timestep(xvector2d_a, dtime, nvalues, ncomponents, &
916  fvalues, error)
917  if (error .GE. 0) then
918  ! write activity array
919  call xf_write_activity_timestep(xvector2d_a, nactive, bactivity, error)
920  end if
921 
922  call tti_test_num_times(xvector2d_a, itimestep, error)
923  enddo
924 
925  ! Write Coordinate file - for Vector2D_A, we will set the coordinate system
926  ! to be Geographic HPGN, with HPGN settings written to the file.
927  call xf_create_coordinate_group(xfileid, xcoordid, status)
928  if (status < 0) then
929  call xf_close_group(xvector2d_a, error)
930  call xf_close_group(xdsetsid, error)
931  call xf_close_file(xfileid, error)
932  error = -1
933  return
934  endif
935 
936  ! set HPGN info for test
937  ihpgnzone = 29 ! Utah
938 
939  call xf_set_horiz_datum(xcoordid, horiz_datum_geographic_hpgn, error)
940  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
941  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
942  call xf_set_vert_units(xcoordid, coord_units_meters, error)
943 
944  ! write additional information
945  call xf_set_hpgn_area(xcoordid, ihpgnzone, error)
946 
947  call xf_close_group(xcoordid, error)
948  xcoordid = 0
949 
950  ! close the dataset
951  call xf_close_group(xvector2d_a, error)
952  call xf_close_group(xdsetsid, error)
953  call xf_close_file(xfileid, error)
954 
955  return
956 END SUBROUTINE
957 ! ttWriteVector2D_A
958 
959 ! --------------------------------------------------------------------------
960 ! FUNCTION TT_WRITE_VECTOR2D_B
961 ! PURPOSE Write scalar Dataset to an HDF5 File
962 ! NOTES This tests dynamic data sets, and activity
963 ! This dataset is dynamic concentrations (mg/L) with output times
964 ! in minutes.
965 ! Dataset is for a mesh and so nActive is the number of elements
966 ! which is not the same as the nValues which would be number of nodes
967 ! reads/writes a reference time in julian days
968 ! --------------------------------------------------------------------------
969 SUBROUTINE tt_write_vector2d_b (a_Filename, a_Compression, a_Overwrite, error)
970  CHARACTER(LEN=*), INTENT(IN) :: a_filename
971  INTEGER, INTENT(IN) :: a_compression
972  LOGICAL, INTENT(IN) :: a_overwrite
973  INTEGER, INTENT(OUT) :: error
974  INTEGER xfileid, xdsetsid, xvector2d_b, xcoordid
975  INTEGER nvalues, ntimes, ncomponents, nactive
976  REAL(DOUBLE) dtime
977  INTEGER itimestep, iactive
978  REAL, DIMENSION(2, 100) :: fvalues
979  INTEGER*1 bactivity(100)
980  INTEGER i, j, status
981 
982  ! initialize the data
983  ncomponents = 2
984  nvalues = 100
985  ntimes = 6
986  nactive = 75
987  dtime = 0.0
988 
989  ! 5th item in data set is always inactive, others active
990  bactivity(1) = 0
991  do iactive = 2, nactive
992  if (mod(iactive-1, 3) == 0) then
993  bactivity(iactive) = 0
994  else
995  bactivity(iactive) = 1
996  endif
997  enddo
998 
999  if (a_overwrite) then
1000  ! open the already-existing file
1001  call xf_open_file(a_filename, .false., xfileid, status)
1002  if (status < 0) then
1003  error = -1
1004  return
1005  endif
1006  ! open the group where we have all the datasets
1007  call xf_open_group(xfileid, datasets_location, xdsetsid, status)
1008  if (status < 0) then
1009  call xf_close_file(xfileid, error)
1010  error = -1
1011  return
1012  endif
1013  else
1014  ! create the file
1015  call xf_create_file(a_filename, .true., xfileid, error)
1016  if (error .LT. 0) then
1017  ! close the dataset
1018  call xf_close_file(xfileid, error)
1019  return
1020  endif
1021 
1022  ! create the group where we will put all the datasets
1023  call xf_create_generic_group(xfileid, datasets_location, xdsetsid, status)
1024  if (status < 0) then
1025  call xf_close_file(xfileid, error)
1026  error = -1
1027  return
1028  endif
1029  endif
1030 
1031  ! Create/Overwrite the vector dataset group
1032  call xf_create_vector_dataset(xdsetsid, vector2d_b_location, 'ft/s', &
1033  ts_seconds, a_compression, xvector2d_b, status)
1034  if (status .LT. 0) then
1035  ! close the dataset
1036  call xf_close_group(xvector2d_b, error)
1037  call xf_close_group(xdsetsid, error)
1038  call xf_close_file(xfileid, error)
1039  error = status
1040  return
1041  endif
1042 
1043  if (.NOT. a_overwrite) then
1044  ! Loop through timesteps adding them to the file
1045  do itimestep = 1, ntimes
1046  ! We will have an 0.5 hour timestep
1047  dtime = itimestep * 0.5
1048  do i = 1, nvalues
1049  do j = 1, ncomponents
1050  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
1051  end do
1052  end do
1053  ! write the dataset array values
1054  call xf_write_vector_timestep(xvector2d_b, dtime, nvalues, ncomponents, &
1055  fvalues, error)
1056  if (error .GE. 0) then
1057  ! write activity array
1058  call xf_write_activity_timestep(xvector2d_b, nactive, bactivity, error)
1059  end if
1060  if (error < 0) then
1061  call xf_close_group(xvector2d_b, error)
1062  call xf_close_group(xdsetsid, error)
1063  call xf_close_file(xfileid, error)
1064  endif
1065  enddo
1066  else
1067  ! Loop through timesteps adding them to the file
1068  do itimestep = 1, ntimes
1069  ! We will have an 1.5 hour timestep
1070  dtime = itimestep * 1.5
1071  do i = 1, nvalues
1072  do j = 1, ncomponents
1073  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
1074  end do
1075  end do
1076  ! write the dataset array values
1077  call xf_write_vector_timestep(xvector2d_b, dtime, nvalues, ncomponents, &
1078  fvalues, error)
1079  if (error .GE. 0) then
1080  ! write activity array
1081  call xf_write_activity_timestep(xvector2d_b, nactive, bactivity, error)
1082  end if
1083  if (error < 0) then
1084  call xf_close_group(xvector2d_b, error)
1085  call xf_close_group(xdsetsid, error)
1086  call xf_close_file(xfileid, error)
1087  endif
1088  enddo
1089  endif
1090 
1091  if (.NOT. a_overwrite) then
1092  ! Write Coordinate file - for ScalarB, we will set the coordinate system
1093  ! to be UTM, with UTM Zone settings written to the file.
1094  call xf_create_coordinate_group(xfileid, xcoordid, status)
1095  if (status < 0) then
1096  call xf_close_group(xvector2d_b, error)
1097  call xf_close_group(xdsetsid, error)
1098  call xf_close_file(xfileid, error)
1099  error = -1
1100  return
1101  endif
1102 
1103  ! write the coordinate data to the file
1104  call xf_set_horiz_datum(xcoordid, horiz_datum_utm, error)
1105  call xf_set_horiz_units(xcoordid, coord_units_meters, error)
1106  call xf_set_vert_datum(xcoordid, vert_datum_local, error)
1107  call xf_set_vert_units(xcoordid, coord_units_meters, error)
1108 
1109  ! write additional information - we'll use the max UTM zone for the test
1110  call xf_set_utm_zone(xcoordid, utm_zone_max, error)
1111 
1112  call xf_close_group(xcoordid, error)
1113  xcoordid = 0
1114  endif
1115 
1116  ! close the dataset
1117  call xf_close_group(xvector2d_b, error)
1118  call xf_close_group(xdsetsid, error)
1119  call xf_close_file(xfileid, error)
1120 
1121  return
1122 END SUBROUTINE
1123 ! ttWriteVector2D_B
1124 
1125 ! --------------------------------------------------------------------------
1126 ! FUNCTION ttWriteVector2D_AToMulti
1127 ! PURPOSE Write scalar Dataset to an HDF5 File
1128 ! NOTES This tests dynamic data sets, and activity
1129 ! This dataset is dynamic concentrations (mg/L) with output times
1130 ! in minutes.
1131 ! Dataset is for a mesh and so nActive is the number of elements
1132 ! which is not the same as the nValues which would be number of nodes
1133 ! reads/writes a reference time in julian days
1134 ! --------------------------------------------------------------------------
1135 SUBROUTINE tt_write_vector2d_a_to_multi (a_FileID, a_GroupID, status)
1136  INTEGER xvector2d_a
1137  INTEGER a_fileid, a_groupid
1138  INTEGER nvalues, ntimes, ncomponents, nactive
1139  REAL(DOUBLE) dtime
1140  INTEGER itimestep, iactive
1141  REAL, DIMENSION(2, 100) :: fvalues ! nComponents, nValues
1142  INTEGER*1 bactivity(100) ! activity
1143  INTEGER i, j, status
1144 
1145  ! initialize the data
1146  ncomponents = 2
1147  nvalues = 100
1148  ntimes = 6
1149  nactive = 75
1150  dtime = 0.0
1151 
1152  ! 5th item in data set is always inactive, others active
1153  bactivity(1) = 0
1154  do iactive = 2, nactive
1155  if (mod(iactive-1, 3) == 0) then
1156  bactivity(iactive) = 0
1157  else
1158  bactivity(iactive) = 1
1159  endif
1160  enddo
1161 
1162  ! Create the vector dataset group
1163  call xf_create_vector_dataset(a_groupid, vector2d_a_location, 'ft/s', &
1164  ts_seconds, none, xvector2d_a, status)
1165  if (status .LT. 0) then
1166  ! close the dataset
1167  call xf_close_group(xvector2d_a, status)
1168  call xf_close_group(a_groupid, status)
1169  call xf_close_file(a_fileid, status)
1170  return
1171  endif
1172 
1173  ! Loop through timesteps adding them to the file
1174  do itimestep = 1, ntimes
1175  ! We will have an 0.5 hour timestep
1176  dtime = itimestep * 0.5
1177 
1178  do i = 1, nvalues
1179  do j = 1, ncomponents
1180  fvalues(j,i) = ((i-1)*ncomponents + j)*dtime
1181  end do
1182  end do
1183 
1184  ! write the dataset array values
1185  call xf_write_vector_timestep(xvector2d_a, dtime, nvalues, ncomponents, &
1186  fvalues, status)
1187  if (status .GE. 0) then
1188  ! write activity array
1189  call xf_write_activity_timestep(xvector2d_a, nactive, bactivity, status)
1190  end if
1191 
1192  call tti_test_num_times(xvector2d_a, itimestep, status)
1193  enddo
1194 
1195  ! close the dataset
1196  call xf_close_group(xvector2d_a, status)
1197  return
1198 END SUBROUTINE
1199 ! ttWriteVector2D_AToMulti
1200 ! --------------------------------------------------------------------------
1201 ! FUNCTION ttiReadScalar
1202 ! PURPOSE Read a scalar from an XMDF file and output information to
1203 ! to a text file
1204 ! NOTES
1205 ! --------------------------------------------------------------------------
1206 SUBROUTINE tti_read_scalar (a_xScalarId, FileUnit, error)
1207  INTEGER, INTENT(IN) :: a_xscalarid
1208  INTEGER, INTENT(IN) :: fileunit
1209  INTEGER, INTENT(OUT) :: error
1210  INTEGER ntimes, nvalues, nactive
1211  LOGICAL*2 busereftime
1212  INTEGER itime
1213  CHARACTER(LEN=100) timeunits
1214  REAL(DOUBLE), ALLOCATABLE :: times(:)
1215  REAL, ALLOCATABLE :: values(:), minimums(:), maximums(:)
1216  INTEGER, ALLOCATABLE :: active(:)
1217  REAL(DOUBLE) reftime
1218 ntimes = none
1219 nvalues = none
1220 nactive = none
1221 
1222  ! read the time units
1223  call xf_get_dataset_time_units(a_xscalarid, timeunits, error)
1224  if (error < 0) return
1225 
1226  WRITE(fileunit,*) 'Time units: ', timeunits(1:len_trim(timeunits))
1227 
1228  ! see if we are using a reftime
1229  call xf_use_reftime(a_xscalarid, busereftime, error)
1230  if (error < 0) then
1231  return
1232  endif
1233  if (busereftime) then
1234  call xf_read_reftime(a_xscalarid, reftime, error)
1235  if (error < 0) then
1236  return
1237  endif
1238  !WRITE(FileUnit,*) 'Reftime: ', Reftime
1239  endif
1240 
1241  ! read in the number of values and number of active values
1242  call xf_get_dataset_numvals(a_xscalarid, nvalues, error)
1243  if (error .GE. 0) then
1244  call xf_get_dataset_numactive(a_xscalarid, nactive, error)
1245  endif
1246  if (error .LT. 0) return
1247 
1248  if (nvalues <= 0) then
1249  WRITE(fileunit, *) 'No data to read in.'
1250  error = -1
1251  return
1252  endif
1253 
1254  ! read in the number of times
1255  call xf_get_dataset_num_times(a_xscalarid, ntimes, error)
1256  if (error < 0) then
1257  return
1258  endif
1259 
1260  ! Read in the individual time values
1261  allocate(times(ntimes))
1262 
1263  call xf_get_dataset_times(a_xscalarid, ntimes, times, error)
1264  if (error < 0) return
1265 
1266  ! Read in the minimum and maximum values
1267  allocate(minimums(ntimes))
1268  allocate(maximums(ntimes))
1269 
1270  call xf_get_dataset_mins(a_xscalarid, ntimes, minimums, error)
1271  if (error >= 0) then
1272  call xf_get_dataset_maxs(a_xscalarid, ntimes, maximums, error)
1273  endif
1274  if (error < 0) then
1275  deallocate(times)
1276  deallocate(minimums)
1277  deallocate(maximums)
1278  return
1279  endif
1280 
1281  allocate(values(nvalues))
1282  if (nactive .GT. 0) then
1283  allocate(active(nactive))
1284  endif
1285 
1286  WRITE(fileunit,*) 'Number Timesteps: ', ntimes
1287  WRITE(fileunit,*) 'Number Values: ', nvalues
1288  WRITE(fileunit,*) 'Number Active: ', nactive
1289  WRITE(fileunit,*) ''
1290 
1291  ! loop through the timesteps, read the values and active values and write
1292  ! them to the text file
1293  do itime = 1, ntimes
1294  call xf_read_scalar_values_timestep(a_xscalarid, itime, nvalues, values, error)
1295  if (error >= 0 .AND. nactive > 0) then
1296  call xf_read_activity_timestep(a_xscalarid, itime, nactive, active, error)
1297  endif
1298 
1299  ! Write the time, min, max, values and active values to the text output
1300  ! file.
1301  WRITE(fileunit,*) 'Timestep at ', times(itime)
1302  WRITE(fileunit,*) 'Min: ', minimums(itime)
1303  WRITE(fileunit,*) 'Max: ', maximums(itime)
1304 
1305  WRITE(fileunit,*) 'Values:'
1306  WRITE(fileunit,*) values(1:nvalues)
1307  WRITE(fileunit,*) ''
1308 
1309  WRITE(fileunit,*) 'Activity:'
1310  WRITE(fileunit,*) active(1:nactive)
1311  WRITE(fileunit,*) ''
1312  end do
1313 
1314  if (allocated(times)) then
1315  deallocate(times)
1316  endif
1317 
1318  if (allocated(minimums)) then
1319  deallocate(minimums)
1320  endif
1321 
1322  if (allocated(maximums)) then
1323  deallocate(maximums)
1324  endif
1325 
1326  if (allocated(values)) then
1327  deallocate(values)
1328  endif
1329 
1330  if (allocated(active)) then
1331  deallocate(active)
1332  endif
1333 
1334  return
1335 END SUBROUTINE
1336 ! ttiReadScalar
1337 
1338 ! --------------------------------------------------------------------------
1339 ! FUNCTION TTI_READ_VECTOR
1340 ! PURPOSE Read a vector from an XMDF file and output information to
1341 ! to a text file
1342 ! NOTES
1343 ! --------------------------------------------------------------------------
1344 SUBROUTINE tti_read_vector (a_xVectorId, FileUnit, error)
1345  INTEGER, INTENT(IN) :: a_xvectorid
1346  INTEGER, INTENT(IN) :: fileunit
1347  INTEGER, INTENT(OUT) :: error
1348  INTEGER ntimes, nvalues, nactive, ncomponents
1349  INTEGER itime, i
1350  LOGICAL*2 busereftime
1351  CHARACTER(LEN=100) timeunits
1352  REAL(DOUBLE), ALLOCATABLE :: times(:)
1353  REAL, ALLOCATABLE, DIMENSION (:, :) :: values
1354  REAL, ALLOCATABLE :: minimums(:), maximums(:)
1355  INTEGER, ALLOCATABLE :: active(:)
1356  REAL(DOUBLE) reftime
1357 
1358 ntimes = none
1359 nvalues = none
1360 nactive = none
1361 ncomponents = none
1362 
1363  ! read the time units
1364  call xf_get_dataset_time_units(a_xvectorid, timeunits, error)
1365  if (error < 0) return
1366 
1367  WRITE(fileunit,*) 'Time units: ', timeunits(1:len_trim(timeunits))
1368 
1369  ! see if we are using a reftime
1370  call xf_use_reftime(a_xvectorid, busereftime, error)
1371  if (error < 0) then
1372  return
1373  endif
1374  if (busereftime) then
1375  call xf_read_reftime(a_xvectorid, reftime, error)
1376  if (error < 0) then
1377  return
1378  endif
1379  !WRITE(FileUnit,*) 'Reftime: ', Reftime
1380  endif
1381 
1382  ! read in the number of values and number of active values
1383  call xf_get_dataset_numvals(a_xvectorid, nvalues, error)
1384  if (error .GE. 0) then
1385  call xf_get_dataset_numcomponents(a_xvectorid, ncomponents, error)
1386  if (error .GE. 0) then
1387  call xf_get_dataset_numactive(a_xvectorid, nactive, error)
1388  endif
1389  endif
1390  if (error .LT. 0) return
1391 
1392  if (nvalues <= 0) then
1393  WRITE(fileunit, *) 'No data to read in.'
1394  error = -1
1395  return
1396  endif
1397 
1398  ! read in the number of times
1399  call xf_get_dataset_num_times(a_xvectorid, ntimes, error)
1400  if (error < 0) then
1401  return
1402  endif
1403 
1404  ! Read in the individual time values
1405  allocate(times(ntimes))
1406 
1407  call xf_get_dataset_times(a_xvectorid, ntimes, times, error)
1408  if (error < 0) return
1409 
1410  ! Read in the minimum and maximum values
1411  allocate(minimums(ntimes))
1412  allocate(maximums(ntimes))
1413 
1414  call xf_get_dataset_mins(a_xvectorid, ntimes, minimums, error)
1415  if (error >= 0) then
1416  call xf_get_dataset_maxs(a_xvectorid, ntimes, maximums, error)
1417  endif
1418  if (error < 0) then
1419  deallocate(times)
1420  deallocate(minimums)
1421  deallocate(maximums)
1422  return
1423  endif
1424 
1425  allocate(values(ncomponents, nvalues))
1426  if (nactive .GT. 0) then
1427  allocate(active(nactive))
1428  endif
1429 
1430  WRITE(fileunit,*) 'Number Timesteps: ', ntimes
1431  WRITE(fileunit,*) 'Number Values: ', nvalues
1432  WRITE(fileunit,*) 'Number Components: ', ncomponents
1433  WRITE(fileunit,*) 'Number Active: ', nactive
1434 
1435  ! loop through the timesteps, read the values and active values and write
1436  ! them to the text file
1437  do itime = 1, ntimes
1438  call xf_read_vector_values_timestep(a_xvectorid, itime, nvalues, &
1439  ncomponents, values, error)
1440  if (error >= 0 .AND. nactive > 0) then
1441  call xf_read_activity_timestep(a_xvectorid, itime, nactive, active, error)
1442  endif
1443 
1444  ! Write the time, min, max, values and active values to the text output
1445  ! file.
1446  WRITE(fileunit,*) ''
1447  WRITE(fileunit,*) 'Timestep at ', times(itime)
1448  WRITE(fileunit,*) 'Min: ', minimums(itime)
1449  WRITE(fileunit,*) 'Max: ', maximums(itime)
1450 
1451  WRITE(fileunit,*) 'Values:'
1452  do i=1, nvalues
1453  WRITE(fileunit,*) values(1:ncomponents,i:i)
1454  enddo
1455  WRITE(fileunit,*) ''
1456 
1457  WRITE(fileunit,*) 'Activity:'
1458  WRITE(fileunit,*) active(1:nactive)
1459  WRITE(fileunit,*) ''
1460  WRITE(fileunit,*) ''
1461 
1462  end do
1463 
1464  if (allocated(times)) then
1465  deallocate(times)
1466  endif
1467 
1468  if (allocated(minimums)) then
1469  deallocate(minimums)
1470  endif
1471 
1472  if (allocated(maximums)) then
1473  deallocate(maximums)
1474  endif
1475 
1476  if (allocated(values)) then
1477  deallocate(values)
1478  endif
1479 
1480  if (allocated(active)) then
1481  deallocate(active)
1482  endif
1483 
1484  return
1485 END SUBROUTINE
1486 ! ttiReadVector
1487 
1488 END MODULE testtimestep

TestTimestep.f90 tests timesteps

1 MODULE testdefs
2 
3 INTEGER, PARAMETER :: numtimes1 = 5
4 INTEGER, PARAMETER :: numvalues1 = 5
5 INTEGER, PARAMETER :: numactive1 = 3
6 INTEGER, PARAMETER :: numtimesadd = 1
7 
8 CHARACTER(LEN=*), PARAMETER :: xmdf_version_out_f = 'XMDF_Version_f.txt'
9 CHARACTER(LEN=*), PARAMETER :: mesh_a_file_f = 'mesh_a_file_f.h5'
10 CHARACTER(LEN=*), PARAMETER :: mesh_b_file_f = 'mesh_b_file_f.h5'
11 CHARACTER(LEN=*), PARAMETER :: mesh_a_out_f = 'mesh_a_file_f.txt'
12 CHARACTER(LEN=*), PARAMETER :: mesh_b_out_f = 'mesh_b_file_f.txt'
13 
14 CHARACTER(LEN=*), PARAMETER :: grid_cart2d_a_file_f = 'grid_cart2d_a_file_f.h5'
15 CHARACTER(LEN=*), PARAMETER :: grid_curv2d_a_file_f = 'grid_curv2d_a_file_f.h5'
16 CHARACTER(LEN=*), PARAMETER :: grid_cart3d_a_file_f = 'grid_cart3d_a_file_f.h5'
17 
18 CHARACTER(LEN=*), PARAMETER :: grid_cart2d_a_out_f = 'grid_cart2d_a_out_f.txt'
19 CHARACTER(LEN=*), PARAMETER :: grid_curv2d_a_out_f = 'grid_curv2d_a_out_f.txt'
20 CHARACTER(LEN=*), PARAMETER :: grid_cart3d_a_out_f = 'grid_cart3d_a_out_f.txt'
21 
22 CHARACTER(LEN=*), PARAMETER :: multidataset_file_f = 'MultiDataSet_f.h5'
23 CHARACTER(LEN=*), PARAMETER :: multidataset_text_f = 'MultiDataSet_f.txt'
24 
25 CHARACTER(LEN=*), PARAMETER :: scalar_a_file_f = 'ScalarA_f.h5'
26 CHARACTER(LEN=*), PARAMETER :: scalar_a_text_f = 'ScalarA_f.txt'
27 CHARACTER(LEN=*), PARAMETER :: scalar_a_pieces_file_f = 'ScalarA_Pieces_f.h5'
28 CHARACTER(LEN=*), PARAMETER :: scalar_a_edited_file_f = 'ScalarA_edited_f.h5'
29 CHARACTER(LEN=*), PARAMETER :: scalar_a_edited_text_f = 'ScalarA_edited_f.txt'
30 CHARACTER(LEN=*), PARAMETER :: scalar_a_pieces_alt_file_f = 'ScalarA_Pieces_alt_f.h5'
31 CHARACTER(LEN=*), PARAMETER :: scalar_b_file_f = 'ScalarB_f.h5'
32 CHARACTER(LEN=*), PARAMETER :: scalar_b_text_f = 'ScalarB_f.txt'
33 CHARACTER(LEN=*), PARAMETER :: vector2d_a_file_f = 'Vector2D_A_f.h5'
34 CHARACTER(LEN=*), PARAMETER :: vector2d_a_text_f = 'Vector2D_A_f.txt'
35 CHARACTER(LEN=*), PARAMETER :: vector2d_a_pieces_file_f = 'Vector2D_A_Pieces_f.h5'
36 CHARACTER(LEN=*), PARAMETER :: vector2d_b_file_f = 'Vector2D_B_f.h5'
37 CHARACTER(LEN=*), PARAMETER :: vector2d_b_text_f = 'Vector2D_B_f.txt'
38 
39 CHARACTER(LEN=*), PARAMETER :: mesh_a_file_c = 'mesh_a_file_c.h5'
40 CHARACTER(LEN=*), PARAMETER :: mesh_b_file_c = 'mesh_b_file_c.h5'
41 CHARACTER(LEN=*), PARAMETER :: mesh_a_out_cf = 'mesh_a_file_cf.txt'
42 CHARACTER(LEN=*), PARAMETER :: mesh_b_out_cf = 'mesh_b_file_cf.txt'
43 
44 CHARACTER(LEN=*), PARAMETER :: grid_cart2d_a_file_c = 'grid_cart2d_a_file_c.h5'
45 CHARACTER(LEN=*), PARAMETER :: grid_curv2d_a_file_c = 'grid_curv2d_a_file_c.h5'
46 CHARACTER(LEN=*), PARAMETER :: grid_cart3d_a_file_c = 'grid_cart3d_a_file_c.h5'
47 
48 CHARACTER(LEN=*), PARAMETER :: grid_cart2d_a_out_cf = 'grid_cart2d_a_out_cf.txt'
49 CHARACTER(LEN=*), PARAMETER :: grid_curv2d_a_out_cf = 'grid_curv2d_a_out_cf.txt'
50 CHARACTER(LEN=*), PARAMETER :: grid_cart3d_a_out_cf = 'grid_cart3d_a_out_cf.txt'
51 
52 CHARACTER(LEN=*), PARAMETER :: scalar_a_file_c = 'ScalarA_c.h5'
53 CHARACTER(LEN=*), PARAMETER :: scalar_a_text_cf = 'ScalarA_cf.txt'
54 CHARACTER(LEN=*), PARAMETER :: scalar_b_file_c = 'ScalarB_c.h5'
55 CHARACTER(LEN=*), PARAMETER :: scalar_b_text_cf = 'ScalarB_cf.txt'
56 CHARACTER(LEN=*), PARAMETER :: vector2d_a_file_c = 'Vector2D_A_c.h5'
57 CHARACTER(LEN=*), PARAMETER :: vector2d_a_text_cf = 'Vector2D_A_cf.txt'
58 CHARACTER(LEN=*), PARAMETER :: vector2d_b_file_c = 'Vector2D_B_c.h5'
59 CHARACTER(LEN=*), PARAMETER :: vector2d_b_text_cf = 'Vector2D_B_cf.txt'
60 
61 CHARACTER(LEN=*), PARAMETER :: calendar_out_f = 'Calendar_f.txt'
62 
63 CHARACTER(LEN=*), PARAMETER :: geompath_a_file_f = 'Geompath_a_file_f.h5'
64 CHARACTER(LEN=*), PARAMETER :: geompath_a_file_f_out = 'Geompath_a_file_f_out.txt'
65 
66 CHARACTER(LEN=*), PARAMETER :: tt_multidataset_file_f = 'TT_MultiDataSet_f.h5'
67 CHARACTER(LEN=*), PARAMETER :: tt_scalar_a_file_f = 'TT_ScalarA_f.h5'
68 CHARACTER(LEN=*), PARAMETER :: tt_scalar_a_text_f = 'TT_ScalarA_f.txt'
69 CHARACTER(LEN=*), PARAMETER :: tt_vector2d_a_file_f = 'TT_Vector2D_A_f.h5'
70 
71  ! Overwrite options in the function xfSetupToWriteDatasets
72 !INTEGER,PARAMETER :: XF_OVERWRITE_CLEAR_FILE = 1
73 !INTEGER,PARAMETER :: XF_OVERWRITE_CLEAR_DATASET_GROUP = 2
74 !INTEGER,PARAMETER :: XF_OVERWRITE_NONE = 3
75 
76 END MODULE testdefs
77 
78 MODULE testsmodule
79 
80 USE xmdf
81 USE xmdfdefs
82 USE testtimestep
83 USE testdatasets
84 USE testmesh
85 USE testgrid
86 USE testdefs
87 USE testgeompaths
88 
89 CONTAINS
90 
91 !**************************
92 !-----------------------------------------------------------------------------
93 ! SUBROUTINE TXI_TEST_TIMESTEPS
94 ! PURPOSE test to see if the code can read timestepfiles
95 ! NOTES
96 !------------------------------------------------------------------------------
97 SUBROUTINE txi_test_timesteps(error)
98 INTEGER, INTENT(OUT) :: error
99 INTEGER status, compression
100 INTEGER multifileid, multigroupid
101 INTEGER numopen
102 
103 CHARACTER(LEN=37) sdoguid
104 
105 sdoguid = '73289C80-6235-4fdc-9649-49E4F5AEB676'
106 
107 status = 1
108 compression = none
109 ! 'Path' should be able to be blank, but now it creates errors if it's blank
110  call xf_setup_to_write_datasets(tt_multidataset_file_f, 'Multidatasets','', &
111  sdoguid, xf_overwrite_clear_file, multifileid, multigroupid, error)
112 
113  ! Write coordinates to multidatasets
114 call tt_write_coords_to_multi(multifileid, error)
115 
116  ! Write scalar A and Vector A to multidatasets.
117 call tt_write_scalar_a_to_multi(multigroupid, error)
118 
119 call tt_write_vector2d_a_to_multi(multifileid, multigroupid, error)
120 
121 call xf_close_group(multigroupid, status)
122 call xf_close_file(multifileid, status)
123 
124 WRITE(*,*) 'Done writing multiple datasets...'
125 
126  ! scalar datasets
127 call tt_write_scalar_a(tt_scalar_a_file_f, compression, status)
128 if (status < 0) then
129  error = status
130  return
131 endif
132 
133 WRITE(*,*) 'Done writing scalar datasets...'
134 
135  ! vector datasets
136 call tt_write_vector2d_a(tt_vector2d_a_file_f, compression, status)
137 if (status < 0) then
138  WRITE(*,*) 'Error writing dataset vector2D_A.'
139  error = status
140  return
141 endif
142 
143 WRITE(*,*) 'Done writing vector datasets...'
144 
145 WRITE(*,*) 'Write edited scalar datasets...'
146 call td_edit_scalar_a_values(scalar_a_edited_file_f, compression, status);
147 if (status < 0) then
148  error = status
149  return
150 endif
151 
152 WRITE(*,*) 'Done writing datasets...'
153 
154  ! Read the files back in
155 call txi_read_x_format_file(tt_scalar_a_file_f, scalar_a_text_f, status)
156 if (status < 0) then
157  WRITE(*,*) 'Error reading SCALAR A File (see TXI_READ_X_FORMAT_FILE)'
158  error = status
159  return
160 endif
161 
162 call txi_read_x_format_file(scalar_a_edited_file_f, scalar_a_edited_text_f, status)
163 if (status < 0) then
164  WRITE(*,*) 'Error reading SCALAR A Edited File (see TXI_READ_X_FORMAT_FILE)'
165  error = status
166  return
167 endif
168 
169 call txi_read_x_format_file(tt_vector2d_a_file_f, vector2d_a_text_f, status)
170 if (status < 0) then
171  WRITE(*,*) 'Error reading VECTOR A Format File'
172  error = status
173  return
174 endif
175 
176 call txi_read_x_format_file(tt_multidataset_file_f, multidataset_text_f, status)
177 if (status < 0) then
178  WRITE(*,*) 'Error reading Multidataset File (see TXI_READ_X_FORMAT_FILE)'
179  error = status
180  return
181 endif
182 
183 WRITE(*,*) 'Done reading datasets...'
184 
185 call xf_get_num_open_identifiers(h5f_obj_all_f, numopen, error)
186 
187 call xfi_close_open_identifiers(h5f_obj_all_f, error)
188 
189 call xf_setup_to_write_datasets(tt_multidataset_file_f, 'Multidatasets','', &
190  sdoguid, xf_overwrite_clear_dataset_grp, multifileid, multigroupid, &
191  error)
192 
193 call tt_write_scalar_a_to_multi(multigroupid, error)
194 
195 call xf_setup_to_write_datasets(tt_multidataset_file_f, 'Multidatasets','', &
196  sdoguid, xf_overwrite_none, multifileid, multigroupid, error)
197 
198 call tt_write_vector2d_a_to_multi(multifileid, multigroupid, error)
199 
200  ! Test reading information at index for multiple timesteps
201 call tt_read_scalar_a_index(tt_scalar_a_file_f, 4, status)
202 if (status < 0) then
203  error = status
204  return
205 endif
206 
207 WRITE(*,*) 'Done reading scalar data at index.'
208 
209 call tt_read_vector2d_a_index(tt_vector2d_a_file_f, 6, status)
210 if (status < 0) then
211  error = status
212  return
213 endif
214 
215 WRITE(*,*) 'Done reading vector data at index.'
216 
217 call tt_read_activity_scalar_a_index(tt_scalar_a_file_f, 6, status)
218 if (status < 0) then
219  error = status
220  return
221 endif
222 
223 error = status
224 return
225 
226 END SUBROUTINE
227 
228 
229 !**************************
230 !-----------------------------------------------------------------------------
231 ! SUBROUTINE TXI_TEST_DATASETS
232 ! PURPOSE test to see if the code can read datasetfiles
233 ! NOTES
234 !------------------------------------------------------------------------------
235 SUBROUTINE txi_test_datasets(error)
236 INTEGER, INTENT(OUT) :: error
237 INTEGER status, compression
238 INTEGER multifileid, multigroupid
239 INTEGER numopen
240 INTEGER, PARAMETER :: nindices = 3
241 INTEGER, DIMENSION(nIndices) :: indices
242 
243 CHARACTER(LEN=37) sdoguid
244 
245 sdoguid = '73289C80-6235-4fdc-9649-49E4F5AEB676'
246 
247 status = 1
248 compression = none
249 ! 'Path' should be able to be blank, but now it creates errors if it's blank
250 call xf_setup_to_write_datasets(multidataset_file_f, 'Multidatasets','', &
251  sdoguid, xf_overwrite_clear_file, multifileid, multigroupid, error)
252 
253  ! Write coordinates to multidatasets
254 call td_write_coords_to_multi(multifileid, error)
255 
256  ! Write scalar A and Vector A to multidatasets.
257 call td_write_scalar_a_to_multi(multigroupid, error)
258 
259 call td_write_vector2d_a_to_multi(multifileid, multigroupid, error)
260 
261 call xf_close_group(multigroupid, status)
262 call xf_close_file(multifileid, status)
263 
264 WRITE(*,*) 'Done writing multiple datasets...'
265 
266  ! scalar datasets
267 call td_write_scalar_a(scalar_a_file_f, compression, status)
268 if (status < 0) then
269  error = status
270  return
271 endif
272 
273 call td_write_scalar_a_pieces(scalar_a_pieces_file_f, compression, status)
274 if (status < 0) then
275  error = status
276  return
277 endif
278 
279 call td_write_scalar_a_pieces_alt_min_max(scalar_a_pieces_alt_file_f, &
280  compression, status)
281 if (status < 0) then
282  error = status
283  return
284 endif
285 
286 WRITE(*,*) 'Done writing scalar datasets...'
287 
288  ! vector datasets
289 call td_write_vector2d_a(vector2d_a_file_f, compression, status)
290 if (status < 0) then
291  WRITE(*,*) 'Error writing dataset vector2D_A.'
292  error = status
293  return
294 endif
295 
296 call td_write_vector2d_a_pieces(vector2d_a_pieces_file_f, compression, status)
297 if (status < 0) then
298  WRITE(*,*) 'Error writing dataset vector2D_A.'
299  error = status
300  return
301 endif
302 
303 WRITE(*,*) 'Done writing vector datasets...'
304 
305 WRITE(*,*) 'Done writing datasets...'
306 
307  ! Read the files back in
308 call txi_read_x_format_file(scalar_a_file_f, scalar_a_text_f, status)
309 if (status < 0) then
310  WRITE(*,*) 'Error reading SCALAR A File (see TXI_READ_X_FORMAT_FILE)'
311  error = status
312  return
313 endif
314 
315 call txi_read_x_format_file(vector2d_a_file_f, vector2d_a_text_f, status)
316 if (status < 0) then
317  WRITE(*,*) 'Error reading VECTOR A Format File'
318  error = status
319  return
320 endif
321 
322 call txi_read_x_format_file(multidataset_file_f, multidataset_text_f, status)
323 if (status < 0) then
324  WRITE(*,*) 'Error reading Multidataset File (see TXI_READ_X_FORMAT_FILE)'
325  error = status
326  return
327 endif
328 
329 WRITE(*,*) 'Done reading datasets...'
330 
331 call xf_get_num_open_identifiers(h5f_obj_all_f, numopen, error)
332 
333 call xfi_close_open_identifiers(h5f_obj_all_f, error)
334 
335 call xf_setup_to_write_datasets(multidataset_file_f, 'Multidatasets','', &
336  sdoguid, xf_overwrite_clear_dataset_grp, multifileid, multigroupid, &
337  error)
338 
339 call td_write_scalar_a_to_multi(multigroupid, error)
340 
341 call xf_setup_to_write_datasets(multidataset_file_f, 'Multidatasets','', &
342  sdoguid, xf_overwrite_none, multifileid, multigroupid, error)
343 
344 call td_write_vector2d_a_to_multi(multifileid, multigroupid, error)
345 
346  ! Test reading information at index for multiple timesteps
347 call td_read_scalar_a_index(scalar_a_file_f, 4, status)
348 if (status < 0) then
349  error = status
350  return
351 endif
352 
353 WRITE(*,*) 'Done reading scalar data at index.'
354 
355 call td_read_vector2d_a_index(vector2d_a_file_f, 6, status)
356 if (status < 0) then
357  error = status
358  return
359 endif
360 
361 WRITE(*,*) 'Done reading vector data at index.'
362 
363 call td_read_activity_scalar_a_index(scalar_a_file_f, 6, status)
364 if (status < 0) then
365  error = status
366  return
367 endif
368 
369  ! Test reading information at multiple indices
370  indices(1) = 2;
371  indices(2) = 3;
372  indices(3) = 5;
373  CALL td_read_scalar_a_indices(scalar_a_file_f, nindices, indices, error)
374 
375 error = status
376 return
377 
378 END SUBROUTINE
379 
380 !------------------------------------------------------------------------------
381 ! FUNCTION TXI_TEST_OVERWRITE_DSETS
382 ! PURPOSE Check to see if already-written datasets can be overwritten
383 ! NOTES
384 !------------------------------------------------------------------------------
385 SUBROUTINE txi_test_overwrite_dsets(error)
386 INTEGER, INTENT(OUT) :: error
387 INTEGER status, compression
388 
389  status = 1
390  compression = none
391 
392  ! scalar datasets
393  call td_write_scalar_b(scalar_b_file_f, compression, .false., status)
394  if (status < 0) then
395  error = status
396  return
397  endif
398  !overwrite scalar datasets
399  call td_write_scalar_b(scalar_b_file_f, compression, .true., status)
400  if (status < 0) then
401  error = status
402  return
403  endif
404 
405  ! vector datasets
406  call td_write_vector2d_b(vector2d_b_file_f, compression, .false., status)
407  if (status < 0) then
408  WRITE(*,*) 'Error writing dataset vector2D_B.'
409  error = status
410  return
411  endif
412  ! overwrite vector datasets
413  call td_write_vector2d_b(vector2d_b_file_f, compression, .true., status)
414  if (status < 0) then
415  WRITE(*,*) 'Error writing dataset vector2D_B.'
416  error = status
417  return
418  endif
419 
420  ! Read the files back in
421  call txi_read_x_format_file(scalar_b_file_f, scalar_b_text_f, status)
422  if (status < 0) then
423  WRITE(*,*) 'Error reading SCALAR B File'
424  error = status
425  return
426  endif
427 
428  call txi_read_x_format_file(vector2d_b_file_f, vector2d_b_text_f, status)
429  if (status < 0) then
430  WRITE(*,*) 'Error reading VECTOR B Format File'
431  error = status
432  return
433  endif
434 
435  error = status
436  return
437 
438 END SUBROUTINE
439 
440 !------------------------------------------------------------------------------
441 ! FUNCTION TXI_TEST_GRIDS
442 ! PURPOSE Test to see if we can read and write grids
443 ! NOTES
444 !------------------------------------------------------------------------------
445 SUBROUTINE txi_test_grids (error)
446 INTEGER, INTENT(OUT) :: error
447 INTEGER compression
448 
449 compression = none
450 
451 WRITE(*,*) ''
452 WRITE(*,*) ''
453 WRITE(*,*) 'Writing grid data.'
454 WRITE(*,*) ''
455 
456 call tg_write_test_grid_cart_2d(grid_cart2d_a_file_f, error)
457 if (error < 0) then
458  WRITE(*,*) 'Error writing grid Cartesian 2D A'
459 endif
460 WRITE(*,*) 'Finished writing grid Cartesian 2D A'
461 
462 call tg_write_test_grid_curv_2d(grid_curv2d_a_file_f, compression, error)
463 if (error < 0) then
464  WRITE(*,*) 'Error writing grid Curvilinear 2D A'
465 endif
466 WRITE(*,*) 'Finished writing grid Curvilinear 2D A'
467 
468 call tg_write_test_grid_cart_3d(grid_cart3d_a_file_f, compression, error)
469 if (error < 0) then
470  WRITE(*,*) 'Error writing grid Cartesian 3D A'
471 endif
472 WRITE(*,*) 'Finished writing grid Cartesian 3D A'
473  ! read the files back in
474 call txi_read_x_format_file(grid_cart2d_a_file_f, grid_cart2d_a_out_f, error)
475 if (error < 0) then
476  WRITE(*,*) 'Error reading grid Cartesian 2D A'
477 endif
478 WRITE(*,*) 'Finished reading grid Cartesian 2D A'
479 
480 call txi_read_x_format_file(grid_curv2d_a_file_f, grid_curv2d_a_out_f, error)
481 if (error < 0) then
482  WRITE(*,*) 'Error reading grid Curvilinear 2D A'
483 endif
484 WRITE(*,*) 'Finished reading grid Curvilinear 2D A'
485 
486 call txi_read_x_format_file(grid_cart3d_a_file_f, grid_cart3d_a_out_f, error)
487 if (error < 0) then
488  WRITE(*,*) 'Error reading grid Cartesian 3D A'
489 endif
490 WRITE(*,*) 'Finished reading grid Cartesian 3D A'
491 
492 END SUBROUTINE
493 
494 !**************************
495 ! ---------------------------------------------------------------------------
496 ! FUNCTION TXI_TEST_MESHS
497 ! PURPOSE test to see if we can read and write meshes
498 ! NOTES
499 ! ---------------------------------------------------------------------------
500 SUBROUTINE txi_test_meshs (error)
501 INTEGER, INTENT(OUT) :: error
502 INTEGER status
503 INTEGER compression
504 
505 status = 1
506 compression = none
507 
508 call tm_write_test_mesh_a(mesh_a_file_f, compression, status)
509 if (status < 0) then
510  WRITE(*,*) 'Error writing TestMeshA'
511  error = status
512  return
513 endif
514 
515 call tm_write_test_mesh_b(mesh_b_file_f, compression, status)
516 if (status < 0) then
517  WRITE(*,*) 'Error writing TestMeshB'
518  error = status
519  return
520 endif
521 
522 WRITE(*,*) 'Finished writing meshes.'
523 
524  ! read the files back in
525 call txi_read_x_format_file(mesh_a_file_f, mesh_a_out_f, status)
526 if (status < 0) then
527  WRITE(*,*) 'Error reading TestMeshA'
528  error = status
529  return
530 endif
531 
532  ! read the files back in
533 call txi_read_x_format_file(mesh_b_file_f, mesh_b_out_f, status)
534 if (status < 0) then
535  WRITE(*,*) 'Error reading TestMeshB'
536  error = status
537  return
538 endif
539 
540 WRITE(*,*) 'Finished reading meshes.'
541 
542 error = status
543 return
544 
545 END SUBROUTINE
546 
547 !**************************
548 
549 !---------------------------------------------------------------------------
550 ! FUNCTION txiTestC
551 ! PURPOSE test to see if fortran code can read file written with C.
552 ! NOTES
553 !---------------------------------------------------------------------------
554 SUBROUTINE txi_test_c (error)
555 INTEGER, INTENT(OUT) :: error
556 INTEGER nstatus
557 INTEGER xfileid
558 
559 error = 1
560 
561  !Check to see if files written with C exist
562  ! Turn off error handling
563 !call H5Eset_auto_f(0, error)
564  ! Try opening a file written with C to see if one exists.
565 call xf_open_file(scalar_a_file_c, .true., xfileid, nstatus)
566  ! If the file written with C doesn't exist, return.
567 if (nstatus < 0) then
568  call xf_close_file(xfileid, error)
569  ! Restore previous error handler
570  !call H5Eset_Auto_f(1, error)
571  error = -1
572  return
573  ! If the file written with C does exist, assume all C files exist.
574 else
575  call xf_close_file(xfileid, error)
576  ! Restore previous error handler
577  !call H5Eset_Auto_f(1, error)
578 endif
579 
580  ! Read the files back in
581 call txi_read_x_format_file(scalar_a_file_c, scalar_a_text_cf, error)
582 if (error < 0) then
583  return
584 endif
585 call txi_read_x_format_file(scalar_b_file_c, scalar_b_text_cf, error)
586 if (error < 0) then
587  return
588 endif
589 
590 call txi_read_x_format_file(vector2d_a_file_c, vector2d_a_text_cf, error)
591 if (error < 0) then
592  return
593 endif
594 call txi_read_x_format_file(vector2d_b_file_c, vector2d_b_text_cf, error)
595 if (error < 0) then
596  return
597 endif
598 
599 WRITE(*,*) 'Done reading C datasets...'
600 
601 call txi_read_x_format_file(grid_cart2d_a_file_c, grid_cart2d_a_out_cf, error)
602 if (error < 0) then
603  WRITE(*,*) 'Error reading C grid Cartesian 2D A'
604 endif
605 WRITE(*,*) 'Finished reading C grid Cartesian 2D A'
606 
607 call txi_read_x_format_file(grid_curv2d_a_file_c, grid_curv2d_a_out_cf, error)
608 if (error < 0) then
609  WRITE(*,*) 'Error reading C grid Curvilinear 2D A'
610 endif
611 WRITE(*,*) 'Finished reading C grid Curvilinear 2D A'
612 
613 call txi_read_x_format_file(grid_cart3d_a_file_c, grid_cart3d_a_out_cf, error)
614 if (error < 0) then
615  WRITE(*,*) 'Error reading C grid Cartesian 3D A'
616 endif
617 WRITE(*,*) 'Finished reading C grid Cartesian 3D A'
618 
619  ! read the files back in
620 call txi_read_x_format_file(mesh_a_file_c, mesh_a_out_cf, error)
621 if (error < 0) then
622  WRITE(*,*) 'Error reading C TestMeshA'
623  return
624 endif
625 
626  ! read the files back in
627 call txi_read_x_format_file(mesh_b_file_c, mesh_b_out_cf, error)
628 if (error < 0) then
629  WRITE(*,*) 'Error reading C TestMeshB'
630  return
631 endif
632 
633 WRITE(*,*) 'Finished reading C meshes.'
634 
635 return
636 
637 END SUBROUTINE
638 
639 !**************************
640 ! --------------------------------------------------------------------------
641 ! FUNCTION txiReadXFormatFile
642 ! PURPOSE Read a file using XMDF and write information about the data
643 ! contained in the file to a output file
644 ! --------------------------------------------------------------------------
645 SUBROUTINE txi_read_x_format_file (a_XmdfFile, a_OutFile, error)
646 CHARACTER(LEN=*), INTENT(IN) :: a_xmdffile
647 CHARACTER(LEN=*), INTENT(IN) :: a_outfile
648 INTEGER, INTENT(OUT) :: error
649 CHARACTER(LEN=BIG_STRING_SIZE) :: individualpath
650 CHARACTER,ALLOCATABLE :: paths(:)
651 INTEGER xfileid, xgroupid
652 INTEGER nmeshgroups, nmaxpathlength, ngridgroups
653 INTEGER fileunit, startloc, nstatus, i, j
654 REAL version
655 
656 xfileid = none
657 xgroupid = none
658 
659  ! Open the XMDF file
660 call xf_open_file(a_xmdffile, .true., xfileid, nstatus)
661 if (nstatus < 0) then
662  call xf_close_file(xfileid, error)
663  error = -1
664  return
665 endif
666 
667  ! open the status file
668 fileunit = 53
669 OPEN(unit=fileunit, file=a_outfile, status='REPLACE', action='WRITE', &
670  iostat = error)
671 if (fileunit == 0) then
672  call xf_close_file(xfileid, error)
673  error = -1
674  return
675 endif
676 
677 WRITE(fileunit,*) 'File ', a_xmdffile, ' opened.'
678 
679  ! write the version number to the file
680 call xf_get_library_version_file(xfileid, version, error)
681 WRITE(fileunit,*) 'XMDF Version: ', version
682 WRITE(fileunit,*) ''
683 
684  ! Read Coordinate System Informatioin to the .txt file if contained in
685  ! file, if not skip
686 call txi_test_coord_system(xfileid, fileunit, nstatus)
687 WRITE(fileunit,*) ''
688 
689  ! read all datasets not beneath a mesh, grid, or cross-sections
690 call td_read_datasets(xfileid,fileunit, nstatus)
691 if (nstatus < 0) then
692  call xf_close_file(xfileid, error)
693  error = -1
694  return
695 endif
696 
697  ! Get the number and paths of datasets in the file.
698 call xf_grp_pths_sz_for_meshes(xfileid, nmeshgroups, &
699  nmaxpathlength, error)
700 if (error >= 0 .AND. nmeshgroups > 0) then
701  allocate (paths(nmaxpathlength*nmeshgroups))
702  call xf_get_group_paths_for_meshes(xfileid, nmeshgroups, nmaxpathlength, &
703  paths, error)
704 endif
705 
706 if (error < 0) then
707  call xf_close_file(xfileid, error)
708  error = -1
709  return
710 endif
711 
712  ! Report the number and paths to individual meshes in the file.
713 WRITE(fileunit,*) 'Number of meshes in file: ', nmeshgroups
714 WRITE(fileunit,*) 'Paths:'
715 do i=1, nmeshgroups
716  do j=1, nmaxpathlength-1
717  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
718  enddo
719  WRITE(fileunit,*) individualpath(1:nmaxpathlength-1)
720 enddo
721 
722 WRITE(fileunit,*) ''
723 
724  ! Open each mesh group
725 !if (nMeshGroups > 0) allocate(IndividualPath(nMaxPathLength + 1))
726 
727 do i=1, nmeshgroups
728  ! copy the portion of the array where a single path is stored
729  startloc = (i-1)*nmaxpathlength + 1
730  individualpath = ''
731  do j = 1, nmaxpathlength - 1
732  individualpath(j:j) = paths(startloc+j-1)
733  enddo
734 
735  WRITE(fileunit,*) 'Reading mesh in group: ', &
736  individualpath(1:nmaxpathlength-1)
737  call xf_open_group(xfileid, individualpath(1:len_trim(individualpath)), &
738  xgroupid, nstatus)
739  if (nstatus >= 0) then
740  call tm_read_mesh(xgroupid, fileunit, nstatus)
741  endif
742  if (nstatus < 0) then
743  WRITE(*,*) 'Error reading mesh..'
744  endif
745 enddo
746 
747 if (allocated(paths)) deallocate(paths)
748 !if (allocated(IndividualPath)) deallocate(IndividualPath)
749 
750  ! Grid stuff
751 call xf_grp_pths_sz_for_grids(xfileid, ngridgroups, &
752  nmaxpathlength, nstatus)
753 if (nstatus >= 0 .AND. ngridgroups > 0) then
754  allocate (paths(nmaxpathlength*ngridgroups))
755  call xf_get_group_paths_for_grids(xfileid, ngridgroups, &
756  nmaxpathlength, paths, nstatus)
757 endif
758 if (nstatus < 0) then
759  call xf_close_file(xfileid, error)
760  error = -1
761  return
762 endif
763 
764  ! Report the number and paths to individual meshes in the file.
765 WRITE(fileunit,*) 'Number of grids in file: ', ngridgroups
766 WRITE(fileunit,*) 'Paths:'
767 do i=1, ngridgroups
768  do j=1, nmaxpathlength-1
769  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
770  enddo
771  WRITE(fileunit,*) individualpath(1:len_trim(individualpath))
772 enddo
773 
774 WRITE(fileunit,*) ''
775 
776 !if (nGridGroups > 0) allocate(IndividualPath(nMaxPathLength + 1))
777 
778  ! Open each grid group
779 do i=1, ngridgroups
780  do j = 1, nmaxpathlength - 1
781  individualpath(j:j) = paths((i-1)*nmaxpathlength+j)
782  enddo
783  WRITE(fileunit,*) 'Reading grid in group: ', &
784  individualpath(1:len_trim(individualpath))
785  call xf_open_group(xfileid, individualpath(1:len_trim(individualpath)), &
786  xgroupid, nstatus)
787  if (nstatus >= 0) then
788  call tg_read_grid(xgroupid, fileunit, nstatus)
789  endif
790  if (nstatus < 0) then
791  WRITE(fileunit,*) 'Error reading grid..'
792  endif
793 enddo
794 
795 if (allocated(paths)) deallocate(paths)
796 !if (allocated(IndividualPath)) deallocate(IndividualPath)
797 
798  ! TODO do grid, and cross-section stuff.
799 
800  ! close the files
801 call xf_close_file(xfileid, error)
802 CLOSE(fileunit)
803 
804 return
805 
806 END SUBROUTINE
807 
808 !-----------------------------------------------------------------------------
809 ! SUBROUTINE TXI_TestCalendar
810 ! PURPOSE Check the Calculations of Julian date from calendar date or
811 ! vice-versa.
812 ! NOTES era is #defined (use #defines): ERA_IS_BCE (BC), ERA_IS_CE (AD)
813 !-----------------------------------------------------------------------------
814 SUBROUTINE txi_test_calendar (error)
815  INTEGER, INTENT(OUT) :: error
816  INTEGER era1, yr1, mo1, day1, hr1, min1, sec1
817  INTEGER era2, yr2, mo2, day2, hr2, min2, sec2
818  INTEGER era3, yr3, mo3, day3, hr3, min3, sec3, fileunit
819  INTEGER era4, yr4, mo4, day4, hr4, min4, sec4, calendarworks
820  DOUBLE PRECISION julian1, julian2, julian3, julian4
821 
822  calendarworks = 0
823  fileunit = 53
824  OPEN(unit=fileunit, file=calendar_out_f, status='REPLACE', action='WRITE', &
825  iostat = error)
826 
827  WRITE(fileunit,*) 'Calendar conversion:'
828 
829  era1 = era_is_bce
830  yr1 = 0
831  mo1 = 0
832  day1 = 0
833  hr1 = 0
834  min1 = 0
835  sec1 = 0
836  julian1 = 2655.5
837  call xf_julian_to_calendar(era1, yr1, mo1, day1, hr1, min1, sec1, julian1, error)
838 
839  era2 = era_is_bce
840  yr2 = 4706
841  mo2 = 4
842  day2 = 10
843  hr2 = 0
844  min2 = 0
845  sec2 = 0
846  julian2 = 0.0
847  call xf_calendar_to_julian(era2, yr2, mo2, day2, hr2, min2, sec2, julian2, error)
848 
849  era3 = era_is_ce
850  yr3 = 2004
851  mo3 = 6
852  day3 = 3
853  hr3 = 2
854  min3 = 8
855  sec3 = 32
856  julian3 = 0.0
857  call xf_calendar_to_julian(era3, yr3, mo3, day3, hr3, min3, sec3, julian3, error)
858 
859  era4 = era_is_bce
860  yr4 = 0
861  mo4 = 0
862  day4 = 0
863  hr4 = 0
864  min4 = 0
865  sec4 = 0
866  julian4 = 2453159.5892592594_double
867  call xf_julian_to_calendar(era4, yr4, mo4, day4, hr4, min4, sec4, julian4, error)
868 
869 WRITE(fileunit,*) ''
870 WRITE(fileunit,*) 'Dates #1 & #2 were calculated with the same date:'
871 WRITE(fileunit,*) ''
872 WRITE(fileunit,*) era1, '/', yr1, '/', mo1, '/', day1
873 WRITE(fileunit,*) '', hr1, ':', min1, ':',sec1, '--- julian =',julian1
874 WRITE(fileunit,*) ''
875 WRITE(fileunit,*) era2, '/', yr2, '/', mo2, '/', day2
876 WRITE(fileunit,*) '', hr2, ':', min2, ':',sec2, '--- julian =',julian2
877 WRITE(fileunit,*) ''
878 WRITE(fileunit,*) 'Dates #3 & #4 were calculated with the same date:'
879 WRITE(fileunit,*) ''
880 WRITE(fileunit,*) era3, '/', yr3, '/', mo3, '/', day3
881 WRITE(fileunit,*) '', hr3, ':', min3, ':',sec3, '--- julian =',julian3
882 WRITE(fileunit,*) ''
883 WRITE(fileunit,*) era4, '/', yr4, '/', mo4, '/', day4
884 WRITE(fileunit,*) '', hr4, ':', min4, ':',sec4, '--- julian =',julian4
885 
886  if (era1==era2 .AND. era3==era4) then
887  if (yr1==yr2 .AND. yr3==yr4) then
888  if (mo1==mo2 .AND. mo3==mo4) then
889  if (day1==day2 .AND. day3==day4) then
890  if (hr1==hr2 .AND. hr3==hr4) then
891  if (min1==min2 .AND. min3==min4) then
892  if (julian1==julian2 .AND. (julian3-.0000001)<julian4 .AND. (julian3+.0000001)>julian4) then
893  WRITE(*,*) ''
894  WRITE(*,*) 'Calendar conversion works correctly.'
895  calendarworks = 1
896  endif
897  endif
898  endif
899  endif
900  endif
901  endif
902  endif
903  if (calendarworks .NE. 1) then
904  WRITE(*,*) 'Calendar Stuff DOES NOT Work Correctly'
905  error = -1
906  else
907  error = 1
908  endif
909  return
910 
911 END SUBROUTINE
912 !------------------------------------------------------------------------------
913 ! FUNCTION TXI_TEST_MULTI_DATASETS
914 ! PURPOSE To test the newly translated functions
915 ! NOTES
916 !------------------------------------------------------------------------------
917 SUBROUTINE txi_test_multi_datasets
918 
919 !call XF_SETUP_TO_WRITE_DATASETS(a_Filename, a_MultiDatasetsGroupPath, &
920 ! a_PathInMultiDatasetsGroup, a_SpatialDataObjectGuid, &
921  ! a_OverwriteOptions, a_FileId, a_GroupId, error)
922 !XF_OPEN_MULTI_DATASETS_GROUP
923 !XF_DELETE_GROUP
924 
925 END SUBROUTINE
926 
927 !------------------------------------------------------------------------------
928 ! SUBROUTINE TXI_WRITE_XMDF_VERSION
929 ! PURPOSE Write the XMDF version number to the screen
930 ! NOTES
931 !------------------------------------------------------------------------------
932 SUBROUTINE txi_test_version ()
933 INTEGER error
934 REAL version
935 
936  WRITE(*,*) ''
937  call xf_get_library_version(version, error)
938  WRITE(*,*) 'The current version of XMDF is: ', version
939  WRITE(*,*) ''
940 
941  return
942 
943 END SUBROUTINE
944 
945 !------------------------------------------------------------------------------
946 ! SUBROUTINE TXI_TEST_COORD_SYSTEM
947 ! PURPOSE Reads a file's Coordinate Group and prints coordinate data out
948 ! to each text file.
949 ! NOTES
950 !------------------------------------------------------------------------------
951 SUBROUTINE txi_test_coord_system (xFileId, a_OutFile, error)
952 INTEGER, INTENT(IN) :: xfileid
953 INTEGER, INTENT(IN) :: a_outfile
954 INTEGER, INTENT(OUT) :: error
955 INTEGER ihorizdatum, ihorizunits, ivertdatum, ivertunits
956 INTEGER ilat, ilon, iutmzone, ispczone, ihpgnarea, iellipse
957 INTEGER bhorizdatum, nstatus
958 REAL(DOUBLE) dcpplat, dcpplon, dmajorr, dminorr
959 INTEGER xcoordid
960 CHARACTER(LEN=BIG_STRING_SIZE) strhorizunits, strvertdatum
961 CHARACTER(LEN=BIG_STRING_SIZE) strvertunits
962 
963  ! Coordinate stuff
964  ! Read Coordinate Info
965  ! Open the Coordinate group
966  !call H5Eset_auto_f(0, error)
967  call xf_open_coordinate_group(xfileid, xcoordid, nstatus)
968  if (nstatus < 0) then
969  WRITE(a_outfile,*) ''
970  WRITE(a_outfile,*) 'Coordinate Group not found'
971  WRITE(a_outfile,*) ''
972  error = nstatus
973  return
974  endif
975  !call H5Eset_auto_f(1, error)
976 
977  WRITE(a_outfile,*) ''
978  WRITE(a_outfile,*) 'Coordinate System:'
979 
980  call xf_get_horiz_datum(xcoordid, ihorizdatum, bhorizdatum)
981  call xf_get_horiz_units(xcoordid, ihorizunits, error)
982  call xf_get_vert_datum(xcoordid, ivertdatum, error)
983  call xf_get_vert_units(xcoordid, ivertunits, error)
984 
985  ! set horizontal units
986  if (ihorizunits == 0) then
987  strhorizunits = 'Horizontal units = US Survey Feet (=0)'
988  else if (ihorizunits == 1) then
989  strhorizunits = 'Horizontal units = International Feet (=1)'
990  else if (ihorizunits == 2) then
991  strhorizunits = 'Horizontal units = Meters (=2)'
992  else
993  strhorizunits = 'ERROR in reading Horizontal units'
994  endif
995 
996  ! set vertical datum
997  if (ivertdatum == 0) then
998  strvertdatum = 'Vertical datum = Local (=0)'
999  else if (ivertdatum == 1) then
1000  strvertdatum = 'Vertical datum = NGVD 29 (=1)'
1001  else if (ivertdatum == 2) then
1002  strvertdatum = 'Vertical datum = NGVD 88 (=2)'
1003  else
1004  strvertdatum = 'ERROR in reading the Vertical datum'
1005  endif
1006 
1007  ! set vertocal units
1008  if (ivertunits == 0) then
1009  strvertunits = 'Vertical units = US Survey Feet (=0)'
1010  else if (ivertunits == 1) then
1011  strvertunits = 'Vertical units = International Feet (=1)'
1012  else if (ivertunits == 2) then
1013  strvertunits = 'Vertical units = Meters (=2)'
1014  else
1015  strvertunits = 'ERROR in reading the Vertical units'
1016  endif
1017 
1018  if (bhorizdatum >= 0) then
1019  SELECT CASE (ihorizdatum)
1020  CASE (horiz_datum_geographic)
1021  call xf_get_ellipse(xcoordid, iellipse, error)
1022  call xf_get_lat(xcoordid, ilat, error)
1023  call xf_get_lon(xcoordid, ilon, error)
1024  ! Write Horizontal and Vertical Info
1025  WRITE(a_outfile,*) 'Horizontal datum = Geographic'
1026  WRITE(a_outfile,*) 'Horizontal units = ', strhorizunits(1:len_trim(strhorizunits))
1027  WRITE(a_outfile,*) 'Vertical datum = ', strvertdatum(1:len_trim(strvertdatum))
1028  WRITE(a_outfile,*) 'Vertical units = ', strvertunits(1:len_trim(strvertunits))
1029  ! Write Latitude data
1030  if (ilat == 0) then
1031  WRITE(a_outfile,*) ' Latitude = North (=0)'
1032  else if (ilat == 1) then
1033  WRITE(a_outfile,*) ' Latitude = South (=1)'
1034  else
1035  WRITE(a_outfile,*) ' LATITUDE INFO INCORRECT'
1036  endif
1037  ! Write Longitude data
1038  if (ilon == 0) then
1039  WRITE(a_outfile,*) ' Longitude = East (=0)'
1040  else if (ilon == 1) then
1041  WRITE(a_outfile,*) ' Longitude = West (=1)'
1042  else
1043  WRITE(a_outfile,*) ' LONGITUDE INFO INCORRECT'
1044  endif
1045  ! Ellipse Information
1046  ! User-defined Ellipse (==32)
1047  if (iellipse == 32) then
1048  WRITE(a_outfile,*) 'Ellipse = User-defined:'
1049  call xf_get_major_r(xcoordid, dmajorr, error)
1050  call xf_get_minor_r(xcoordid, dminorr, error)
1051  WRITE(a_outfile,*) ' MajorR = ', dmajorr
1052  WRITE(a_outfile,*) ' MinorR = ', dminorr
1053  else
1054  WRITE(a_outfile,*) 'Ellipse = ', iellipse
1055  endif
1056  WRITE(a_outfile,*) ''
1057  CASE (horiz_datum_utm)
1058  call xf_get_utm_zone(xcoordid, iutmzone, error)
1059  ! output info to text file
1060  if (ihorizdatum == horiz_datum_utm) then
1061  WRITE(a_outfile,*) 'Horizontal datum = UTM'
1062  else if (ihorizdatum == horiz_datum_utm_nad27) then
1063  WRITE(a_outfile,*) 'Horizontal datum = UTM NAD27 (US)'
1064  else
1065  WRITE(a_outfile,*) 'Horizontal datum = UTM NAD83 (US)'
1066  endif
1067  WRITE(a_outfile,*) 'Horizontal units = ', &
1068  strhorizunits(1:len_trim(strhorizunits))
1069  WRITE(a_outfile,*) 'Vertical datum = ', &
1070  strvertdatum(1:len_trim(strvertdatum))
1071  WRITE(a_outfile,*) 'Vertical units = ', &
1072  strvertunits(1:len_trim(strvertunits))
1073  WRITE(a_outfile,*) 'UTM Zone = ', iutmzone
1074  WRITE(a_outfile,*) ''
1075 
1076  CASE (horiz_datum_utm_nad27, horiz_datum_utm_nad83)
1077  call xf_get_utm_zone(xcoordid, iutmzone, error)
1078  ! output info to text file
1079  if (ihorizdatum == horiz_datum_utm) then
1080  WRITE(a_outfile,*) 'Horizontal datum = UTM'
1081  else if (ihorizdatum == horiz_datum_utm_nad27) then
1082  WRITE(a_outfile,*) 'Horizontal datum = UTM NAD27 (US)'
1083  else
1084  WRITE(a_outfile,*) 'Horizontal datum = UTM NAD83 (US)'
1085  endif
1086  WRITE(a_outfile,*) 'Horizontal units = ', &
1087  strhorizunits(1:len_trim(strhorizunits))
1088  WRITE(a_outfile,*) 'Vertical datum = ', &
1089  strvertdatum(1:len_trim(strvertdatum))
1090  WRITE(a_outfile,*) 'Vertical units = ', &
1091  strvertunits(1:len_trim(strvertunits))
1092  WRITE(a_outfile,*) 'UTM Zone = ', iutmzone
1093  WRITE(a_outfile,*) ''
1094  CASE (horiz_datum_state_plane_nad27, horiz_datum_state_plane_nad83)
1095  call xf_get_spc_zone(xcoordid, ispczone, error)
1096  ! output info to text file
1097  if (ihorizdatum == horiz_datum_state_plane_nad27) then
1098  WRITE(a_outfile,*) 'Horizontal datum = State Plane NAD27 (US)'
1099  else
1100  WRITE(a_outfile,*) 'Horizontal datum = State Plane NAD83 (US)'
1101  endif
1102  WRITE(a_outfile,*) 'Horizontal units = ', &
1103  strhorizunits(1:len_trim(strhorizunits))
1104  WRITE(a_outfile,*) 'Vertical datum = ', &
1105  strvertdatum(1:len_trim(strvertdatum))
1106  WRITE(a_outfile,*) 'Vertical units = ', &
1107  strvertunits(1:len_trim(strvertunits))
1108  WRITE(a_outfile,*) 'SPC Zone = ', ispczone
1109  WRITE(a_outfile,*) ''
1110  CASE (horiz_datum_utm_hpgn, horiz_datum_state_plane_hpgn, &
1111  horiz_datum_geographic_hpgn)
1112  call xf_get_hpgn_area(xcoordid, ihpgnarea, error)
1113  if (ihorizdatum == horiz_datum_utm_hpgn) then
1114  WRITE(a_outfile,*) 'Horizontal datum = UTM HPGN (US)'
1115  else if (ihorizdatum == horiz_datum_state_plane_hpgn) then
1116  WRITE(a_outfile,*) 'Horizontal datum = State Plane HPGN (US)'
1117  else
1118  WRITE(a_outfile,*) 'Horizontal datum = Geographic HPGN (US)'
1119  endif
1120  WRITE(a_outfile,*) 'Horizontal units = ', &
1121  strhorizunits(1:len_trim(strhorizunits))
1122  WRITE(a_outfile,*) 'Vertical datum = ', &
1123  strvertdatum(1:len_trim(strvertdatum))
1124  WRITE(a_outfile,*) 'Vertical units = ', &
1125  strvertunits(1:len_trim(strvertunits))
1126  WRITE(a_outfile,*) 'HPGN Area = ', ihpgnarea
1127  WRITE(a_outfile,*) ''
1128  CASE (horiz_datum_cpp)
1129  call xf_get_cpp_lat(xcoordid, dcpplat, error)
1130  call xf_get_cpp_lon(xcoordid, dcpplon, error)
1131  WRITE(a_outfile,*) 'Horizontal datum = CPP (Carte Parallelo-Grammatique Projection)'
1132  WRITE(a_outfile,*) 'Horizontal units = ', &
1133  strhorizunits(1:len_trim(strhorizunits))
1134  WRITE(a_outfile,*) 'Vertical datum = ', &
1135  strvertdatum(1:len_trim(strvertdatum))
1136  WRITE(a_outfile,*) 'Vertical units = ', &
1137  strvertunits(1:len_trim(strvertunits))
1138  WRITE(a_outfile,*) 'CPP Latitude = ', dcpplat
1139  WRITE(a_outfile,*) 'CPP Longitude = ', dcpplon
1140  WRITE(a_outfile,*) ''
1141  CASE (horiz_datum_local, horiz_datum_geographic_nad27, &
1142  horiz_datum_geographic_nad83)
1143  ! do other systems
1144  if (ihorizdatum == horiz_datum_local) then
1145  WRITE(a_outfile,*) 'Horizontal datum = Local'
1146  else if (ihorizdatum == horiz_datum_geographic_nad27) then
1147  WRITE(a_outfile,*) 'Horizontal datum = Geographic NAD27 (US)'
1148  else
1149  WRITE(a_outfile,*) 'Horizontal datum = Geographic NAD83 (US)'
1150  endif
1151  WRITE(a_outfile,*) 'Horizontal units = ', &
1152  strhorizunits(1:len_trim(strhorizunits))
1153  WRITE(a_outfile,*) 'Vertical datum = ', &
1154  strvertdatum(1:len_trim(strvertdatum))
1155  WRITE(a_outfile,*) 'Vertical units = ', &
1156  strvertunits(1:len_trim(strvertunits))
1157  WRITE(a_outfile,*) ''
1158  CASE default
1159  WRITE(a_outfile,*) 'ERROR: The coordinate information is not found in the .h5 file'
1160  error = -1
1161  return
1162  END SELECT
1163  else
1164  WRITE(a_outfile,*) 'Coordinate information in HDF5 file is incomplete.'
1165  WRITE(a_outfile,*) ''
1166  endif
1167 
1168  call xf_close_group(xcoordid, error)
1169  xcoordid = 0
1170 
1171  return
1172 
1173 END SUBROUTINE
1174 
1175 SUBROUTINE txi_test_geometric_paths(error)
1176  INTEGER, INTENT(OUT) :: error
1177  INTEGER compression
1178 
1179  compression = -1
1180 
1181  ! test writing a geometric path file */
1182  WRITE(*,*)''
1183  WRITE(*,*)'Writing geometric path data'
1184  WRITE(*,*)''
1185 
1186  call tm_write_test_paths(geompath_a_file_f, compression, error);
1187  if (error < 0) then
1188  WRITE(*,*) 'Error writing geometric path data A'
1189  return
1190  endif
1191  WRITE(*,*) 'Finished writing geometric path data A'
1192 
1193  ! test reading a geometric path file */
1194  call tm_read_test_paths(geompath_a_file_f, geompath_a_file_f_out, error)
1195 
1196  return
1197 
1198 END SUBROUTINE txi_test_geometric_paths
1199 
1200 !****************************
1201 
1202 END MODULE testsmodule
1203 
1204 PROGRAM tests
1205 
1206 USE testdatasets
1207 USE xmdf
1208 USE testmesh
1209 USE testdefs
1210 USE testsmodule
1211 
1212 INTEGER error
1213 ! Assign pp to force crash on error
1214 
1215 call xf_initialize(error)
1216 
1217  ! test the dataset routines
1218 call txi_test_timesteps(error)
1219 if (error < 0) then
1220  WRITE(*,*) 'Error in writing timesteps!'
1221  pause 'Press ENTER to exit...'
1222 endif
1223 
1224  ! test the dataset routines
1225 call txi_test_datasets(error)
1226 if (error < 0) then
1227  WRITE(*,*) 'Error in writing datasets!'
1228  pause 'Press ENTER to exit...'
1229 endif
1230 
1231  ! test overwriting datasets
1232 call txi_test_overwrite_dsets(error)
1233 if (error < 0) then
1234  WRITE(*,*) 'Error in overwriting datasets!'
1235  pause 'Press ENTER to exit...'
1236 else
1237  WRITE(*,*) 'Finished writing datasets...'
1238 endif
1239 
1240  ! test mesh stuff
1241 call txi_test_meshs(error)
1242 if (error < 0) then
1243  WRITE(*,*) 'Error in TXI_TEST_MESHS!'
1244  pause 'Press ENTER to exit...'
1245 endif
1246 
1247  ! test grid stuff
1248 call txi_test_grids(error)
1249 if (error < 0) then
1250  WRITE(*,*) 'Error in TXI_TEST_GRIDS!'
1251  pause 'Press ENTER to exit...'
1252 endif
1253 
1254  ! test c in fortran
1255 call txi_test_c(error)
1256 ! Ignore error because fortran will be run once again with file there
1257 
1258  ! test calendar stuff
1259 call txi_test_calendar(error)
1260 if (error < 0) then
1261  WRITE(*,*) 'Error in TXI_TEST_CALENDAR!'
1262  pause 'Press ENTER to exit...'
1263 endif
1264 
1265  ! test version
1266 call txi_test_version
1267 
1268  ! test geometric paths
1269 call txi_test_geometric_paths(error)
1270 if (error < 0) then
1271  WRITE(*,*) 'Error in TXI_TEST_GEOMETRIC_PATHS!'
1272  pause 'Press ENTER to exit...'
1273 endif
1274 
1275 call xf_close(error)
1276 if (error < 0) then
1277  WRITE(*,*) 'Error in XF_CLOSE!'
1278  pause 'Press ENTER to exit...'
1279 endif
1280 
1281 ! PAUSE 'Press ENTER to exit...'
1282 
1283 END PROGRAM
Generated on Fri Jan 20 2017 10:10:01 for XMDF by   doxygen 1.8.6