column_diagnostics/include/column_diagnostics.inc

   1 !***********************************************************************
   2 !*                   GNU Lesser General Public License
   3 !*
   4 !* This file is part of the GFDL Flexible Modeling System (FMS).
   5 !*
   6 !* FMS is free software: you can redistribute it and/or modify it under
   7 !* the terms of the GNU Lesser General Public License as published by
   8 !* the Free Software Foundation, either version 3 of the License, or (at
   9 !* your option) any later version.
  10 !*
  11 !* FMS is distributed in the hope that it will be useful, but WITHOUT
  12 !* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  13 !* FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  14 !* for more details.
  15 !*
  16 !* You should have received a copy of the GNU Lesser General Public
  17 !* License along with FMS.  If not, see <http://www.gnu.org/licenses/>.
  18 !***********************************************************************
  19 !> @file
  20
  21 !> @brief initialize_diagnostic_columns returns the (i, j, lat, lon) coord-
  22 !!    inates of any diagnostic columns that are located on the current
  23 !!    processor.
  24 subroutine INITIALIZE_DIAGNOSTIC_COLUMNS_     &
  25                    (module, num_diag_pts_latlon, num_diag_pts_ij,  &
  26                     global_i , global_j , global_lat_latlon,   &
  27                     global_lon_latlon, lonb_in, latb_in,  &
  28                     do_column_diagnostics,  &
  29                     diag_lon, diag_lat, diag_i, diag_j, diag_units)
  30
  31 !---------------------------------------------------------------------
  32 !    initialize_diagnostic_columns returns the (i, j, lat, lon) coord-
  33 !    inates of any diagnostic columns that are located on the current
  34 !    processor.
  35 !----------------------------------------------------------------------
  36
  37 !---------------------------------------------------------------------
  38 character(len=*),      intent(in)    :: module                !< module calling this subroutine
  39 integer,               intent(in)    :: num_diag_pts_latlon   !< number of diagnostic columns specified
  40                                                               !! by lat-lon  coordinates
  41 integer,               intent(in)    :: num_diag_pts_ij       !< number of diagnostic columns specified
  42                                                               !! by global (i,j) coordinates
  43 integer, dimension(:), intent(in)    :: global_i              !< specified global i coordinates
  44 integer, dimension(:), intent(in)    :: global_j              !< specified global j coordinates
  45 real(FMS_CD_KIND_), dimension(:), intent(in)   :: global_lat_latlon     !< specified global lat coordinates
  46 real(FMS_CD_KIND_), dimension(:), intent(in)   :: global_lon_latlon     !< specified global lon coordinates
  47 real(FMS_CD_KIND_), dimension(:,:), intent(in) :: lonb_in, latb_in
  48 logical, dimension(:,:), intent(out) :: do_column_diagnostics !< is a diagnostic column in this jrow ?
  49 integer, dimension(:), intent(inout) :: diag_i                !< processor i indices of diagnstic columns
  50 integer, dimension(:), intent(inout) :: diag_j                !< processor j indices of diagnstic columns
  51 real(FMS_CD_KIND_), dimension(:), intent(out) :: diag_lat              !< latitudes of diagnostic columns [ degrees ]
  52 real(FMS_CD_KIND_), dimension(:), intent(out) :: diag_lon              !< longitudes of diagnostic columns [ degrees ]
  53 integer, dimension(:), intent(out)   :: diag_units            !< unit number for each diagnostic column
  54 !---------------------------------------------------------------------
  55
  56 !---------------------------------------------------------------------
  57 !    intent(in) variables:
  58 !
  59 !       module                module calling this subroutine
  60 !       num_diag_pts_latlon   number of diagnostic columns specified
  61 !                             by lat-lon  coordinates
  62 !       num_diag_pts_ij       number of diagnostic columns specified
  63 !                             by global (i,j) coordinates
  64 !       global_i              specified global i coordinates
  65 !       global_j              specified global j coordinates
  66 !       global_lat_latlon     specified global lat coordinates
  67 !       global_lon_latlon     specified global lon coordinates
  68 !
  69 !    intent(out) variables:
  70 !
  71 !       do_column_diagnostics is a diagnostic column in this jrow ?
  72 !       diag_i                processor i indices of diagnstic columns
  73 !       diag_j                processor j indices of diagnstic columns
  74 !       diag_lat              latitudes of diagnostic columns
  75 !                             [ degrees ]
  76 !       diag_lon              longitudes of diagnostic columns
  77 !                             [ degrees ]
  78 !       diag_units            unit number for each diagnostic column
  79 !
  80 !---------------------------------------------------------------------
  81
  82 !--------------------------------------------------------------------
  83 !    local variables:
  84
  85       real(FMS_CD_KIND_), dimension(size(diag_i,1)) :: global_lat !< latitudes for all diagnostic columns [ degrees ]
  86       real(FMS_CD_KIND_), dimension(size(diag_i,1)) :: global_lon !< longitudes for all diagnostic columns [ degrees ]
  87       real(FMS_CD_KIND_), dimension(size(latb_in,1)-1, size(latb_in,2)-1) ::  &
  88                                   distance, distance_x, distance_y, &
  89                                   distance_x2, distance2
  90       real(FMS_CD_KIND_), dimension(size(latb_in,1), size(latb_in,2)) :: latb_deg
  91       real(FMS_CD_KIND_), dimension(size(lonb_in,1), size(lonb_in,2)) :: lonb_deg
  92       real(FMS_CD_KIND_) :: dellat, dellon
  93       real(FMS_CD_KIND_) :: latb_max, latb_min, lonb_max, lonb_min
  94
  95       integer            ::  num_diag_pts !< total number of diagnostic columns
  96       integer            ::  i            !< do loop indices
  97       integer            ::  j            !< do loop indices
  98       integer            ::  nn           !< do loop indices
  99       real(FMS_CD_KIND_) ::  ref_lat
 100       real(FMS_CD_KIND_) ::  current_distance
 101       character(len=8)   ::  char         !< character string for diaganostic column index
 102       character(len=32)  ::  filename     !< filename for output file for diagnostic column
 103       logical            ::  allow_ij_input
 104       logical            ::  open_file
 105       integer            ::  io
 106
 107       integer, parameter :: lkind=FMS_CD_KIND_
 108       real(FMS_CD_KIND_) :: tmp
 109 !--------------------------------------------------------------------
 110 !    local variables:
 111 !
 112 !       global_lat      latitudes for all diagnostic columns [ degrees ]
 113 !       global_lon      longitudes for all diagnostic columns
 114 !                       [ degrees ]
 115 !       num_diag_pts    total number of diagnostic columns
 116 !       i, j, nn        do loop indices
 117 !       char            character string for diaganostic column index
 118 !       filename        filename for output file for diagnostic column
 119 !
 120 !---------------------------------------------------------------------
 121
 122       if (.not. module_is_initialized) call column_diagnostics_init
 123
 124 !--------------------------------------------------------------------
 125 !    save the input lat and lon fields. define the delta of latitude
 126 !    and longitude.
 127 !--------------------------------------------------------------------
 128       latb_deg = real( real(latb_in,r8_kind)*RADIAN, FMS_CD_KIND_)  !< unit conversion in r8_kind
 129       lonb_deg = real( real(lonb_in,r8_kind)*RADIAN, FMS_CD_KIND_ ) !< unit conversion in r8_kind
 130       dellat = latb_in(1,2) - latb_in(1,1)
 131       dellon = lonb_in(2,1) - lonb_in(1,1)
 132       latb_max = MAXVAL (latb_deg(:,:))
 133       latb_min = MINVAL (latb_deg(:,:))
 134       lonb_max = MAXVAL (lonb_deg(:,:))
 135       lonb_min = MINVAL (lonb_deg(:,:))
 136       if (lonb_min < 10.0_lkind .or. lonb_max > 350.0_lkind) then
 137         lonb_min = 0.0_lkind
 138         lonb_max = 360.0_lkind
 139       endif
 140
 141       allow_ij_input = .true.
 142       ref_lat = latb_in(1,1)
 143       do i =2,size(latb_in,1)
 144         if (latb_in(i,1) /= ref_lat) then
 145           allow_ij_input = .false.
 146           exit
 147         endif
 148       end do
 149
 150       if ( .not. allow_ij_input .and. num_diag_pts_ij /= 0) then
 151         call error_mesg ('column_diagnostics_mod', &
 152         'cannot specify column diagnostics column with (i,j) &
 153            &coordinates when using cubed sphere -- must specify &
 154                                     & lat/lon coordinates', FATAL)
 155       endif
 156
 157 !----------------------------------------------------------------------
 158
 159 !----------------------------------------------------------------------
 160 !    initialize column_diagnostics flag and diag unit numbers. define
 161 !    total number of diagnostic columns.
 162 !----------------------------------------------------------------------
 163       do_column_diagnostics = .false.
 164       diag_units(:) = -1
 165       diag_i(:) = -99
 166       diag_j(:) = -99
 167       diag_lat(:) = -999.0_lkind
 168       diag_lon(:) = -999.0_lkind
 169       num_diag_pts = size(diag_i(:))
 170
 171 !--------------------------------------------------------------------
 172 !    define an array of lat-lon values for all diagnostic columns.
 173 !--------------------------------------------------------------------
 174       do nn = 1, num_diag_pts_latlon
 175         global_lat(nn) = global_lat_latlon(nn)
 176         global_lon(nn) = global_lon_latlon(nn)
 177       end do
 178
 179       do nn = 1, num_diag_pts_ij
 180         tmp = (-0.5_lkind*acos(-1.0_lkind) + 0.5_lkind*dellat) + real(global_j(nn)-1,FMS_CD_KIND_)*dellat
 181         global_lat(nn+num_diag_pts_latlon) = real( real(tmp,r8_kind)*RADIAN, FMS_CD_KIND_ )
 182         tmp = 0.5_lkind*dellon + real(global_i(nn)-1,FMS_CD_KIND_)*dellon
 183         global_lon(nn+num_diag_pts_latlon) = real( real(tmp,r8_kind)*RADIAN, FMS_CD_KIND_ )
 184      end do
 185
 186 !----------------------------------------------------------------------
 187 !    loop over all diagnostic points to check for their presence on
 188 !    this processor.
 189 !----------------------------------------------------------------------
 190       do nn=1,num_diag_pts
 191         open_file = .false.
 192
 193 !----------------------------------------------------------------------
 194 !    verify that the values of lat and lon are valid.
 195 !----------------------------------------------------------------------
 196         if (global_lon(nn) >= 0.0_lkind .and. &
 197             global_lon(nn) <= 360.0_lkind) then
 198         else
 199           call error_mesg ('column_diagnostics_mod', &
 200                ' invalid longitude', FATAL)
 201         endif
 202         if (global_lat(nn) >= -90.0_lkind .and. &
 203             global_lat(nn) <=  90.0_lkind) then
 204         else
 205           call error_mesg ('column_diagnostics_mod', &
 206                ' invalid latitude', FATAL)
 207         endif
 208
 209 !--------------------------------------------------------------------
 210 !    if the desired diagnostics column is within the current
 211 !    processor's domain, define the total and coordinate distances from
 212 !    each of the processor's grid points to the diagnostics point.
 213 !--------------------------------------------------------------------
 214
 215         if (global_lat(nn) .ge. latb_min .and.  &
 216             global_lat(nn) .le. latb_max) then
 217           if (global_lon(nn) .ge. lonb_min     .and.&
 218               global_lon(nn) .le. lonb_max)  then
 219             do j=1,size(latb_deg,2) - 1
 220               do i=1,size(lonb_deg,1) - 1
 221                 distance_y(i,j) = ABS(global_lat(nn) - latb_deg(i,j))
 222                 distance_x(i,j) = ABS(global_lon(nn) - lonb_deg(i,j))
 223                 distance_x2(i,j) = ABS((global_lon(nn)-360.0_lkind) - lonb_deg(i,j))
 224                 distance(i,j) = (global_lat(nn)-latb_deg(i,j))**2 + (global_lon(nn)-lonb_deg(i,j))**2
 225                 distance2(i,j) = (global_lat(nn)-latb_deg(i,j))**2 + ((global_lon(nn)-360.0_lkind) - lonb_deg(i,j))**2
 226               end do
 227             end do
 228
 229 !--------------------------------------------------------------------
 230 !    find the grid point on the processor that is within the specified
 231 !    critical distance and also closest to the requested diagnostics
 232 !    column. save the (i,j) coordinates and (lon,lat) of this model
 233 !    grid point. set a flag indicating that a disgnostics file should
 234 !    be opened on this processor for this diagnostic point.
 235 !--------------------------------------------------------------------
 236             current_distance = distance(1,1)
 237             do j=1,size(latb_deg,2) - 1
 238               do i=1,size(lonb_deg,1) - 1
 239                 if (distance_x(i,j) <= real(crit_xdistance,FMS_CD_KIND_) .and. &
 240                     distance_y(i,j) <= real(crit_ydistance,FMS_CD_KIND_)) then
 241                   if (distance(i,j) < current_distance) then
 242                     current_distance = distance(i,j)
 243                     do_column_diagnostics(i,j) = .true.
 244                     diag_j(nn) = j
 245                     diag_i(nn) = i
 246                     diag_lon(nn) = lonb_deg(i,j)
 247                     diag_lat(nn) = latb_deg(i,j)
 248                     open_file = .true.
 249                   endif
 250                 endif
 251
 252 !---------------------------------------------------------------------
 253 !    check needed because of the 0.0 / 360.0 longitude periodicity.
 254 !---------------------------------------------------------------------
 255                 if (distance_x2(i,j)<= real(crit_xdistance,FMS_CD_KIND_) .and. &
 256                     distance_y(i,j) <= real(crit_ydistance,FMS_CD_KIND_)) then
 257                   if (distance2(i,j) < current_distance) then
 258                     current_distance = distance2(i,j)
 259                     do_column_diagnostics(i,j) = .true.
 260                     diag_j(nn) = j
 261                     diag_i(nn) = i
 262                     diag_lon(nn) = lonb_deg(i,j)
 263                     diag_lat(nn) = latb_deg(i,j)
 264                     open_file = .true.
 265                   endif
 266                 endif
 267               end do
 268             end do
 269
 270 !--------------------------------------------------------------------
 271 !    if the point has been found on this processor, open a diagnostics
 272 !    file.
 273 !--------------------------------------------------------------------
 274             if (open_file) then
 275               write (char, '(i2)') nn
 276               filename = trim(module) // '_point' //    &
 277                          trim(adjustl(char)) // '.out'
 278               if(mpp_npes() > 10000) then
 279                  write( filename,'(a,i6.6)' )trim(filename)//'.', mpp_pe()-mpp_root_pe()
 280               else
 281                  write( filename,'(a,i4.4)' )trim(filename)//'.', mpp_pe()-mpp_root_pe()
 282               endif
 283               open(newunit=diag_units(nn), file=trim(filename), action='WRITE', position='rewind', iostat=io)
 284               if(io/=0) call error_mesg ('column_diagnostics_mod', 'Error in opening file '//trim(filename), FATAL)
 285             endif  ! (open_file)
 286           endif
 287         endif
 288       end do
 289
 290 !---------------------------------------------------------------------
 291
 292
 293 end subroutine INITIALIZE_DIAGNOSTIC_COLUMNS_
 294
 295
 296
 297
 298 !####################################################################
 299 !> @brief column_diagnostics_header writes out information concerning
 300 !!    time and location of following data into the column_diagnostics
 301 !!    output file.
 302 subroutine COLUMN_DIAGNOSTICS_HEADER_     &
 303                               (module, diag_unit, Time, nn, diag_lon, &
 304                                diag_lat, diag_i, diag_j)
 305
 306 !--------------------------------------------------------------------
 307 !    column_diagnostics_header writes out information concerning
 308 !    time and location of following data into the column_diagnostics
 309 !    output file.
 310 !--------------------------------------------------------------------
 311
 312 !--------------------------------------------------------------------
 313 character(len=*),      intent(in)  :: module    !< module name calling this subroutine
 314 type(time_type),       intent(in)  :: Time      !< current model time [ time_type ]
 315 integer,               intent(in)  :: diag_unit !< unit number for column_diagnostics output
 316 integer,               intent(in)  :: nn        !< index of diagnostic column currently active
 317 real(FMS_CD_KIND_),    dimension(:), intent(in)  :: diag_lon  !< longitude of current diagnostic column [ degrees ]
 318 real(FMS_CD_KIND_),    dimension(:), intent(in)  :: diag_lat  !< latitude of current diagnostic column [ degrees ]
 319 integer, dimension(:), intent(in)  :: diag_i    !< i coordinate of current diagnostic column
 320 integer, dimension(:), intent(in)  :: diag_j    !< j coordinate of current diagnostic column
 321
 322 !--------------------------------------------------------------------
 323 !    intent(in) variables
 324 !
 325 !       module     module name calling this subroutine
 326 !       Time       current model time [ time_type ]
 327 !       diag_unit  unit number for column_diagnostics output
 328 !       nn         index of diagnostic column currently active
 329 !       diag_lon   longitude of current diagnostic column [ degrees ]
 330 !       diag_lat   latitude of current diagnostic column [ degrees ]
 331 !       diag_i     i coordinate of current diagnostic column
 332 !       diag_j     j coordinate of current diagnostic column
 333 !
 334 !---------------------------------------------------------------------
 335
 336
 337 !--------------------------------------------------------------------
 338 !     local variables:
 339
 340       integer           :: year   !< integers defining the current time
 341       integer           :: month  !< integers defining the current time
 342       integer           :: day    !< integers defining the current time
 343       integer           :: hour   !< integers defining the current time
 344       integer           :: minute !< integers defining the current time
 345       integer           :: second !< integers defining the current time
 346       character(len=9)  :: mon    !< character string for the current month
 347       character(len=64) :: header !< title for the output
 348
 349 !--------------------------------------------------------------------
 350 !     local variables:
 351 !
 352 !       year, month, day, hour, minute, seconds
 353 !                      integers defining the current time
 354 !       mon            character string for the current month
 355 !       header         title for the output
 356 !
 357 !--------------------------------------------------------------------
 358
 359       if (.not. module_is_initialized) call column_diagnostics_init
 360
 361 !--------------------------------------------------------------------
 362 !    convert the time type to a date and time for printing. convert
 363 !    month to a character string.
 364 !--------------------------------------------------------------------
 365       call get_date (Time, year, month, day, hour, minute, second)
 366       mon = month_name(month)
 367
 368 !---------------------------------------------------------------------
 369 !    write timestamp and column location information to the diagnostic
 370 !    columns output unit.
 371 !---------------------------------------------------------------------
 372       write (diag_unit,'(a)')  ' '
 373       write (diag_unit,'(a)')  ' '
 374       write (diag_unit,'(a)')   &
 375               '======================================================'
 376       write (diag_unit,'(a)')  ' '
 377       header = '               PRINTING ' // module // '  DIAGNOSTICS'
 378       write (diag_unit,'(a)')  header
 379       write (diag_unit,'(a)')  ' '
 380       write (diag_unit,'(a, i6,2x, a,i4,i4,i4,i4)')  ' time stamp:',  &
 381                                            year, trim(mon), day, &
 382                                            hour, minute, second
 383       write (diag_unit,'(a, i4)')      &
 384             ' DIAGNOSTIC POINT COORDINATES, point #', nn
 385       write (diag_unit,'(a)')  ' '
 386       write (diag_unit,'(a,f8.3,a,f8.3)') ' longitude = ',    &
 387                    diag_lon(nn), ' latitude  = ', diag_lat(nn)
 388       write (diag_unit,'(a, i6, a,i6,a,i6)')    &
 389                                ' on processor # ', mpp_pe(),   &
 390                                ' :   processor i =', diag_i(nn),     &
 391                                ' ,   processor j =', diag_j(nn)
 392       write (diag_unit,'(a)')  ' '
 393
 394 !---------------------------------------------------------------------
 395
 396 end subroutine COLUMN_DIAGNOSTICS_HEADER_
 397 !@}