Ignore:
Timestamp:
Jul 27, 2018 1:36:03 PM (6 years ago)
Author:
suehring
Message:

New Inifor features: grid stretching, improved command-interface, support start dates in different formats in both YYYYMMDD and YYYYMMDDHH, Ability to manually control input file prefixes (--radiation-prefix, --soil-preifx, --flow-prefix, --soilmoisture-prefix) for compatiblity with DWD forcast naming scheme; GNU-style short and long option; Prepared output of large-scale forcing profiles (no computation yet); Added preprocessor flag netcdf4 to switch output format between netCDF 3 and 4; Updated netCDF variable names and attributes to comply with PIDS v1.9; Inifor bugfixes: Improved compatibility with older Intel Intel compilers by avoiding implicit array allocation; Added origin_lon/_lat values and correct reference time in dynamic driver global attributes; corresponding PALM changes: adjustments to revised Inifor; variables names in dynamic driver adjusted; enable geostrophic forcing also in offline nested mode; variable names in LES-LES and COSMO offline nesting changed; lateral boundary flags for nesting, in- and outflow conditions renamed

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/UTIL/inifor/src/transform.f90

    r2718 r3182  
    2121! Current revisions:
    2222! -----------------
     23! Introduced new PALM grid stretching
     24! Removed unnecessary subroutine parameters
     25! Renamed kcur to k_intermediate
    2326!
    2427!
     
    8083       TYPE(grid_definition), INTENT(IN) ::  outgrid
    8184       REAL(dp), INTENT(IN)              ::  in_arr(0:,0:,0:)
    82        REAL(dp), INTENT(OUT)             ::  out_arr(0:,0:,0:)
    83 
    84        INTEGER :: i, j, k, l, nx, ny, nz
    85 
    86        nx = UBOUND(out_arr, 1)
    87        ny = UBOUND(out_arr, 2)
     85       REAL(dp), INTENT(OUT)             ::  out_arr(0:,0:,:)
     86
     87       INTEGER :: i, j, k, l, nz
     88
    8889       nz = UBOUND(out_arr, 3)
    8990
    90        DO j = 0, ny
    91        DO i = 0, nx
    92        DO k = nz, 0, -1
     91       DO j = LBOUND(out_arr, 2), UBOUND(out_arr, 2)
     92       DO i = LBOUND(out_arr, 1), UBOUND(out_arr, 1)
     93       DO k = nz, LBOUND(out_arr, 3), -1
    9394
    9495          ! TODO: Remove IF clause and extrapolate based on a critical vertical
     
    101102             out_arr(i,j,k) = 0.0_dp
    102103             DO l = 1, 2
    103                 out_arr(i,j,k) = out_arr(i,j,k) +                                 &
    104                     outgrid % w_verti(i,j,k,l) *                                  &
     104                out_arr(i,j,k) = out_arr(i,j,k) +                              &
     105                    outgrid % w_verti(i,j,k,l) *                               &
    105106                    in_arr(i,j,outgrid % kk(i,j,k, l) )
    106107             END DO
     
    139140    ! I index 0-based for the indices of the outvar to be consistent with the
    140141    ! outgrid indices and interpolation weights.
    141        TYPE(grid_definition), INTENT(IN) ::  outgrid
    142        REAL(dp), INTENT(IN)              ::  invar(0:,0:,0:)
    143        REAL(dp), INTENT(OUT)             ::  outvar(0:,0:,0:)
     142       TYPE(grid_definition), INTENT(IN)  ::  outgrid
     143       REAL(dp), INTENT(IN)               ::  invar(0:,0:,0:)
     144       REAL(dp), INTENT(OUT)              ::  outvar(0:,0:,0:)
    144145       TYPE(nc_var), INTENT(IN), OPTIONAL ::  ncvar
    145146
     
    413414
    414415    END SUBROUTINE rotate_to_cosmo
     416       
    415417
    416418
     
    427429!>                     -------------
    428430!>           jj, lat
    429 !>              ^      j
    430 !>              |       \          i
     431!>              ^        j
     432!>              |         \          i
    431433!>  jj(i,j,2/3) + ... 2 ---\--------/------ 3
    432434!>              |     | ^   \      /        |
     
    459461!>
    460462!------------------------------------------------------------------------------!
    461     SUBROUTINE find_horizontal_neighbours(cosmo_lat, cosmo_lon, cosmo_dxi,     &
    462        cosmo_dyi, palm_clat, palm_clon, palm_ii, palm_jj)
     463    SUBROUTINE find_horizontal_neighbours(cosmo_lat, cosmo_lon,                &
     464                                          palm_clat, palm_clon,                &
     465                                          palm_ii, palm_jj)
    463466
    464467       REAL(dp), DIMENSION(0:), INTENT(IN)        ::  cosmo_lat, cosmo_lon
    465468       REAL(dp), DIMENSION(0:,0:), INTENT(IN)     ::  palm_clat, palm_clon
    466        REAL(dp), INTENT(IN)                       ::  cosmo_dxi, cosmo_dyi
     469       REAL(dp)                                   ::  cosmo_dxi, cosmo_dyi
    467470       INTEGER, DIMENSION(0:,0:,1:), INTENT(OUT)  ::  palm_ii, palm_jj
    468471
     
    472475       lon0 = cosmo_lon(0)
    473476       lat0 = cosmo_lat(0)
     477       cosmo_dxi = 1.0_dp / (cosmo_lon(1) - cosmo_lon(0))
     478       cosmo_dyi = 1.0_dp / (cosmo_lat(1) - cosmo_lat(0))
    474479
    475480       DO j = 0, UBOUND(palm_clon, 2)!palm_grid % ny
     
    508513       TYPE(grid_definition), INTENT(IN)    ::  palm_intermediate
    509514
    510        INTEGER  ::  i, j, k, nx, ny, nz, nlev, kcur
     515       INTEGER  ::  i, j, k, nx, ny, nz, nlev, k_intermediate
    511516       LOGICAL  ::  point_is_below_grid, point_is_above_grid,                  &
    512517                    point_is_in_current_cell
     
    523528       DO j = 0, ny
    524529
    525           kcur = 0
     530          k_intermediate = 0
    526531
    527532          column_base = palm_intermediate % h(i,j,0)
     
    532537          ! cell, or above column_top. Keep increasing current cell index until
    533538          ! the current cell overlaps with the current_height.
    534           DO k = 0, nz
     539          DO k = 1, nz
    535540
    536541             ! Memorize the top and bottom boundaries of the coarse cell and the
    537542             ! current height within it
    538543             current_height = palm_grid % z(k) + palm_grid % z0
    539              h_top    = palm_intermediate % h(i,j,kcur+1)
    540              h_bottom = palm_intermediate % h(i,j,kcur)
     544             h_top    = palm_intermediate % h(i,j,k_intermediate+1)
     545             h_bottom = palm_intermediate % h(i,j,k_intermediate)
    541546
    542547             point_is_above_grid = (current_height > column_top) !22000m, very unlikely
     
    556561                palm_grid % w_verti(i,j,k,1:2) = - 2.0_dp
    557562
     563                message = "PALM-4U grid extends above COSMO-DE model top."
     564                CALL abort('find_vertical_neighbours_and_weights', message)
     565
    558566             ELSE IF (point_is_below_grid)  THEN
    559567
     
    564572                ! cycle through intermediate levels until current
    565573                ! intermediate-grid cell overlaps with current_height
    566                 DO WHILE (.NOT. point_is_in_current_cell .AND. kcur <= nlev-1)
    567                    kcur = kcur + 1
    568 
    569                    h_top    = palm_intermediate % h(i,j,kcur+1)
    570                    h_bottom = palm_intermediate % h(i,j,kcur)
     574                DO WHILE (.NOT. point_is_in_current_cell .AND. k_intermediate <= nlev-1)
     575                   k_intermediate = k_intermediate + 1
     576
     577                   h_top    = palm_intermediate % h(i,j,k_intermediate+1)
     578                   h_bottom = palm_intermediate % h(i,j,k_intermediate)
    571579                   point_is_in_current_cell = (                                &
    572580                      current_height >= h_bottom .AND.                         &
     
    575583                END DO
    576584
    577                 ! kcur = 48 indicates the last section (indices 48 and 49), i.e.
    578                 ! kcur = 49 is not the beginning of a valid cell.
    579                 IF (kcur > nlev-1)  THEN
    580                    message = "Index " // TRIM(str(kcur)) // " is above intermediate grid range."
     585                ! k_intermediate = 48 indicates the last section (indices 48 and 49), i.e.
     586                ! k_intermediate = 49 is not the beginning of a valid cell.
     587                IF (k_intermediate > nlev-1)  THEN
     588                   message = "Index " // TRIM(str(k_intermediate)) //          &
     589                             " is above intermediate grid range."
    581590                   CALL abort('find_vertical_neighbours', message)
    582591                END IF
    583592   
    584                 palm_grid % kk(i,j,k,1) = kcur
    585                 palm_grid % kk(i,j,k,2) = kcur + 1
     593                palm_grid % kk(i,j,k,1) = k_intermediate
     594                palm_grid % kk(i,j,k,2) = k_intermediate + 1
    586595
    587596                ! copmute vertical weights
     
    643652!         
    644653    SUBROUTINE compute_horizontal_interp_weights(cosmo_lat, cosmo_lon,         &
    645        cosmo_dxi, cosmo_dyi, palm_clat, palm_clon, palm_ii, palm_jj, palm_w_horiz)
     654       palm_clat, palm_clon, palm_ii, palm_jj, palm_w_horiz)
    646655       
    647656       REAL(dp), DIMENSION(0:), INTENT(IN)        ::  cosmo_lat, cosmo_lon
    648        REAL(dp), INTENT(IN)                       ::  cosmo_dxi, cosmo_dyi
     657       REAL(dp)                                   ::  cosmo_dxi, cosmo_dyi
    649658       REAL(dp), DIMENSION(0:,0:), INTENT(IN)     ::  palm_clat, palm_clon
    650659       INTEGER, DIMENSION(0:,0:,1:), INTENT(IN)   ::  palm_ii, palm_jj
     
    654663       REAL(dp) ::  wl, wp
    655664       INTEGER  ::  i, j
     665
     666       cosmo_dxi = 1.0_dp / (cosmo_lon(1) - cosmo_lon(0))
     667       cosmo_dyi = 1.0_dp / (cosmo_lat(1) - cosmo_lat(0))
    656668
    657669       DO j = 0, UBOUND(palm_clon, 2)
Note: See TracChangeset for help on using the changeset viewer.