Changeset 3866 for palm/trunk/UTIL/inifor/src/inifor_util.f90
- Timestamp:
- Apr 5, 2019 2:25:01 PM (5 years ago)
- File:
-
- 1 edited
-
palm/trunk/UTIL/inifor/src/inifor_util.f90 (modified) (17 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/UTIL/inifor/src/inifor_util.f90
r3785 r3866 26 26 ! ----------------- 27 27 ! $Id$ 28 ! Use PALM's working precision 29 ! Improved coding style 30 ! 31 ! 32 ! 3785 2019-03-06 10:41:14Z eckhard 28 33 ! Prefixed all INIFOR modules with inifor_ 29 34 ! … … 58 63 !> The util module provides miscellaneous utility routines for INIFOR. 59 64 !------------------------------------------------------------------------------! 60 #if defined ( __netcdf )61 65 MODULE inifor_util 62 66 63 67 USE inifor_defs, & 64 ONLY : dp, PI, DATE, SNAME68 ONLY : PI, DATE, SNAME, wp 65 69 USE inifor_types, & 66 70 ONLY : grid_definition … … 87 91 END TYPE 88 92 89 93 INTERFACE 90 94 91 95 !------------------------------------------------------------------------------! … … 95 99 !> structure. 96 100 !------------------------------------------------------------------------------! 97 98 99 100 101 102 103 104 105 106 101 FUNCTION strptime(string, format, timeinfo) BIND(c, NAME='strptime') 102 IMPORT :: C_CHAR, C_SIZE_T, tm_struct 103 104 IMPLICIT NONE 105 106 CHARACTER(KIND=C_CHAR), DIMENSION(*), INTENT(IN) :: string, format 107 TYPE(tm_struct), INTENT(OUT) :: timeinfo 108 109 INTEGER(C_SIZE_T) :: strptime 110 END FUNCTION 107 111 108 112 … … 113 117 !> structure to a string in the given 'format'. 114 118 !------------------------------------------------------------------------------! 115 116 117 118 119 120 121 122 123 124 125 126 119 FUNCTION strftime(string, string_len, format, timeinfo) BIND(c, NAME='strftime') 120 IMPORT :: C_CHAR, C_SIZE_T, tm_struct 121 122 IMPLICIT NONE 123 124 CHARACTER(KIND=C_CHAR), DIMENSION(*), INTENT(OUT) :: string 125 CHARACTER(KIND=C_CHAR), DIMENSION(*), INTENT(IN) :: format 126 INTEGER(C_SIZE_T), INTENT(IN) :: string_len 127 TYPE(tm_struct), INTENT(IN) :: timeinfo 128 129 INTEGER(C_SIZE_T) :: strftime 130 END FUNCTION 127 131 128 132 … … 135 139 !> e.g. increments the date if hours overfow 24. 136 140 !------------------------------------------------------------------------------! 137 138 139 140 141 142 143 144 145 146 147 141 FUNCTION mktime(timeinfo) BIND(c, NAME='mktime') 142 IMPORT :: C_PTR, tm_struct 143 144 IMPLICIT NONE 145 146 TYPE(tm_struct), INTENT(IN) :: timeinfo 147 148 TYPE(C_PTR) :: mktime 149 END FUNCTION 150 151 END INTERFACE 148 152 149 153 CONTAINS … … 156 160 !> format 157 161 !------------------------------------------------------------------------------! 158 CHARACTER(LEN=DATE) FUNCTION add_hours_to(date_string, hours) 159 CHARACTER(LEN=DATE), INTENT(IN) :: date_string 160 INTEGER, INTENT(IN) :: hours 161 162 CHARACTER(KIND=C_CHAR, LEN=*), PARAMETER :: format_string = "%Y%m%d%H" 163 CHARACTER(KIND=C_CHAR, LEN=DATE) :: c_date_string 164 TYPE(C_PTR) :: c_pointer 165 TYPE(tm_struct) :: time_info 166 INTEGER :: err 167 168 c_date_string = date_string 169 170 ! Convert C string to C tm struct 171 CALL init_tm(time_info) 172 err = strptime(c_date_string, format_string, time_info) 162 CHARACTER(LEN=DATE) FUNCTION add_hours_to(date_string, hours) 163 CHARACTER(LEN=DATE), INTENT(IN) :: date_string 164 INTEGER, INTENT(IN) :: hours 165 166 CHARACTER(KIND=C_CHAR, LEN=*), PARAMETER :: format_string = "%Y%m%d%H" 167 CHARACTER(KIND=C_CHAR, LEN=DATE) :: c_date_string 168 TYPE(C_PTR) :: c_pointer 169 TYPE(tm_struct) :: time_info 170 INTEGER :: err 171 172 c_date_string = date_string 173 174 ! Convert C string to C tm struct 175 CALL init_tm(time_info) 176 err = strptime(c_date_string, format_string, time_info) 177 178 ! Manipulate and normalize C tm struct 179 time_info%tm_hour = time_info%tm_hour + hours 180 c_pointer = mktime(time_info) 181 182 ! Convert back to C string 183 err = strftime(c_date_string, INT(DATE, KIND=C_SIZE_T), & 184 format_string, time_info) 185 186 add_hours_to = c_date_string 187 END FUNCTION 188 189 190 !------------------------------------------------------------------------------! 191 ! Description: 192 ! ------------ 193 !> Print all members of the given tm structure 194 !------------------------------------------------------------------------------! 195 SUBROUTINE print_tm(timeinfo) 196 TYPE(tm_struct), INTENT(IN) :: timeinfo 197 198 PRINT *, "sec: ", timeinfo%tm_sec, & !< seconds after the minute [0, 61] 199 "min: ", timeinfo%tm_min, & !< minutes after the hour [0, 59] 200 "hr: ", timeinfo%tm_hour, & !< hours since midnight [0, 23] 201 "day: ", timeinfo%tm_mday, & !< day of the month [1, 31] 202 "mon: ", timeinfo%tm_mon, & !< month since January [0, 11] 203 "yr: ", timeinfo%tm_year, & !< years since 1900 204 "wday:", timeinfo%tm_wday, & !< days since Sunday [0, 6] 205 "yday:", timeinfo%tm_yday, & !< days since January 1st [0, 356] 206 "dst: ", timeinfo%tm_isdst !< Daylight Saving time flag 207 END SUBROUTINE print_tm 208 173 209 174 ! Manipulate and normalize C tm struct175 time_info % tm_hour = time_info % tm_hour + hours176 c_pointer = mktime(time_info)177 178 ! Convert back to C string179 err = strftime(c_date_string, INT(DATE, KIND=C_SIZE_T), &180 format_string, time_info)181 182 add_hours_to = c_date_string183 END FUNCTION184 185 186 !------------------------------------------------------------------------------!187 ! Description:188 ! ------------189 !> Print all members of the given tm structure190 !------------------------------------------------------------------------------!191 SUBROUTINE print_tm(timeinfo)192 TYPE(tm_struct), INTENT(IN) :: timeinfo193 194 PRINT *, "sec: ", timeinfo % tm_sec, & !< seconds after the minute [0, 61]195 "min: ", timeinfo % tm_min, & !< minutes after the hour [0, 59]196 "hr: ", timeinfo % tm_hour, & !< hours since midnight [0, 23]197 "day: ", timeinfo % tm_mday, & !< day of the month [1, 31]198 "mon: ", timeinfo % tm_mon, & !< month since January [0, 11]199 "yr: ", timeinfo % tm_year, & !< years since 1900200 "wday:", timeinfo % tm_wday, & !< days since Sunday [0, 6]201 "yday:", timeinfo % tm_yday, & !< days since January 1st [0, 356]202 "dst: ", timeinfo % tm_isdst !< Daylight Saving time flag203 END SUBROUTINE print_tm204 205 206 210 !------------------------------------------------------------------------------! 207 211 ! Description: … … 209 213 !> Initialize the given tm structure with zero values 210 214 !------------------------------------------------------------------------------! 211 212 213 214 timeinfo %tm_sec = 0215 timeinfo %tm_min = 0216 timeinfo %tm_hour = 0217 timeinfo %tm_mday = 0218 timeinfo %tm_mon = 0219 timeinfo %tm_year = 0220 timeinfo %tm_wday = 0221 timeinfo %tm_yday = 0222 223 224 225 226 timeinfo %tm_isdst = -1227 215 SUBROUTINE init_tm(timeinfo) 216 TYPE(tm_struct), INTENT(INOUT) :: timeinfo 217 218 timeinfo%tm_sec = 0 219 timeinfo%tm_min = 0 220 timeinfo%tm_hour = 0 221 timeinfo%tm_mday = 0 222 timeinfo%tm_mon = 0 223 timeinfo%tm_year = 0 224 timeinfo%tm_wday = 0 225 timeinfo%tm_yday = 0 226 227 ! We use UTC times, so marking Daylight Saving Time (DST) 'not available' 228 ! (< 0). If this is set to 0, mktime will convert the timeinfo to DST and 229 ! add one hour. 230 timeinfo%tm_isdst = -1 231 END SUBROUTINE init_tm 228 232 229 233 … … 234 238 !> and stop 235 239 !------------------------------------------------------------------------------! 236 237 238 REAL(dp), INTENT(IN) :: start, stop239 REAL(dp), INTENT(INOUT) :: array(0:)240 241 242 243 244 245 246 247 248 249 250 251 array(i) = start + REAL(i, dp) / n * (stop - start)252 253 254 255 256 240 SUBROUTINE linspace(start, stop, array) 241 242 REAL(wp), INTENT(IN) :: start, stop 243 REAL(wp), INTENT(INOUT) :: array(0:) 244 INTEGER :: i, n 245 246 n = UBOUND(array, 1) 247 248 IF (n .EQ. 0) THEN 249 250 array(0) = start 251 252 ELSE 253 254 DO i = 0, n 255 array(i) = start + REAL(i, wp) / n * (stop - start) 256 ENDDO 257 258 ENDIF 259 260 END SUBROUTINE linspace 257 261 258 262 … … 263 267 !> (COSMO) to bottom-up (PALM) 264 268 !------------------------------------------------------------------------------! 265 266 267 REAL(dp), INTENT(INOUT) :: input_arr(:,:,:)268 269 270 271 269 SUBROUTINE reverse(input_arr) 270 271 REAL(wp), INTENT(INOUT) :: input_arr(:,:,:) 272 273 input_arr = input_arr(:,:,size(input_arr, 3):1:-1) 274 275 END SUBROUTINE reverse 272 276 273 277 … … 277 281 !> 278 282 !------------------------------------------------------------------------------! 279 280 281 REAL(dp), DIMENSION(:,:,:), INTENT(IN) :: avg_1, avg_2282 REAL(dp), INTENT(IN) :: t1, t2, t3283 REAL(dp), DIMENSION(:,:,:), INTENT(OUT) :: avg_3284 285 REAL(dp) :: ti283 SUBROUTINE deaverage(avg_1, t1, avg_2, t2, avg_3, t3) 284 285 REAL(wp), DIMENSION(:,:,:), INTENT(IN) :: avg_1, avg_2 286 REAL(wp), INTENT(IN) :: t1, t2, t3 287 REAL(wp), DIMENSION(:,:,:), INTENT(OUT) :: avg_3 288 289 REAL(wp) :: ti 286 290 287 ti = 1.0_dp / t3288 289 290 291 291 ti = 1.0_wp / t3 292 293 avg_3(:,:,:) = ti * ( t2 * avg_2(:,:,:) - t1 * avg_1(:,:,:) ) 294 295 END SUBROUTINE deaverage 292 296 293 297 … … 297 301 !> Compute the COSMO-DE/-D2 basic state pressure profile 298 302 !------------------------------------------------------------------------------! 299 300 301 REAL(dp), INTENT(IN) :: z(1:) !< height [m]302 REAL(dp), INTENT(IN) :: beta !< logarithmic lapse rate, dT / d ln(p) [K]303 REAL(dp), INTENT(IN) :: p_sl !< reference pressure [Pa]304 REAL(dp), INTENT(IN) :: t_sl !< reference tempereature [K]305 REAL(dp), INTENT(IN) :: rd !< ideal gas constant of dry air [J/kg/K]306 REAL(dp), INTENT(IN) :: g !< acceleration of Earth's gravity [m/s^2]307 REAL(dp), INTENT(OUT) :: p0(1:) !< COSMO basic state pressure [Pa]308 REAL(dp) :: root_frac, factor !< precomputed factors309 310 311 root_frac = (2.0_dp * beta * g) / (rd * t_sl*t_sl)312 313 314 factor * ( 1.0_dp - SQRT( 1.0_dp - root_frac * z(:) ) ) &315 316 317 303 SUBROUTINE get_basic_state(z, beta, p_sl, t_sl, rd, g, p0) 304 305 REAL(wp), INTENT(IN) :: z(1:) !< height [m] 306 REAL(wp), INTENT(IN) :: beta !< logarithmic lapse rate, dT / d ln(p) [K] 307 REAL(wp), INTENT(IN) :: p_sl !< reference pressure [Pa] 308 REAL(wp), INTENT(IN) :: t_sl !< reference tempereature [K] 309 REAL(wp), INTENT(IN) :: rd !< ideal gas constant of dry air [J/kg/K] 310 REAL(wp), INTENT(IN) :: g !< acceleration of Earth's gravity [m/s^2] 311 REAL(wp), INTENT(OUT) :: p0(1:) !< COSMO basic state pressure [Pa] 312 REAL(wp) :: root_frac, factor !< precomputed factors 313 314 factor = - t_sl / beta 315 root_frac = (2.0_wp * beta * g) / (rd * t_sl*t_sl) 316 317 p0(:) = p_sl * EXP( & 318 factor * ( 1.0_wp - SQRT( 1.0_wp - root_frac * z(:) ) ) & 319 ) 320 321 END SUBROUTINE get_basic_state 318 322 319 323 … … 326 330 !> theta = T * (p_ref/p)^(R/c_p) = T * e^( R/c_p * ln(p_ref/p) ) 327 331 !------------------------------------------------------------------------------! 328 329 REAL(dp), DIMENSION(:,:,:), INTENT(INOUT) :: t330 REAL(dp), DIMENSION(:,:,:), INTENT(IN) :: p331 REAL(dp), INTENT(IN) :: p_ref, r, cp332 REAL(dp) :: rcp333 334 335 336 337 332 SUBROUTINE potential_temperature(t, p, p_ref, r, cp) 333 REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) :: t 334 REAL(wp), DIMENSION(:,:,:), INTENT(IN) :: p 335 REAL(wp), INTENT(IN) :: p_ref, r, cp 336 REAL(wp) :: rcp 337 338 rcp = r/cp 339 t(:,:,:) = t(:,:,:) * EXP( rcp * LOG(p_ref / p(:,:,:)) ) 340 341 END SUBROUTINE potential_temperature 338 342 339 343 … … 343 347 !> Compute the density in place of the given temperature (t_rho). 344 348 !------------------------------------------------------------------------------! 345 SUBROUTINE moist_density(t_rho, p, qv, rd, rv) 346 REAL(dp), DIMENSION(:,:,:), INTENT(INOUT) :: t_rho 347 REAL(dp), DIMENSION(:,:,:), INTENT(IN) :: p, qv 348 REAL(dp), INTENT(IN) :: rd, rv 349 350 t_rho(:,:,:) = p(:,:,:) / ( & 351 (rv * qv(:,:,:) + rd * (1.0_dp - qv(:,:,:))) * t_rho(:,:,:) & 352 ) 353 354 END SUBROUTINE moist_density 355 356 357 ! Convert a real number to a string in scientific notation 358 ! showing four significant digits. 359 CHARACTER(LEN=SNAME) FUNCTION real_to_str(val, format) 360 361 REAL(dp), INTENT(IN) :: val 362 CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: format 363 364 IF (PRESENT(format)) THEN 365 WRITE(real_to_str, format) val 366 ELSE 367 WRITE(real_to_str, '(E11.4)') val 368 ENDIF 369 real_to_str = ADJUSTL(real_to_str) 370 371 END FUNCTION real_to_str 349 SUBROUTINE moist_density(t_rho, p, qv, rd, rv) 350 REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) :: t_rho 351 REAL(wp), DIMENSION(:,:,:), INTENT(IN) :: p, qv 352 REAL(wp), INTENT(IN) :: rd, rv 353 354 t_rho(:,:,:) = p(:,:,:) / ( & 355 (rv * qv(:,:,:) + rd * (1.0_wp - qv(:,:,:))) * t_rho(:,:,:) & 356 ) 357 358 END SUBROUTINE moist_density 359 360 !------------------------------------------------------------------------------! 361 ! Description: 362 ! ------------ 363 ! Convert a real number to a string in scientific notation showing four 364 ! significant digits. 365 !------------------------------------------------------------------------------! 366 CHARACTER(LEN=SNAME) FUNCTION real_to_str(val, format) 367 368 REAL(wp), INTENT(IN) :: val 369 CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: format 370 371 IF (PRESENT( format ) ) THEN 372 WRITE( real_to_str, format ) val 373 ELSE 374 WRITE( real_to_str, '(E11.4)' ) val 375 ENDIF 376 real_to_str = ADJUSTL( real_to_str ) 377 378 END FUNCTION real_to_str 372 379 373 380 … … 377 384 !> Converts the given real value to a string 378 385 !------------------------------------------------------------------------------! 379 380 381 REAL(dp), INTENT(IN) :: val382 383 384 385 386 386 CHARACTER(LEN=16) FUNCTION real_to_str_f(val) 387 388 REAL(wp), INTENT(IN) :: val 389 390 WRITE(real_to_str_f, '(F16.8)') val 391 real_to_str_f = ADJUSTL(real_to_str_f) 392 393 END FUNCTION real_to_str_f 387 394 388 395 … … 392 399 !> Converts the given integer value to a string 393 400 !------------------------------------------------------------------------------! 394 395 396 397 398 399 400 401 401 CHARACTER(LEN=10) FUNCTION str(val) 402 403 INTEGER, INTENT(IN) :: val 404 405 WRITE(str, '(i10)') val 406 str = ADJUSTL(str) 407 408 END FUNCTION str 402 409 403 410 … … 407 414 !> If the given path is not conlcuded by a slash, add one. 408 415 !------------------------------------------------------------------------------! 409 410 411 412 413 414 415 416 417 418 419 420 416 SUBROUTINE normalize_path(path) 417 418 CHARACTER(LEN=*), INTENT(INOUT) :: path 419 INTEGER :: n 420 421 n = LEN_TRIM(path) 422 423 IF (path(n:n) .NE. '/') THEN 424 path = TRIM(path) // '/' 425 ENDIF 426 427 END SUBROUTINE 421 428 422 429 END MODULE inifor_util 423 #endif424
Note: See TracChangeset
for help on using the changeset viewer.
