[Dart-dev] [5673] DART/branches/development/location: added a [0, 1] periodic 3d locations mod, and an X, Y, Z 3d Cartesian one.
nancy at ucar.edu
nancy at ucar.edu
Mon Apr 9 11:11:37 MDT 2012
Revision: 5673
Author: nancy
Date: 2012-04-09 11:11:37 -0600 (Mon, 09 Apr 2012)
Log Message:
-----------
added a [0,1] periodic 3d locations mod, and an X,Y,Z 3d Cartesian one.
Added Paths:
-----------
DART/branches/development/location/location_test3.f90
DART/branches/development/location/threed/
DART/branches/development/location/threed/location_mod.f90
DART/branches/development/location/threed/test/
DART/branches/development/location/threed/test/input.nml
DART/branches/development/location/threed/test/mkmf_location_test
DART/branches/development/location/threed/test/path_names_location_test
DART/branches/development/location/threed/test/test.in
DART/branches/development/location/threed_cartesian/
DART/branches/development/location/threed_cartesian/location_mod.f90
DART/branches/development/location/threed_cartesian/location_mod.html
DART/branches/development/location/threed_cartesian/location_mod.nml
DART/branches/development/location/threed_cartesian/test/
DART/branches/development/location/threed_cartesian/test/input.nml
DART/branches/development/location/threed_cartesian/test/mkmf_location_test
DART/branches/development/location/threed_cartesian/test/mkmf_location_test3
DART/branches/development/location/threed_cartesian/test/path_names_location_test
DART/branches/development/location/threed_cartesian/test/path_names_location_test3
DART/branches/development/location/threed_cartesian/test/test.in
Removed Paths:
-------------
DART/branches/development/location/threed/location_mod.f90
DART/branches/development/location/threed/test/
DART/branches/development/location/threed/test/input.nml
DART/branches/development/location/threed/test/mkmf_location_test
DART/branches/development/location/threed/test/path_names_location_test
DART/branches/development/location/threed/test/test.in
DART/branches/development/location/threed_cartesian/location_mod.f90
DART/branches/development/location/threed_cartesian/location_mod.html
DART/branches/development/location/threed_cartesian/location_mod.nml
DART/branches/development/location/threed_cartesian/test/
DART/branches/development/location/threed_cartesian/test/input.nml
DART/branches/development/location/threed_cartesian/test/mkmf_location_test
DART/branches/development/location/threed_cartesian/test/mkmf_location_test3
DART/branches/development/location/threed_cartesian/test/path_names_location_test
DART/branches/development/location/threed_cartesian/test/path_names_location_test3
DART/branches/development/location/threed_cartesian/test/test.in
-------------- next part --------------
Copied: DART/branches/development/location/location_test3.f90 (from rev 5672, DART/branches/mpas/location/location_test3.f90)
===================================================================
--- DART/branches/development/location/location_test3.f90 (rev 0)
+++ DART/branches/development/location/location_test3.f90 2012-04-09 17:11:37 UTC (rev 5673)
@@ -0,0 +1,161 @@
+! DART software - Copyright 2004 - 2011 UCAR. This open source software is
+! provided by UCAR, "as is", without charge, subject to all terms of use at
+! http://www.image.ucar.edu/DAReS/DART/DART_download
+
+program location_test3
+
+! <next few lines under version control, do not edit>
+! $URL$
+! $Id$
+! $Revision$
+! $Date$
+
+! Simple test program to exercise oned location module.
+
+use location_mod
+use types_mod, only : r8
+use utilities_mod, only : get_unit, error_handler, E_ERR
+use random_seq_mod, only : random_seq_type, init_random_seq, random_uniform
+
+
+implicit none
+
+! version controlled file description for error handling, do not edit
+character(len=128), parameter :: &
+ source = "$URL$", &
+ revision = "$Revision$", &
+ revdate = "$Date$"
+
+integer, parameter :: nl = 20
+type(location_type) :: loc1(nl), loc2, near1
+integer :: iunit, i, j, dummy(nl), num_close, close_ind(nl), rc, near1index
+real(r8) :: loc2_val(3), dist(nl)
+real(r8) :: minv(3), maxv(3), v(3), minr(3), maxr(3), maxdist
+type(random_seq_type) :: ran_seq
+logical :: ran_seq_init = .false.
+type(get_close_type) :: cc_gc
+character(len=100) :: buf
+
+
+
+! Open an output file
+iunit = get_unit()
+open(iunit, file = 'location_test_file')
+
+if(.not. ran_seq_init) then
+ call init_random_seq(ran_seq)
+ ran_seq_init = .TRUE.
+endif
+
+minv(1) = -100.0_r8
+minv(2) = -200.0_r8
+minv(3) = -500.0_r8
+
+maxv(1) = 10000.0_r8
+maxv(2) = 20000.0_r8
+maxv(3) = 50000.0_r8
+
+minr(:) = 1e38_r8
+maxr(:) = -1e38_r8
+
+maxdist = 10000.0_r8
+
+! set a known min/max at the box edges so we know any
+! random points generated will be inside for now.
+loc1(1) = set_location(minv(1), minv(2), minv(3))
+call write_location(0, loc1(1), charstring=buf)
+write(*,*) 'location ', 1, trim(buf)
+
+loc1(2) = set_location(maxv(1), maxv(2), maxv(3))
+call write_location(0, loc1(2), charstring=buf)
+write(*,*) 'location ', 2, trim(buf)
+
+! Set the location list
+do i = 3, nl
+ do j=1, 3
+ v(j) = random_uniform(ran_seq) * (maxv(j)-minv(j)) + minv(j)
+ if (v(j) < minr(j)) minr(j) = v(j)
+ if (v(j) > maxr(j)) maxr(j) = v(j)
+ enddo
+
+ loc1(i) = set_location(v(1), v(2), v(3))
+
+ call write_location(0, loc1(i), charstring=buf)
+ write(*,*) 'location ', i, trim(buf)
+enddo
+
+! Write this location to the file
+do i = 1, nl
+ call write_location(iunit, loc1(i))
+enddo
+
+call get_close_maxdist_init(cc_gc, maxdist)
+call get_close_obs_init(cc_gc, nl, loc1)
+
+call print_get_close_type(cc_gc)
+
+dummy = 0
+
+! generate another random and find nearest from list
+do i = 1, nl
+ do j=1, 3
+ v(j) = random_uniform(ran_seq) * (maxr(j)-minr(j)) + minr(j)
+ enddo
+
+ loc2 = set_location(v(1), v(2), v(3))
+ call write_location(iunit, loc2)
+ call write_location(0, loc2, charstring=buf)
+ print *, 'generated a random point at ', trim(buf)
+
+ do j=1, nl
+ dist(j) = get_dist(loc1(j), loc2)
+ if (dist(j) <= maxdist) then
+ print *, 'dist to point ', j, 'is less than maxdist', dist(j)
+ endif
+ enddo
+
+
+ call get_close_obs(cc_gc, loc2, 0, loc1, dummy, num_close, close_ind, dist)
+ if (num_close > 0) then
+ print *, 'num close = ', num_close
+ do j=1, min(num_close, nl)
+ print *, j, close_ind(j)
+ if (close_ind(j) >= 1 .and. close_ind(j) <= nl) then
+ call write_location(0, loc1(close_ind(j)), charstring=buf)
+ write(*,*) 'close box loc ', trim(buf), dist(j)
+ endif
+ enddo
+ endif
+
+ call find_nearest(cc_gc, loc2, loc1, near1index, rc)
+ if (rc /= 0) then
+ print *, 'bad return from find nearest, ', rc
+ cycle
+ endif
+ if (near1index < 1) then
+ print *, 'near1index < 1, ', near1index
+ cycle
+ endif
+ print *, 'nearest location index = ', near1index
+ call write_location(0, loc1(near1index), charstring=buf)
+ write(*,*) 'near loc ', trim(buf)
+
+enddo
+
+close(iunit)
+
+! Now read them back in and compute the distances from loc1
+open(iunit, file = 'location_test_file')
+
+do i = 1, nl
+ loc1(i) = read_location(iunit)
+enddo
+do i = 1, nl
+ loc2 = read_location(iunit)
+ write(*, *) 'distance ', i, ' is ', get_dist(loc1(i), loc2)
+enddo
+
+close(iunit)
+
+end program location_test3
+
Deleted: DART/branches/development/location/threed/location_mod.f90
===================================================================
--- DART/branches/mpas/location/threed/location_mod.f90 2012-04-09 16:48:53 UTC (rev 5672)
+++ DART/branches/development/location/threed/location_mod.f90 2012-04-09 17:11:37 UTC (rev 5673)
@@ -1,729 +0,0 @@
-! DART software - Copyright 2004 - 2011 UCAR. This open source software is
-! provided by UCAR, "as is", without charge, subject to all terms of use at
-! http://www.image.ucar.edu/DAReS/DART/DART_download
-
-module location_mod
-
-! <next few lines under version control, do not edit>
-! $URL$
-! $Id$
-! $Revision$
-! $Date$
-
-! Implements location interfaces for a three dimensional cyclic region.
-! The internal representation of the location is currently implemented
-! as (x, y) from 0.0 to 1.0 in both dimensions.
-!
-! If you are looking for a geophysical locations module, look at either
-! the threed_sphere or threed_cartesian versions of this file.
-
-use types_mod, only : r8, MISSING_R8
-use utilities_mod, only : register_module, error_handler, E_ERR, ascii_file_format, &
- nc_check
-use random_seq_mod, only : random_seq_type, init_random_seq, random_uniform
-
-implicit none
-private
-
-public :: location_type, get_location, set_location, &
- set_location_missing, is_location_in_region, &
- write_location, read_location, interactive_location, query_location, &
- LocationDims, LocationName, LocationLName, get_close_obs, &
- get_close_maxdist_init, get_close_obs_init, get_close_type, &
- operator(==), operator(/=), get_dist, get_close_obs_destroy, &
- nc_write_location_atts, nc_get_location_varids, nc_write_location, &
- vert_is_height, vert_is_pressure, vert_is_undef, vert_is_level, &
- vert_is_surface, has_vertical_localization
-
-! version controlled file description for error handling, do not edit
-character(len=128), parameter :: &
- source = "$URL$", &
- revision = "$Revision$", &
- revdate = "$Date$"
-
-type location_type
- private
- real(r8) :: x, y, z
-end type location_type
-
-! Needed as stub but not used in this low-order model
-type get_close_type
- private
- integer :: num
- real(r8) :: maxdist
-end type get_close_type
-
-type(random_seq_type) :: ran_seq
-logical :: ran_seq_init = .false.
-logical, save :: module_initialized = .false.
-
-integer, parameter :: LocationDims = 3
-character(len = 129), parameter :: LocationName = "loc3D"
-character(len = 129), parameter :: LocationLName = "threed cyclic locations: x, y, z"
-
-character(len = 129) :: errstring
-
-interface operator(==); module procedure loc_eq; end interface
-interface operator(/=); module procedure loc_ne; end interface
-
-interface set_location
- module procedure set_location_single
- module procedure set_location_array
-end interface set_location
-
-contains
-
-!----------------------------------------------------------------------------
-
-subroutine initialize_module
-
-if (module_initialized) return
-
-call register_module(source, revision, revdate)
-module_initialized = .true.
-
-end subroutine initialize_module
-
-!----------------------------------------------------------------------------
-
-function get_dist(loc1, loc2, kind1, kind2)
-
-! Return the distance between 2 locations. Since this is a periodic
-! domain, the shortest distance may wrap around.
-
-type(location_type), intent(in) :: loc1, loc2
-integer, optional, intent(in) :: kind1, kind2
-real(r8) :: get_dist
-
-real(r8) :: x_dif, y_dif, z_dif
-
-if ( .not. module_initialized ) call initialize_module
-
-! Periodic domain, if distance is greater than half wraparound the other way.
-x_dif = abs(loc1%x - loc2%x)
-if (x_dif > 0.5_r8) x_dif = 1.0_r8 - x_dif
-y_dif = abs(loc1%y - loc2%y)
-if (y_dif > 0.5_r8) y_dif = 1.0_r8 - y_dif
-z_dif = abs(loc1%z - loc2%z)
-if (z_dif > 0.5_r8) z_dif = 1.0_r8 - z_dif
-
-get_dist = sqrt ( x_dif * x_dif + y_dif * y_dif + z_dif * z_dif)
-
-end function get_dist
-
-!---------------------------------------------------------------------------
-
-function loc_eq(loc1,loc2)
-
-! interface operator used to compare two locations.
-! Returns true only if all components are 'the same' to within machine
-! precision.
-
-type(location_type), intent(in) :: loc1, loc2
-logical :: loc_eq
-
-if ( .not. module_initialized ) call initialize_module
-
-loc_eq = .false.
-
-if ( abs(loc1%x - loc2%x ) > epsilon(loc1%x ) ) return
-if ( abs(loc1%y - loc2%y ) > epsilon(loc1%y ) ) return
-if ( abs(loc1%z - loc2%z ) > epsilon(loc1%z ) ) return
-
-loc_eq = .true.
-
-end function loc_eq
-
-!---------------------------------------------------------------------------
-
-function loc_ne(loc1,loc2)
-
-! interface operator used to compare two locations.
-! Returns true if locations are not identical to machine precision.
-
-type(location_type), intent(in) :: loc1, loc2
-logical :: loc_ne
-
-if ( .not. module_initialized ) call initialize_module
-
-loc_ne = (.not. loc_eq(loc1,loc2))
-
-end function loc_ne
-
-!---------------------------------------------------------------------------
-
-function get_location(loc)
-
-! Given a location type, return the x, y, z
-
-type(location_type), intent(in) :: loc
-real(r8), dimension(3) :: get_location
-
-if ( .not. module_initialized ) call initialize_module
-
-get_location(1) = loc%x
-get_location(2) = loc%y
-get_location(3) = loc%z
-
-end function get_location
-
-!----------------------------------------------------------------------------
-
-function set_location_single(x, y, z)
-
-! Given an x, y, z triplet, put the values into the location.
-
-real(r8), intent(in) :: x, y, z
-type (location_type) :: set_location_single
-
-if ( .not. module_initialized ) call initialize_module
-
-if(x < 0.0_r8 .or. x > 1.0_r8) then
- write(errstring,*)'x (',x,') is not within range [0,1]'
- call error_handler(E_ERR, 'set_location', errstring, source, revision, revdate)
-endif
-
-if(y < 0.0_r8 .or. y > 1.0_r8) then
- write(errstring,*)'y (',y,') is not within range [0,1]'
- call error_handler(E_ERR, 'set_location', errstring, source, revision, revdate)
-endif
-
-if(z < 0.0_r8 .or. z > 1.0_r8) then
- write(errstring,*)'z (',z,') is not within range [0,1]'
- call error_handler(E_ERR, 'set_location', errstring, source, revision, revdate)
-endif
-
-set_location_single%x = x
-set_location_single%y = y
-set_location_single%z = z
-
-end function set_location_single
-
-!----------------------------------------------------------------------------
-
-function set_location_array(list)
-
-! location semi-independent interface routine
-! given 3 float numbers, call the underlying set_location routine
-
-real(r8), intent(in) :: list(:)
-type (location_type) :: set_location_array
-
-if ( .not. module_initialized ) call initialize_module
-
-if (size(list) < 3) then
- write(errstring,*)'requires 3 input values'
- call error_handler(E_ERR, 'set_location', errstring, source, revision, revdate)
-endif
-
-set_location_array = set_location_single(list(1), list(2), list(3))
-
-end function set_location_array
-
-!----------------------------------------------------------------------------
-
-function set_location_missing()
-
-! fill in the contents to a known value.
-
-type (location_type) :: set_location_missing
-
-if ( .not. module_initialized ) call initialize_module
-
-set_location_missing%x = MISSING_R8
-set_location_missing%y = MISSING_R8
-set_location_missing%z = MISSING_R8
-
-end function set_location_missing
-
-!---------------------------------------------------------------------------
-
-function query_location(loc,attr)
-
-! Returns the value of the attribute
-!
-
-type(location_type), intent(in) :: loc
-character(len=*), optional, intent(in) :: attr
-real(r8) :: query_location
-
-if ( .not. module_initialized ) call initialize_module
-
-! see the long comment in this routine in the threed_sphere
-! module for warnings about compiler bugs before you change
-! this code.
-
-query_location = loc%x
-
-if (.not. present(attr)) return
-
-select case(attr)
- case ('x','X')
- query_location = loc%x
- case ('y','Y')
- query_location = loc%y
- case ('z','Z')
- query_location = loc%z
- case default
- call error_handler(E_ERR, 'query_location; threed', &
- 'Only x, y, or z are legal attributes to request from location', source, revision, revdate)
-end select
-
-end function query_location
-
-!----------------------------------------------------------------------------
-
-subroutine write_location(locfile, loc, fform, charstring)
-
-! Writes a 3D location to the file.
-! additional functionality: if optional argument charstring is specified,
-! it must be long enough to hold the string, and the location information is
-! written into it instead of to a file. fform must be ascii (which is the
-! default if not specified) to use this option.
-
-integer, intent(in) :: locfile
-type(location_type), intent(in) :: loc
-character(len = *), intent(in), optional :: fform
-character(len = *), intent(out), optional :: charstring
-
-integer :: charlength
-logical :: writebuf
-
-! 10 format(1x,2(f22.14,1x)) ! old
-10 format(1X,3(F20.16,1X))
-
-if ( .not. module_initialized ) call initialize_module
-
-! writing to a file (normal use) or to a character buffer?
-writebuf = present(charstring)
-
-! output file; test for ascii or binary, write what's asked, and return
-if (.not. writebuf) then
- if (ascii_file_format(fform)) then
- write(locfile, '(''loc3D'')' )
- write(locfile, 10) loc%x, loc%y, loc%z
- else
- write(locfile) loc%x, loc%y, loc%z
- endif
- return
-endif
-
-! you only get here if you're writing to a buffer and not
-! to a file, and you can't have binary format set.
-if (.not. ascii_file_format(fform)) then
- call error_handler(E_ERR, 'write_location', &
- 'Cannot use string buffer with binary format', &
- source, revision, revdate)
-endif
-
-! format the location to be more human-friendly; which in
-! this case doesn't change the value.
-
-! this must be the sum of the formats below.
-charlength = 38
-
-if (len(charstring) < charlength) then
- write(errstring, *) 'charstring buffer must be at least ', charlength, ' chars long'
- call error_handler(E_ERR, 'write_location', errstring, source, revision, revdate)
-endif
-
-write(charstring, '(A,F9.7,2(2X,F9.7))') 'X/Y/Z: ', loc%x, loc%y, loc%z
-
-
-end subroutine write_location
-
-!----------------------------------------------------------------------------
-
-function read_location(locfile, fform)
-
-! Reads a 3D location from locfile that was written by write_location.
-! See write_location for additional discussion.
-
-integer, intent(in) :: locfile
-type(location_type) :: read_location
-character(len = *), intent(in), optional :: fform
-
-character(len=5) :: header
-
-if ( .not. module_initialized ) call initialize_module
-
-if (ascii_file_format(fform)) then
- read(locfile, '(a5)' ) header
- if(header /= 'loc3D') then
- write(errstring,*)'Expected location header "loc3D" in input file, got ', header
- call error_handler(E_ERR, 'read_location', errstring, source, revision, revdate)
- endif
- ! Now read the location data value
- read(locfile, *) read_location%x, read_location%y, read_location%z
-else
- read(locfile) read_location%x, read_location%y, read_location%z
-endif
-
-end function read_location
-
-!--------------------------------------------------------------------------
-
-subroutine interactive_location(location, set_to_default)
-
-! Allows for interactive input of a location. Also gives option of selecting
-! a uniformly distributed random location.
-
-type(location_type), intent(out) :: location
-logical, intent(in), optional :: set_to_default
-
-real(r8) :: v(3)
-character(len=1) :: l(3)
-integer :: i
-
-if ( .not. module_initialized ) call initialize_module
-
-! If set_to_default is true, then just zero out and return
-if(present(set_to_default)) then
- if(set_to_default) then
- location%x = 0.0
- location%y = 0.0
- location%z = 0.0
- return
- endif
-endif
-
-l(1) = 'X'
-l(2) = 'Y'
-l(3) = 'Z'
-
-do i=1, 3
- write(*, *) 'Input ', l(i), ' location for this obs: value 0 to 1 or a negative number for '
- write(*, *) 'Uniformly distributed random location'
- read(*, *) v(i)
-
- do while(v(i) > 1.0_r8)
- write(*, *) 'Input value greater than 1.0 is illegal, please try again'
- read(*, *) v(i)
- end do
-
- if(v(i) < 0.0_r8) then
-
- ! Need to make sure random sequence is initialized
-
- if(.not. ran_seq_init) then
- call init_random_seq(ran_seq)
- ran_seq_init = .TRUE.
- endif
-
- ! Uniform location from 0 to 1 for this location type
-
- v(i) = random_uniform(ran_seq)
- write(*, *) 'random ',l(i),' location is ', v(i)
-
- endif
-enddo
-
-location%x = v(1)
-location%y = v(2)
-location%z = v(3)
-
-end subroutine interactive_location
-
-!----------------------------------------------------------------------------
-
-function nc_write_location_atts( ncFileID, fname, ObsNumDimID ) result (ierr)
-
-! Writes the "location module" -specific attributes to a netCDF file.
-
-use typeSizes
-use netcdf
-
-integer, intent(in) :: ncFileID ! handle to the netcdf file
-character(len=*), intent(in) :: fname ! file name (for printing purposes)
-integer, intent(in) :: ObsNumDimID ! handle to the dimension that grows
-integer :: ierr
-
-integer :: LocDimID
-integer :: VarID
-
-if ( .not. module_initialized ) call initialize_module
-
-ierr = -1 ! assume things will fail ...
-
-! define the rank/dimension of the location information
-call nc_check(nf90_def_dim(ncid=ncFileID, name='location', len=LocationDims, &
- dimid = LocDimID), 'nc_write_location_atts', 'def_dim:location '//trim(fname))
-
-! Define the observation location variable and attributes
-
-call nc_check(nf90_def_var(ncid=ncFileID, name='location', xtype=nf90_double, &
- dimids=(/ LocDimID, ObsNumDimID /), varid=VarID), &
- 'nc_write_location_atts', 'location:def_var')
-
-call nc_check(nf90_put_att(ncFileID, VarID, 'description', &
- 'location coordinates'), 'nc_write_location_atts', 'location:description')
-call nc_check(nf90_put_att(ncFileID, VarID, 'location_type', &
- trim(LocationName)), 'nc_write_location_atts', 'location:location_type')
-call nc_check(nf90_put_att(ncFileID, VarID, 'long_name', &
- trim(LocationLName)), 'nc_write_location_atts', 'location:long_name')
-call nc_check(nf90_put_att(ncFileID, VarID, 'storage_order', &
- 'X Y'), 'nc_write_location_atts', 'location:storage_order')
-call nc_check(nf90_put_att(ncFileID, VarID, 'units', &
- 'none none'), 'nc_write_location_atts', 'location:units')
-
-! no vertical array here.
-
-! If we made it to here without error-ing out ... we're good.
-
-ierr = 0
-
-end function nc_write_location_atts
-
-!----------------------------------------------------------------------------
-
-subroutine nc_get_location_varids( ncFileID, fname, LocationVarID, WhichVertVarID )
-
-! Return the LocationVarID and WhichVertVarID variables from a given netCDF file.
-!
-! ncFileId the netcdf file descriptor
-! fname the name of the netcdf file (for error messages only)
-! LocationVarID the integer ID of the 'location' variable in the netCDF file
-! WhichVertVarID the integer ID of the 'which_vert' variable in the netCDF file
-!
-! In this instance, WhichVertVarID will never be defined, ... set to a bogus value
-
-use typeSizes
-use netcdf
-
-integer, intent(in) :: ncFileID ! handle to the netcdf file
-character(len=*), intent(in) :: fname ! file name (for printing purposes)
-integer, intent(out) :: LocationVarID, WhichVertVarID
-
-if ( .not. module_initialized ) call initialize_module
-
-call nc_check(nf90_inq_varid(ncFileID, 'location', varid=LocationVarID), &
- 'nc_get_location_varids', 'inq_varid:location '//trim(fname))
-
-WhichVertVarID = -99
-
-end subroutine nc_get_location_varids
-
-!----------------------------------------------------------------------------
-
-subroutine nc_write_location(ncFileID, LocationVarID, loc, obsindex, WhichVertVarID)
-
-! Writes a SINGLE location to the specified netCDF variable and file.
-! The LocationVarID and WhichVertVarID must be the values returned from
-! the nc_get_location_varids call.
-
-use typeSizes
-use netcdf
-
-integer, intent(in) :: ncFileID, LocationVarID
-type(location_type), intent(in) :: loc
-integer, intent(in) :: obsindex
-integer, intent(in) :: WhichVertVarID
-
-real(r8), dimension(LocationDims) :: locations
-
-if ( .not. module_initialized ) call initialize_module
-
-locations = get_location( loc )
-
-call nc_check(nf90_put_var(ncFileID, LocationVarId, locations, &
- start=(/ 1, obsindex /), count=(/ LocationDims, 1 /) ), &
- 'nc_write_location', 'put_var:location')
-
-if ( WhichVertVarID >= 0 ) then
- write(errstring,*)'WhichVertVarID supposed to be negative ... is ',WhichVertVarID
- call error_handler(E_ERR, 'nc_write_location', errstring, source, revision, revdate)
-endif ! if less than zero (as it should be) ... just ignore
-
-end subroutine nc_write_location
-
-!----------------------------------------------------------------------------
-
-subroutine get_close_obs_init(gc, num, obs)
-
-! Initializes part of get_close accelerator that depends on the particular obs
-
-type(get_close_type), intent(inout) :: gc
-integer, intent(in) :: num
-type(location_type), intent(in) :: obs(num)
-
-! Set the value of num_obs in the structure
-gc%num = num
-
-end subroutine get_close_obs_init
-
-!----------------------------------------------------------------------------
-
-subroutine get_close_obs_destroy(gc)
-
-type(get_close_type), intent(inout) :: gc
-
-end subroutine get_close_obs_destroy
-
-!----------------------------------------------------------------------------
-
-subroutine get_close_maxdist_init(gc, maxdist, maxdist_list)
-
-type(get_close_type), intent(inout) :: gc
-real(r8), intent(in) :: maxdist
-real(r8), intent(in), optional :: maxdist_list(:)
-
-! Set the maximum distance in the structure
-gc%maxdist = maxdist
-
-end subroutine get_close_maxdist_init
-
-!----------------------------------------------------------------------------
-
-subroutine get_close_obs(gc, base_obs_loc, base_obs_kind, obs, obs_kind, &
- num_close, close_ind, dist)
-
-! Default version with no smarts; no need to be smart in 1D
-! Kinds are available here if one wanted to do more refined distances.
-
-type(get_close_type), intent(in) :: gc
-type(location_type), intent(in) :: base_obs_loc, obs(:)
-integer, intent(in) :: base_obs_kind, obs_kind(:)
-integer, intent(out) :: num_close, close_ind(:)
-real(r8), optional, intent(out) :: dist(:)
-
-integer :: i
-real(r8) :: this_dist
-
-! the list of locations in the obs() argument must be the same
-! as the list of locations passed into get_close_obs_init(), so
-! gc%num and size(obs) better be the same. if the list changes,
-! you have to destroy the old gc and init a new one.
-if (size(obs) /= gc%num) then
- write(errstring,*)'obs() array must match one passed to get_close_obs_init()'
- call error_handler(E_ERR, 'get_close_obs', errstring, source, revision, revdate)
-endif
-
-! Return list of obs that are within maxdist and their distances
-num_close = 0
-do i = 1, gc%num
- this_dist = get_dist(base_obs_loc, obs(i), base_obs_kind, obs_kind(i))
- if(this_dist <= gc%maxdist) then
- ! Add this ob to the list
- num_close = num_close + 1
- close_ind(num_close) = i
- if (present(dist)) dist(num_close) = this_dist
- endif
-end do
-
-end subroutine get_close_obs
-
-!----------------------------------------------------------------------------
-
-function is_location_in_region(loc, minl, maxl)
-
-! Returns true if the given location is between the other two.
-
-logical :: is_location_in_region
-type(location_type), intent(in) :: loc, minl, maxl
-
-if ( .not. module_initialized ) call initialize_module
-
-! assume failure and return as soon as we are confirmed right.
-! set to success only at the bottom after all tests have passed.
-is_location_in_region = .false.
-
-! FIXME: this is a triply cyclic domain. check if min
-! limit > max; if so, then wrap around.
-!if (minl%x <= maxl%x) .and. ...
-if ((loc%x < minl%x) .or. (loc%x > maxl%x)) return
-if ((loc%y < minl%y) .or. (loc%y > maxl%y)) return
-if ((loc%z < minl%z) .or. (loc%z > maxl%z)) return
-
-is_location_in_region = .true.
-
-end function is_location_in_region
-
-!----------------------------------------------------------------------------
-! stubs - always say no, but allow this code to be compiled with
-! common code that sometimes needs vertical info.
-!----------------------------------------------------------------------------
-
-function vert_is_undef(loc)
-
-! Stub, always returns false.
-
-logical :: vert_is_undef
-type(location_type), intent(in) :: loc
-
-vert_is_undef = .false.
-
-end function vert_is_undef
-
-!----------------------------------------------------------------------------
-
-function vert_is_surface(loc)
-
-! Stub, always returns false.
-
-logical :: vert_is_surface
-type(location_type), intent(in) :: loc
-
-vert_is_surface = .false.
-
-end function vert_is_surface
-
-!----------------------------------------------------------------------------
-
-function vert_is_pressure(loc)
-
-! Stub, always returns false.
-
-logical :: vert_is_pressure
-type(location_type), intent(in) :: loc
-
-vert_is_pressure = .false.
-
-end function vert_is_pressure
-
-!----------------------------------------------------------------------------
-
-function vert_is_height(loc)
-
-! Stub, always returns false.
-
-logical :: vert_is_height
-type(location_type), intent(in) :: loc
-
-vert_is_height = .false.
-
-end function vert_is_height
-
-!----------------------------------------------------------------------------
-
-function vert_is_level(loc)
-
-! Stub, always returns false.
-
-logical :: vert_is_level
-type(location_type), intent(in) :: loc
-
-vert_is_level = .false.
-
-end function vert_is_level
-
-!---------------------------------------------------------------------------
-
-function has_vertical_localization()
-
-! Always returns false since this type of location doesn't support
-! vertical localization.
-
-logical :: has_vertical_localization
-
-if ( .not. module_initialized ) call initialize_module
-
-has_vertical_localization = .false.
-
-end function has_vertical_localization
-
-
-!----------------------------------------------------------------------------
-! end of location/threed/location_mod.f90
-!----------------------------------------------------------------------------
-
-end module location_mod
Copied: DART/branches/development/location/threed/location_mod.f90 (from rev 5672, DART/branches/mpas/location/threed/location_mod.f90)
===================================================================
--- DART/branches/development/location/threed/location_mod.f90 (rev 0)
+++ DART/branches/development/location/threed/location_mod.f90 2012-04-09 17:11:37 UTC (rev 5673)
@@ -0,0 +1,729 @@
+! DART software - Copyright 2004 - 2011 UCAR. This open source software is
+! provided by UCAR, "as is", without charge, subject to all terms of use at
+! http://www.image.ucar.edu/DAReS/DART/DART_download
+
+module location_mod
+
+! <next few lines under version control, do not edit>
+! $URL$
+! $Id$
+! $Revision$
+! $Date$
+
+! Implements location interfaces for a three dimensional cyclic region.
+! The internal representation of the location is currently implemented
+! as (x, y) from 0.0 to 1.0 in both dimensions.
+!
+! If you are looking for a geophysical locations module, look at either
+! the threed_sphere or threed_cartesian versions of this file.
+
+use types_mod, only : r8, MISSING_R8
+use utilities_mod, only : register_module, error_handler, E_ERR, ascii_file_format, &
+ nc_check
+use random_seq_mod, only : random_seq_type, init_random_seq, random_uniform
+
+implicit none
+private
+
+public :: location_type, get_location, set_location, &
+ set_location_missing, is_location_in_region, &
+ write_location, read_location, interactive_location, query_location, &
+ LocationDims, LocationName, LocationLName, get_close_obs, &
+ get_close_maxdist_init, get_close_obs_init, get_close_type, &
+ operator(==), operator(/=), get_dist, get_close_obs_destroy, &
+ nc_write_location_atts, nc_get_location_varids, nc_write_location, &
+ vert_is_height, vert_is_pressure, vert_is_undef, vert_is_level, &
+ vert_is_surface, has_vertical_localization
+
+! version controlled file description for error handling, do not edit
+character(len=128), parameter :: &
+ source = "$URL$", &
+ revision = "$Revision$", &
+ revdate = "$Date$"
+
+type location_type
+ private
+ real(r8) :: x, y, z
+end type location_type
+
+! Needed as stub but not used in this low-order model
+type get_close_type
+ private
+ integer :: num
+ real(r8) :: maxdist
+end type get_close_type
+
+type(random_seq_type) :: ran_seq
+logical :: ran_seq_init = .false.
+logical, save :: module_initialized = .false.
+
+integer, parameter :: LocationDims = 3
+character(len = 129), parameter :: LocationName = "loc3D"
+character(len = 129), parameter :: LocationLName = "threed cyclic locations: x, y, z"
+
+character(len = 129) :: errstring
+
+interface operator(==); module procedure loc_eq; end interface
+interface operator(/=); module procedure loc_ne; end interface
+
+interface set_location
+ module procedure set_location_single
+ module procedure set_location_array
+end interface set_location
+
+contains
+
+!----------------------------------------------------------------------------
+
+subroutine initialize_module
+
+if (module_initialized) return
+
+call register_module(source, revision, revdate)
+module_initialized = .true.
+
+end subroutine initialize_module
+
+!----------------------------------------------------------------------------
+
+function get_dist(loc1, loc2, kind1, kind2)
+
+! Return the distance between 2 locations. Since this is a periodic
+! domain, the shortest distance may wrap around.
+
+type(location_type), intent(in) :: loc1, loc2
+integer, optional, intent(in) :: kind1, kind2
+real(r8) :: get_dist
+
+real(r8) :: x_dif, y_dif, z_dif
+
+if ( .not. module_initialized ) call initialize_module
+
+! Periodic domain, if distance is greater than half wraparound the other way.
+x_dif = abs(loc1%x - loc2%x)
+if (x_dif > 0.5_r8) x_dif = 1.0_r8 - x_dif
+y_dif = abs(loc1%y - loc2%y)
+if (y_dif > 0.5_r8) y_dif = 1.0_r8 - y_dif
+z_dif = abs(loc1%z - loc2%z)
+if (z_dif > 0.5_r8) z_dif = 1.0_r8 - z_dif
+
+get_dist = sqrt ( x_dif * x_dif + y_dif * y_dif + z_dif * z_dif)
+
+end function get_dist
+
+!---------------------------------------------------------------------------
+
+function loc_eq(loc1,loc2)
+
+! interface operator used to compare two locations.
+! Returns true only if all components are 'the same' to within machine
+! precision.
+
+type(location_type), intent(in) :: loc1, loc2
+logical :: loc_eq
+
+if ( .not. module_initialized ) call initialize_module
+
+loc_eq = .false.
+
+if ( abs(loc1%x - loc2%x ) > epsilon(loc1%x ) ) return
+if ( abs(loc1%y - loc2%y ) > epsilon(loc1%y ) ) return
+if ( abs(loc1%z - loc2%z ) > epsilon(loc1%z ) ) return
+
+loc_eq = .true.
+
+end function loc_eq
+
+!---------------------------------------------------------------------------
+
+function loc_ne(loc1,loc2)
+
+! interface operator used to compare two locations.
+! Returns true if locations are not identical to machine precision.
+
+type(location_type), intent(in) :: loc1, loc2
+logical :: loc_ne
+
+if ( .not. module_initialized ) call initialize_module
+
+loc_ne = (.not. loc_eq(loc1,loc2))
+
+end function loc_ne
+
+!---------------------------------------------------------------------------
+
+function get_location(loc)
+
+! Given a location type, return the x, y, z
+
+type(location_type), intent(in) :: loc
+real(r8), dimension(3) :: get_location
+
+if ( .not. module_initialized ) call initialize_module
+
+get_location(1) = loc%x
+get_location(2) = loc%y
+get_location(3) = loc%z
+
+end function get_location
+
+!----------------------------------------------------------------------------
+
+function set_location_single(x, y, z)
+
+! Given an x, y, z triplet, put the values into the location.
+
+real(r8), intent(in) :: x, y, z
+type (location_type) :: set_location_single
+
+if ( .not. module_initialized ) call initialize_module
+
+if(x < 0.0_r8 .or. x > 1.0_r8) then
+ write(errstring,*)'x (',x,') is not within range [0,1]'
+ call error_handler(E_ERR, 'set_location', errstring, source, revision, revdate)
+endif
+
+if(y < 0.0_r8 .or. y > 1.0_r8) then
+ write(errstring,*)'y (',y,') is not within range [0,1]'
+ call error_handler(E_ERR, 'set_location', errstring, source, revision, revdate)
+endif
+
+if(z < 0.0_r8 .or. z > 1.0_r8) then
+ write(errstring,*)'z (',z,') is not within range [0,1]'
+ call error_handler(E_ERR, 'set_location', errstring, source, revision, revdate)
+endif
+
+set_location_single%x = x
+set_location_single%y = y
+set_location_single%z = z
+
+end function set_location_single
+
+!----------------------------------------------------------------------------
+
+function set_location_array(list)
+
+! location semi-independent interface routine
+! given 3 float numbers, call the underlying set_location routine
+
+real(r8), intent(in) :: list(:)
+type (location_type) :: set_location_array
+
+if ( .not. module_initialized ) call initialize_module
+
+if (size(list) < 3) then
+ write(errstring,*)'requires 3 input values'
+ call error_handler(E_ERR, 'set_location', errstring, source, revision, revdate)
+endif
+
+set_location_array = set_location_single(list(1), list(2), list(3))
+
+end function set_location_array
+
+!----------------------------------------------------------------------------
+
+function set_location_missing()
+
+! fill in the contents to a known value.
+
+type (location_type) :: set_location_missing
+
+if ( .not. module_initialized ) call initialize_module
+
+set_location_missing%x = MISSING_R8
+set_location_missing%y = MISSING_R8
+set_location_missing%z = MISSING_R8
+
+end function set_location_missing
+
+!---------------------------------------------------------------------------
+
+function query_location(loc,attr)
+
+! Returns the value of the attribute
+!
+
+type(location_type), intent(in) :: loc
+character(len=*), optional, intent(in) :: attr
+real(r8) :: query_location
+
+if ( .not. module_initialized ) call initialize_module
+
+! see the long comment in this routine in the threed_sphere
+! module for warnings about compiler bugs before you change
+! this code.
+
+query_location = loc%x
+
+if (.not. present(attr)) return
+
+select case(attr)
+ case ('x','X')
+ query_location = loc%x
+ case ('y','Y')
+ query_location = loc%y
+ case ('z','Z')
+ query_location = loc%z
+ case default
+ call error_handler(E_ERR, 'query_location; threed', &
+ 'Only x, y, or z are legal attributes to request from location', source, revision, revdate)
+end select
+
+end function query_location
+
+!----------------------------------------------------------------------------
+
+subroutine write_location(locfile, loc, fform, charstring)
+
+! Writes a 3D location to the file.
+! additional functionality: if optional argument charstring is specified,
+! it must be long enough to hold the string, and the location information is
+! written into it instead of to a file. fform must be ascii (which is the
+! default if not specified) to use this option.
+
+integer, intent(in) :: locfile
+type(location_type), intent(in) :: loc
+character(len = *), intent(in), optional :: fform
+character(len = *), intent(out), optional :: charstring
+
+integer :: charlength
+logical :: writebuf
+
+! 10 format(1x,2(f22.14,1x)) ! old
+10 format(1X,3(F20.16,1X))
+
+if ( .not. module_initialized ) call initialize_module
+
+! writing to a file (normal use) or to a character buffer?
+writebuf = present(charstring)
+
+! output file; test for ascii or binary, write what's asked, and return
+if (.not. writebuf) then
+ if (ascii_file_format(fform)) then
+ write(locfile, '(''loc3D'')' )
+ write(locfile, 10) loc%x, loc%y, loc%z
+ else
+ write(locfile) loc%x, loc%y, loc%z
+ endif
+ return
+endif
+
+! you only get here if you're writing to a buffer and not
+! to a file, and you can't have binary format set.
+if (.not. ascii_file_format(fform)) then
+ call error_handler(E_ERR, 'write_location', &
+ 'Cannot use string buffer with binary format', &
+ source, revision, revdate)
+endif
+
+! format the location to be more human-friendly; which in
+! this case doesn't change the value.
+
+! this must be the sum of the formats below.
+charlength = 38
+
+if (len(charstring) < charlength) then
+ write(errstring, *) 'charstring buffer must be at least ', charlength, ' chars long'
+ call error_handler(E_ERR, 'write_location', errstring, source, revision, revdate)
+endif
+
+write(charstring, '(A,F9.7,2(2X,F9.7))') 'X/Y/Z: ', loc%x, loc%y, loc%z
+
+
+end subroutine write_location
+
+!----------------------------------------------------------------------------
+
+function read_location(locfile, fform)
+
+! Reads a 3D location from locfile that was written by write_location.
+! See write_location for additional discussion.
+
+integer, intent(in) :: locfile
+type(location_type) :: read_location
+character(len = *), intent(in), optional :: fform
+
+character(len=5) :: header
+
+if ( .not. module_initialized ) call initialize_module
+
+if (ascii_file_format(fform)) then
+ read(locfile, '(a5)' ) header
+ if(header /= 'loc3D') then
+ write(errstring,*)'Expected location header "loc3D" in input file, got ', header
+ call error_handler(E_ERR, 'read_location', errstring, source, revision, revdate)
+ endif
+ ! Now read the location data value
+ read(locfile, *) read_location%x, read_location%y, read_location%z
+else
+ read(locfile) read_location%x, read_location%y, read_location%z
+endif
+
+end function read_location
+
+!--------------------------------------------------------------------------
+
+subroutine interactive_location(location, set_to_default)
+
+! Allows for interactive input of a location. Also gives option of selecting
+! a uniformly distributed random location.
+
+type(location_type), intent(out) :: location
+logical, intent(in), optional :: set_to_default
+
+real(r8) :: v(3)
+character(len=1) :: l(3)
+integer :: i
+
+if ( .not. module_initialized ) call initialize_module
+
+! If set_to_default is true, then just zero out and return
+if(present(set_to_default)) then
+ if(set_to_default) then
+ location%x = 0.0
+ location%y = 0.0
+ location%z = 0.0
+ return
+ endif
+endif
+
+l(1) = 'X'
+l(2) = 'Y'
+l(3) = 'Z'
+
+do i=1, 3
+ write(*, *) 'Input ', l(i), ' location for this obs: value 0 to 1 or a negative number for '
+ write(*, *) 'Uniformly distributed random location'
+ read(*, *) v(i)
+
+ do while(v(i) > 1.0_r8)
+ write(*, *) 'Input value greater than 1.0 is illegal, please try again'
+ read(*, *) v(i)
+ end do
+
+ if(v(i) < 0.0_r8) then
+
+ ! Need to make sure random sequence is initialized
+
+ if(.not. ran_seq_init) then
+ call init_random_seq(ran_seq)
+ ran_seq_init = .TRUE.
+ endif
+
@@ Diff output truncated at 40000 characters. @@
More information about the Dart-dev
mailing list