psmile_neigh_nearx_irr2_3d_real.F90

Go to the documentation of this file.
00001 !-----------------------------------------------------------------------
00002 ! Copyright 2006-2010, NEC Europe Ltd., London, UK.
00003 ! All rights reserved. Use is subject to OASIS4 license terms.
00004 !-----------------------------------------------------------------------
00005 !BOP
00006 !
00007 ! !ROUTINE: PSMILe_Neigh_nearx_irr2_3d_real
00008 !
00009 ! !INTERFACE:
00010 
00011       subroutine psmile_neigh_nearx_irr2_3d_real (grid_id,        &
00012                       coords1, coords2, coords3,                  &
00013                       x_coords, y_coords, z_coords, coords_shape, &
00014                       mask_array, mask_shape, mask_available,     &
00015                       sin_values, cos_values, grid_valid_shape,   &
00016                       search_mode, srcloc, nloc, nprev, nsearch,  &
00017                       neighbors_3d, num_neigh,                    &
00018                       extra_search, ierror)
00019 !
00020 ! !USES:
00021 !
00022       use PRISM_constants
00023 !
00024       use PSMILe, dummy_interface => PSMILe_Neigh_nearx_irr2_3d_real
00025 
00026       Implicit none
00027 !
00028 ! !INPUT PARAMETERS:
00029 !
00030       Integer, Intent (In)            :: grid_id
00031 
00032 !     Link to the info on the component in which the donor cells
00033 !     should be searched.
00034 
00035       Integer, Intent (In)            :: nloc
00036 !
00037 !     Total number of locations to be transferred
00038 !
00039       Integer, Intent (In)            :: nprev
00040 !
00041 !     Number of locations for which neighbours were already searched.
00042 !
00043       Integer, Intent (In)            :: nsearch
00044 !
00045 !     Number of locations for which neighbours should be searched in this
00046 !     sweep.
00047 !
00048       Integer, Intent (In)            :: srcloc (ndim_3d, nloc)
00049 
00050 !     Indices of the grid cell in which the point was found.
00051 !
00052       Real, Intent (In)               :: coords1 (nloc)
00053       Real, Intent (In)               :: coords2 (nloc)
00054       Real, Intent (In)               :: coords3 (nloc)
00055 
00056 !     Coordinates of the target points (which were searched and found)
00057 
00058       Integer, Intent (In)            :: coords_shape (2, ndim_3d)
00059 
00060 !     Dimension of coordinates (method) x_coords, ...
00061 
00062       Real, Intent (In)               ::               
00063          x_coords(coords_shape(1,1):coords_shape(2,1), 
00064                   coords_shape(1,2):coords_shape(2,2))
00065       Real, Intent (In)               ::               
00066          y_coords(coords_shape(1,1):coords_shape(2,1), 
00067                   coords_shape(1,2):coords_shape(2,2))
00068       Real, Intent (In)               ::               
00069          z_coords(coords_shape(1,3):coords_shape(2,3))
00070 
00071 !     Coordinates of the method
00072 
00073       Integer, Intent (In)            :: mask_shape (2, ndim_3d)
00074 
00075 !     Dimension of (method) mask array
00076 
00077       Logical, Intent (In)            :: 
00078          mask_array (mask_shape (1,1):mask_shape (2,1), 
00079                      mask_shape (1,2):mask_shape (2,2), 
00080                      mask_shape (1,3):mask_shape (2,3))
00081 !     Mask of the method
00082 
00083       Logical, Intent (In)            :: mask_available
00084 
00085 !     Is mask specified in array "mask_array" ?
00086 !     mask_available = .false. : Mask is not available
00087 !                      .true.  : Mask is     available
00088 
00089       Integer, Intent (In)            :: grid_valid_shape (2, ndim_3d)
00090 
00091 !     Specifies the valid block shape for the "ndim_3d"-dimensional block
00092 
00093       Real, Intent (In)               :: 
00094          sin_values (grid_valid_shape(1,1):grid_valid_shape(2,1), 
00095                      grid_valid_shape(1,2):grid_valid_shape(2,2), 2)
00096 !
00097 !     Sin values of Longitudes and Latitudes (x_coords, y_coords)
00098 !
00099       Real, Intent (In)               :: 
00100          cos_values (grid_valid_shape(1,1):grid_valid_shape(2,1), 
00101                      grid_valid_shape(1,2):grid_valid_shape(2,2), 2)
00102 !    
00103 !     Cos values of Longitudes and Latitudes (x_coords, y_coords)
00104 !
00105       Integer, Intent (In)            :: search_mode
00106 !
00107 !     Specifies the search mode for nearest neigbours with
00108 !        search_mode = 3 : Full 3d-search
00109 !        search_mode = 2 : Search in 2d-hyperplane (1st 2 direction lon, lat)
00110 
00111       Integer, Intent (In)            :: num_neigh
00112 !
00113 !     Last dimension of neighbors array "neighbors_3d" and
00114 !     number of neighbors to be searched.
00115 !
00116 ! !INPUT/OUTPUT PARAMETERS:
00117 !
00118       Type (Extra_search_info)        :: extra_search
00119 !
00120 !     Number of locations where
00121 !       (*) required mask values were not true
00122 !
00123 ! !OUTPUT PARAMETERS:
00124 !
00125       Integer, Intent (Out)           :: neighbors_3d (ndim_3d, nloc, num_neigh)
00126 
00127 !     Indices of neighbor locations (interpolation bases)
00128 !
00129       Integer, Intent (Out)           :: ierror
00130 
00131 !     Returns the error code of PSMILE_Neigh_nearx_irr2_3d_real;
00132 !             ierror = 0 : No error
00133 !             ierror > 0 : Severe error
00134 !
00135 ! !DEFINED PARAMETERS:
00136 !
00137 !  lon   = Index of Longitudes in arrays "sin_values" and "cos_values"
00138 !  lat   = Index of Latitudes  in arrays "sin_values" and "cos_values"
00139 !
00140       Integer, Parameter              :: lon = 1
00141       Integer, Parameter              :: lat = 2
00142 !
00143 ! !LOCAL VARIABLES
00144 !
00145 !     ... For locations controlled
00146 !
00147       Real, Pointer                   :: sin_search (:, :)
00148       Real, Pointer                   :: cos_search (:, :)
00149       Real, Pointer                   :: z_search (:)
00150 !
00151 !     ... (1d) Indices in (entire) data to be extra searched
00152 !
00153       Integer                         :: i, ibeg, iend
00154       Integer                         :: j, jbeg, jend
00155 !
00156       Integer, Pointer                :: indices (:)
00157 !
00158 !     ... 3d Indices in local data to be extra searched
00159 !
00160       Integer                         :: sd
00161       Integer                         :: ijk0 (ndim_3d)
00162       Integer, Pointer                :: ijk (:, :)
00163       Integer                         :: width (ndim_3d)
00164 !
00165 !     ... for error handling
00166 !
00167       Integer, Parameter              :: nerrp = 2
00168       Integer                         :: ierrp (nerrp)
00169 !
00170 !
00171 ! !DESCRIPTION:
00172 !
00173 ! Subroutine "PSMILe_Neigh_nearx_irr2_3d_real" searches the next "num_neigh"
00174 ! nearest neighbours on the method-grid (x_coords, y_coords, z_coords)
00175 ! for the subgrid coords sent by the destination process.
00176 !
00177 ! There are 2 different of search modes supported by the routine:
00178 ! search_mode = 3: The nearest neighbours are searched are searched
00179 !                  within the 3d region.
00180 ! search_mode = 2: The nearest neighbours are searched only in the
00181 !                  2d hyperplane (lon- and lat-direction).
00182 !
00183 ! Subroutine "PSMILe_Neigh_nearx_irr2_3d_real" is designed for source grids of
00184 ! of type "PRISM_Irrlonlat_regvrt" and destination grids of type
00185 ! "PRISM_Irrlonlatvrt".
00186 !
00187 ! !REVISION HISTORY:
00188 !
00189 !   Date      Programmer   Description
00190 ! ----------  ----------   -----------
00191 ! 03.07.06    H. Ritzdorf  created
00192 !
00193 !EOP
00194 !----------------------------------------------------------------------
00195 !
00196 !  $Id: psmile_neigh_nearx_irr2_3d_real.F90 3018 2011-03-16 12:56:53Z redler $
00197 !  $Author: redler $
00198 !
00199    Character(len=len_cvs_string), save :: mycvs = 
00200        '$Id: psmile_neigh_nearx_irr2_3d_real.F90 3018 2011-03-16 12:56:53Z redler $'
00201 !                            
00202 !----------------------------------------------------------------------
00203 !
00204 !  Initialization
00205 !
00206 #ifdef VERBOSE
00207       print 9990, trim(ch_id)
00208 
00209       call psmile_flushstd
00210 #endif /* VERBOSE */
00211 !
00212       ierror = 0
00213 !
00214       ibeg = nprev + 1
00215       iend = nprev + nsearch
00216 !
00217       indices => extra_search%indices
00218 !
00219 !===> Get range in indices to be controlled
00220 !     Note: indices are locations in the full 3d cube and
00221 !           contains the previous numbers
00222 !           TODO: DAS SOLLTE in lokale Indices aendern
00223 !
00224       do jbeg = 1, extra_search%n_extra
00225          if (indices(jbeg) >= ibeg) exit
00226       end do
00227 !
00228       if (jbeg > extra_search%n_extra) then
00229 #ifdef VERBOSE
00230          print 9980, trim(ch_id), ierror
00231          call psmile_flushstd
00232 #endif /* VERBOSE */
00233          return
00234       endif
00235 !
00236       do jend = jbeg, extra_search%n_extra
00237          if (indices(jend) > iend) exit
00238       end do
00239       jend = jend - 1
00240 !
00241       if (jbeg > jend) then
00242 #ifdef VERBOSE
00243          print 9980, trim(ch_id), ierror
00244          call psmile_flushstd
00245 #endif /* VERBOSE */
00246          return
00247       endif
00248 !
00249       ibeg = indices(jbeg)
00250       iend = indices(jend)
00251 !
00252 #ifdef PRISM_ASSERTION
00253 !
00254 !===> Internal control
00255 !
00256       if (search_mode < 2 .or. search_mode > ndim_3d) then
00257          print *, trim(ch_id), " search_mode ", search_mode
00258          call psmile_assert (__FILE__, __LINE__, &
00259                              'search_mode is out of range ')
00260       endif
00261 !
00262       if (nloc < nprev+nsearch ) then
00263          print *, trim(ch_id), " nloc, nprev, nsearch ", nloc, nprev, nsearch
00264          call psmile_assert (__FILE__, __LINE__, &
00265                              'nloc < nprev + PRODUCT (nlocs) ')
00266       endif
00267 !
00268 !===> Are the locations within the correct shape ?
00269 !
00270 !cdir vector
00271          do j = jbeg, jend
00272          i = indices(j)
00273 !
00274          if (srcloc(1,i) < coords_shape(1,1) -1 .or. &
00275              srcloc(1,i) > coords_shape(2,1)    .or. &
00276              srcloc(2,i) < coords_shape(1,2) -1 .or. &
00277              srcloc(2,i) > coords_shape(2,2)    .or. &
00278              srcloc(3,i) < coords_shape(1,3) -1 .or. &
00279              srcloc(3,i) > coords_shape(2,3)) exit
00280          end do
00281 !
00282       if (j <= jend) then
00283          print *, "i, j, jbeg, jend", i, j, jbeg, jend
00284          print *, "srcloc ", srcloc(:, i)
00285          print *, "shape  ", coords_shape
00286          call psmile_assert (__FILE__, __LINE__, &
00287                              'wrong index')
00288       endif
00289 #endif
00290 !
00291 !-----------------------------------------------------------------------
00292 !     Allocate and set temporary array which is used to transform
00293 !     the coordinates of points to be searched.
00294 !     (*) Longitudes and Latitudes are transformed
00295 !         from degrees into radients.
00296 !     (*) The z-values are currently not transformed
00297 !         ??? Whats about PRISM_Irrlonlat_sigmavrt, ...
00298 !
00299 !     If global search is required, the values are stored in
00300 !     structure "extra_seach" in order to be reused.
00301 !-----------------------------------------------------------------------
00302 !
00303       if (extra_search%global_search) then
00304 
00305          sin_search => extra_search%sin_search_real
00306          cos_search => extra_search%cos_search_real
00307            z_search => extra_search%z_search_real
00308 
00309       else
00310 
00311          Allocate (sin_search(jbeg:jend, lat), &
00312                    cos_search(jbeg:jend, lat), &
00313                      z_search(jbeg:jend), STAT = ierror)
00314 
00315          if ( ierror > 0 ) then
00316             ierrp (1) = ierror
00317             ierrp (2) = (jend-jbeg+1) * (lat * 2 +  1)
00318 
00319             ierror = PRISM_Error_Alloc
00320             call psmile_error ( ierror, 'sin_search, ..., z_search', &
00321                                 ierrp, 2, __FILE__, __LINE__ )
00322             return
00323          endif
00324 
00325       endif
00326 !
00327       Allocate (ijk(ndim_3d, jbeg:jend), STAT = ierror)
00328 
00329       if ( ierror > 0 ) then
00330          ierrp (1) = ierror
00331          ierrp (2) = (jend-jbeg+1) * ndim_3d
00332 
00333          ierror = PRISM_Error_Alloc
00334          call psmile_error ( ierror, 'ijk', &
00335                              ierrp, 2, __FILE__, __LINE__ )
00336          return
00337       endif
00338 !
00339 !===> Compute values required to compute distances
00340 !
00341 !cdir vector
00342          do j = jbeg, jend
00343          sin_search (j, lon) = coords1 (indices(j)) * real_deg2rad
00344          end do
00345 !
00346 !cdir vector
00347          do j = jbeg, jend
00348          sin_search (j, lat) = coords2 (indices(j)) * real_deg2rad
00349          end do
00350 !
00351       cos_search (jbeg:jend, :) = cos (sin_search (jbeg:jend, :))
00352       sin_search (jbeg:jend, :) = sin (sin_search (jbeg:jend, :))
00353 !
00354 !cdir vector
00355          do j = jbeg, jend
00356          z_search (j) = coords3 (indices(j))
00357          end do
00358 !
00359 !-----------------------------------------------------------------------
00360 !     Compute coarse grid indices
00361 !-----------------------------------------------------------------------
00362 !
00363       ijk0 = Grids(grid_id)%ijk0
00364 !
00365 !===> Get coarse grid index on level nlev-2.
00366 !     By coarsening with a faktor of 4
00367 !
00368       sd = search_mode
00369       width (1:sd) = 3
00370 !
00371       if (search_mode == 2) then
00372          ijk (ndim_3d, jbeg:jend) = srcloc (ndim_3d, indices(jbeg:jend))
00373 !
00374          width (ndim_3d) = 0
00375       endif
00376 !
00377          do i = 1, sd
00378          ijk (i, jbeg:jend) = ijk0(i) + &
00379                               ((srcloc(i, indices(jbeg:jend)) - &
00380                                ijk0(i)) / 4) * 4
00381          end do ! i
00382 !
00383 !===> Compute nearest neigbours for all points
00384 !     specified in "indices(jbeg:jend)"
00385 !     TODO: sub-arrays "sin_search(jbeg:jend,:)" vermeiden
00386 !
00387       call psmile_neigh_nearx_sub_irr_real (grid_id,                  &
00388               x_coords, y_coords, z_coords, coords_shape,             &
00389               mask_array, mask_shape, mask_available,                 &
00390               sin_values, cos_values, grid_valid_shape, search_mode,  &
00391               neighbors_3d, num_neigh, nloc, extra_search,            &
00392               ijk, sin_search(jbeg:jend,:), cos_search(jbeg:jend,:),  &
00393                      z_search(jbeg:jend), jbeg, jend,                 &
00394               width, ierror)
00395 !
00396 !===> Deallocate memory
00397 !
00398       Deallocate (ijk)
00399 !
00400       if (.not. extra_search%global_search) then
00401          Deallocate (sin_search, cos_search, z_search)
00402       endif
00403 !
00404 !===> All done
00405 !
00406 #ifdef VERBOSE
00407       print 9980, trim(ch_id), ierror
00408 
00409       call psmile_flushstd
00410 #endif /* VERBOSE */
00411 !
00412 !  Formats:
00413 !
00414 9990 format (1x, a, ': psmile_neigh_nearx_irr2_3d_real')
00415 9980 format (1x, a, ': psmile_neigh_nearx_irr2_3d_real: eof, ierror =', i3)
00416 
00417       end subroutine PSMILe_Neigh_nearx_irr2_3d_real