Changeset 4190 for palm/trunk/SOURCE
- Timestamp:
- Aug 27, 2019 3:42:37 PM (5 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 2 edited
-
netcdf_data_input_mod.f90 (modified) (4 diffs)
-
radiation_model_mod.f90 (modified) (16 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/netcdf_data_input_mod.f90
r4186 r4190 25 25 ! ----------------- 26 26 ! $Id$ 27 ! type real_1d changed to real_1d_3d 28 ! 29 ! 4186 2019-08-23 16:06:14Z suehring 27 30 ! Minor formatting adjustments 28 31 ! … … 396 399 END TYPE int_2d_32bit 397 400 ! 398 !-- Define data type to read 1D real variables399 TYPE real_1d 401 !-- Define data type to read 1D or 3D real variables. 402 TYPE real_1d_3d 400 403 LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used 401 404 405 INTEGER(iwp) :: lod = -1 !< level-of-detail 406 402 407 REAL(wp) :: fill = -9999.9_wp !< fill value 403 408 404 REAL(wp), DIMENSION(:), ALLOCATABLE :: var !< respective variable 405 END TYPE real_1d 409 REAL(wp), DIMENSION(:), ALLOCATABLE :: var1d !< respective 1D variable 410 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: var3d !< respective 3D variable 411 END TYPE real_1d_3d 406 412 ! 407 413 !-- Define data type to read 2D real variables … … 409 415 LOGICAL :: from_file = .FALSE. !< flag indicating whether an input variable is available and read from file or default values are used 410 416 417 INTEGER(iwp) :: lod !< level-of-detail 418 411 419 REAL(wp) :: fill = -9999.9_wp !< fill value 412 420 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: var !< respective variable … … 747 755 ! 748 756 !-- Public data structures 749 PUBLIC real_1d ,&757 PUBLIC real_1d_3d, & 750 758 real_2d 751 759 ! 752 760 !-- Public variables 753 761 PUBLIC albedo_pars_f, albedo_type_f, basal_area_density_f, buildings_f, & 754 building_id_f, building_pars_f, building_type_f, char_fill, & 762 building_id_f, building_pars_f, building_type_f, & 763 char_fill, & 764 char_lod, & 755 765 chem_emis, chem_emis_att, chem_emis_att_type, chem_emis_val_type, & 756 766 coord_ref_sys, & -
palm/trunk/SOURCE/radiation_model_mod.f90
r4188 r4190 28 28 ! ----------------- 29 29 ! $Id$ 30 ! Implement external radiation forcing also for level-of-detail = 2 31 ! (horizontally 2D radiation) 32 ! 33 ! 4188 2019-08-26 14:15:47Z suehring 30 34 ! Minor adjustment in error message 31 35 ! … … 235 239 water_type_f, & 236 240 char_fill, & 241 char_lod, & 237 242 check_existence, & 238 243 close_input_file, & … … 247 252 pids_id, & 248 253 open_read_file, & 249 real_1d ,&254 real_1d_3d, & 250 255 vars_pids 251 256 … … 780 785 REAL(wp), DIMENSION(:), ALLOCATABLE :: emiss_surf !< emissivity of the wall surface 781 786 ! 782 !-- External radiation 783 TYPE( real_1d ) :: rad_lw_in_f !< external incoming longwave radiation, from observation or model 784 TYPE( real_1d ) :: rad_sw_in_f !< external incoming shortwave radiation, from observation or model785 TYPE( real_1d ) :: rad_sw_in_dif_f !< external incoming shortwave radiation, diffuse part, from observation or model786 TYPE( real_1d ) :: time_rad_f !< time dimension for external radiation, from observation or model787 787 !-- External radiation. Depending on the given level of detail either a 1D or 788 !-- a 3D array will be allocated. 789 TYPE( real_1d_3d ) :: rad_lw_in_f !< external incoming longwave radiation, from observation or model 790 TYPE( real_1d_3d ) :: rad_sw_in_f !< external incoming shortwave radiation, from observation or model 791 TYPE( real_1d_3d ) :: rad_sw_in_dif_f !< external incoming shortwave radiation, diffuse part, from observation or model 792 TYPE( real_1d_3d ) :: time_rad_f !< time dimension for external radiation, from observation or model 788 793 789 794 INTERFACE radiation_check_data_output … … 1458 1463 INTEGER(iwp) :: l !< running index for orientation of vertical surfaces 1459 1464 INTEGER(iwp) :: m !< running index for surface elements 1460 INTEGER(iwp) :: ntime = 0 !< number of available external radiation timesteps 1465 INTEGER(iwp) :: ntime = 0 !< number of available external radiation timesteps 1461 1466 #if defined( __rrtmg ) 1462 1467 INTEGER(iwp) :: ind_type !< running index for subgrid-surface tiles … … 2501 2506 IF ( check_existence( vars_pids, 'time_rad' ) ) THEN 2502 2507 CALL netcdf_data_input_get_dimension_length( pids_id, & 2503 ntime, 'time_rad' ) 2508 ntime, & 2509 'time_rad' ) 2504 2510 2505 ALLOCATE( time_rad_f%var (0:ntime-1) )2511 ALLOCATE( time_rad_f%var1d(0:ntime-1) ) 2506 2512 ! 2507 2513 !-- Read variable 2508 CALL get_variable( pids_id, 'time_rad', time_rad_f%var )2514 CALL get_variable( pids_id, 'time_rad', time_rad_f%var1d ) 2509 2515 2510 2516 time_rad_f%from_file = .TRUE. … … 2513 2519 !-- Input shortwave downwelling. 2514 2520 IF ( check_existence( vars_pids, 'rad_sw_in' ) ) THEN 2515 ALLOCATE( rad_sw_in_f%var(0:ntime-1) )2516 2521 ! 2517 !-- Read_FillValue attribute2522 !-- Get _FillValue attribute 2518 2523 CALL get_attribute( pids_id, char_fill, rad_sw_in_f%fill, & 2519 2524 .FALSE., 'rad_sw_in' ) 2520 ! 2521 !-- Read variable 2522 CALL get_variable( pids_id, 'rad_sw_in', rad_sw_in_f%var ) 2523 2524 rad_sw_in_f%from_file = .TRUE. 2525 ! 2526 !-- Get level-of-detail 2527 CALL get_attribute( pids_id, char_lod, rad_sw_in_f%lod, & 2528 .FALSE., 'rad_sw_in' ) 2529 ! 2530 !-- Level-of-detail 1 - radiation depends only on time_rad 2531 IF ( rad_sw_in_f%lod == 1 ) THEN 2532 ALLOCATE( rad_sw_in_f%var1d(0:ntime-1) ) 2533 CALL get_variable( pids_id, 'rad_sw_in', rad_sw_in_f%var1d ) 2534 rad_sw_in_f%from_file = .TRUE. 2535 ! 2536 !-- Level-of-detail 2 - radiation depends on time_rad, y, x 2537 ELSEIF ( rad_sw_in_f%lod == 2 ) THEN 2538 ALLOCATE( rad_sw_in_f%var3d(0:ntime-1,nys:nyn,nxl:nxr) ) 2539 2540 CALL get_variable( pids_id, 'rad_sw_in', rad_sw_in_f%var3d, & 2541 nxl, nxr, nys, nyn, 0, ntime-1 ) 2542 2543 rad_sw_in_f%from_file = .TRUE. 2544 ELSE 2545 message_string = '"rad_sw_in" has no valid lod attribute' 2546 CALL message( 'radiation_init', 'PA0646', 1, 2, 0, 6, 0 ) 2547 ENDIF 2525 2548 ENDIF 2526 2549 ! 2527 2550 !-- Input longwave downwelling. 2528 2551 IF ( check_existence( vars_pids, 'rad_lw_in' ) ) THEN 2529 ALLOCATE( rad_lw_in_f%var(0:ntime-1) )2530 2552 ! 2531 !-- Read_FillValue attribute2553 !-- Get _FillValue attribute 2532 2554 CALL get_attribute( pids_id, char_fill, rad_lw_in_f%fill, & 2533 2555 .FALSE., 'rad_lw_in' ) 2534 ! 2535 !-- Read variable 2536 CALL get_variable( pids_id, 'rad_lw_in', rad_lw_in_f%var ) 2537 2538 rad_lw_in_f%from_file = .TRUE. 2556 ! 2557 !-- Get level-of-detail 2558 CALL get_attribute( pids_id, char_lod, rad_lw_in_f%lod, & 2559 .FALSE., 'rad_lw_in' ) 2560 ! 2561 !-- Level-of-detail 1 - radiation depends only on time_rad 2562 IF ( rad_lw_in_f%lod == 1 ) THEN 2563 ALLOCATE( rad_lw_in_f%var1d(0:ntime-1) ) 2564 CALL get_variable( pids_id, 'rad_lw_in', rad_lw_in_f%var1d ) 2565 rad_lw_in_f%from_file = .TRUE. 2566 ! 2567 !-- Level-of-detail 2 - radiation depends on time_rad, y, x 2568 ELSEIF ( rad_lw_in_f%lod == 2 ) THEN 2569 ALLOCATE( rad_lw_in_f%var3d(0:ntime-1,nys:nyn,nxl:nxr) ) 2570 2571 CALL get_variable( pids_id, 'rad_lw_in', rad_lw_in_f%var3d, & 2572 nxl, nxr, nys, nyn, 0, ntime-1 ) 2573 2574 rad_lw_in_f%from_file = .TRUE. 2575 ELSE 2576 message_string = '"rad_lw_in" has no valid lod attribute' 2577 CALL message( 'radiation_init', 'PA0646', 1, 2, 0, 6, 0 ) 2578 ENDIF 2539 2579 ENDIF 2540 2580 ! 2541 2581 !-- Input shortwave downwelling, diffuse part. 2542 2582 IF ( check_existence( vars_pids, 'rad_sw_in_dif' ) ) THEN 2543 ALLOCATE( rad_sw_in_dif_f%var(0:ntime-1) )2544 2583 ! 2545 2584 !-- Read _FillValue attribute 2546 2585 CALL get_attribute( pids_id, char_fill, rad_sw_in_dif_f%fill, & 2547 2586 .FALSE., 'rad_sw_in_dif' ) 2548 ! 2549 !-- Read variable 2550 CALL get_variable( pids_id, 'rad_sw_in_dif', rad_sw_in_dif_f%var ) 2551 2552 rad_sw_in_dif_f%from_file = .TRUE. 2587 ! 2588 !-- Get level-of-detail 2589 CALL get_attribute( pids_id, char_lod, rad_sw_in_dif_f%lod, & 2590 .FALSE., 'rad_sw_in_dif' ) 2591 ! 2592 !-- Level-of-detail 1 - radiation depends only on time_rad 2593 IF ( rad_sw_in_dif_f%lod == 1 ) THEN 2594 ALLOCATE( rad_sw_in_dif_f%var1d(0:ntime-1) ) 2595 CALL get_variable( pids_id, 'rad_sw_in_dif', & 2596 rad_sw_in_dif_f%var1d ) 2597 rad_sw_in_dif_f%from_file = .TRUE. 2598 ! 2599 !-- Level-of-detail 2 - radiation depends on time_rad, y, x 2600 ELSEIF ( rad_sw_in_dif_f%lod == 2 ) THEN 2601 ALLOCATE( rad_sw_in_dif_f%var3d(0:ntime-1,nys:nyn,nxl:nxr) ) 2602 2603 CALL get_variable( pids_id, 'rad_sw_in_dif', & 2604 rad_sw_in_dif_f%var3d, & 2605 nxl, nxr, nys, nyn, 0, ntime-1 ) 2606 2607 rad_sw_in_dif_f%from_file = .TRUE. 2608 ELSE 2609 message_string = '"rad_sw_in_dif" has no valid lod attribute' 2610 CALL message( 'radiation_init', 'PA0646', 1, 2, 0, 6, 0 ) 2611 ENDIF 2553 2612 ENDIF 2554 2613 ! … … 2573 2632 ENDIF 2574 2633 2575 IF ( ANY( rad_sw_in_f%var == rad_sw_in_f%fill ) ) THEN 2576 message_string = 'External radiation forcing: rad_sw_in must ' // & 2577 'not contain any fill values' 2578 CALL message( 'radiation_init', 'PA0197', 1, 2, 0, 6, 0 ) 2579 ENDIF 2580 2581 IF ( ANY( rad_lw_in_f%var == rad_lw_in_f%fill ) ) THEN 2582 message_string = 'External radiation forcing: rad_lw_in must ' // & 2583 'not contain any fill values' 2584 CALL message( 'radiation_init', 'PA0198', 1, 2, 0, 6, 0 ) 2585 ENDIF 2586 2587 IF ( rad_sw_in_dif_f%from_file ) THEN 2588 IF ( ANY( rad_sw_in_dif_f%var == rad_sw_in_dif_f%fill ) ) THEN 2589 message_string = 'External radiation forcing: ' // & 2590 'rad_lw_in_dif must not contain any fill ' // & 2591 'values' 2592 CALL message( 'radiation_init', 'PA0199', 1, 2, 0, 6, 0 ) 2593 ENDIF 2594 ENDIF 2595 2596 IF ( time_rad_f%var(0) /= 0.0_wp ) THEN 2634 IF ( time_rad_f%var1d(0) /= 0.0_wp ) THEN 2597 2635 message_string = 'External radiation forcing: first point in ' // & 2598 2636 'time is /= 0.0.' … … 2600 2638 ENDIF 2601 2639 2602 IF ( end_time - spinup_time > time_rad_f%var (ntime-1) ) THEN2640 IF ( end_time - spinup_time > time_rad_f%var1d(ntime-1) ) THEN 2603 2641 message_string = 'External radiation forcing does not cover ' // & 2604 2642 'the entire simulation time.' 2605 2643 CALL message( 'radiation_init', 'PA0314', 1, 2, 0, 6, 0 ) 2606 2644 ENDIF 2645 ! 2646 !-- Check for fill values in radiation 2647 IF ( ALLOCATED( rad_sw_in_f%var1d ) ) THEN 2648 IF ( ANY( rad_sw_in_f%var1d == rad_sw_in_f%fill ) ) THEN 2649 message_string = 'External radiation array "rad_sw_in" ' // & 2650 'must not contain any fill values.' 2651 CALL message( 'radiation_init', 'PA0197', 1, 2, 0, 6, 0 ) 2652 ENDIF 2653 ENDIF 2607 2654 2655 IF ( ALLOCATED( rad_lw_in_f%var1d ) ) THEN 2656 IF ( ANY( rad_lw_in_f%var1d == rad_lw_in_f%fill ) ) THEN 2657 message_string = 'External radiation array "rad_lw_in" ' // & 2658 'must not contain any fill values.' 2659 CALL message( 'radiation_init', 'PA0198', 1, 2, 0, 6, 0 ) 2660 ENDIF 2661 ENDIF 2662 2663 IF ( ALLOCATED( rad_sw_in_dif_f%var1d ) ) THEN 2664 IF ( ANY( rad_sw_in_dif_f%var1d == rad_sw_in_dif_f%fill ) ) THEN 2665 message_string = 'External radiation array "rad_sw_in_dif" ' //& 2666 'must not contain any fill values.' 2667 CALL message( 'radiation_init', 'PA0199', 1, 2, 0, 6, 0 ) 2668 ENDIF 2669 ENDIF 2670 2671 IF ( ALLOCATED( rad_sw_in_f%var3d ) ) THEN 2672 IF ( ANY( rad_sw_in_f%var3d == rad_sw_in_f%fill ) ) THEN 2673 message_string = 'External radiation array "rad_sw_in" ' // & 2674 'must not contain any fill values.' 2675 CALL message( 'radiation_init', 'PA0197', 1, 2, 0, 6, 0 ) 2676 ENDIF 2677 ENDIF 2678 2679 IF ( ALLOCATED( rad_lw_in_f%var3d ) ) THEN 2680 IF ( ANY( rad_lw_in_f%var3d == rad_lw_in_f%fill ) ) THEN 2681 message_string = 'External radiation array "rad_lw_in" ' // & 2682 'must not contain any fill values.' 2683 CALL message( 'radiation_init', 'PA0198', 1, 2, 0, 6, 0 ) 2684 ENDIF 2685 ENDIF 2686 2687 IF ( ALLOCATED( rad_sw_in_dif_f%var3d ) ) THEN 2688 IF ( ANY( rad_sw_in_dif_f%var3d == rad_sw_in_dif_f%fill ) ) THEN 2689 message_string = 'External radiation array "rad_sw_in_dif" ' //& 2690 'must not contain any fill values.' 2691 CALL message( 'radiation_init', 'PA0199', 1, 2, 0, 6, 0 ) 2692 ENDIF 2693 ENDIF 2694 ! 2695 !-- Currently, 2D external radiation input is not possible in 2696 !-- combination with topography where average radiation is used. 2697 IF ( ( rad_sw_in_f%lod == 2 .OR. rad_sw_in_f%lod == 2 .OR. & 2698 rad_sw_in_dif_f%lod == 2 ) .AND. average_radiation ) THEN 2699 message_string = 'External radiation with lod = 2 is currently '//& 2700 'not possible with average_radiation = .T..' 2701 CALL message( 'radiation_init', 'PA0670', 1, 2, 0, 6, 0 ) 2702 ENDIF 2703 ! 2704 !-- All radiation input should have the same level of detail. The sum 2705 !-- of lods divided by the number of available radiation arrays must be 2706 !-- 1 (if all are lod = 1) or 2 (if all are lod = 2). 2707 IF ( REAL( MERGE( rad_sw_in_f%lod, 0, rad_sw_in_f%from_file ) + & 2708 MERGE( rad_sw_in_f%lod, 0, rad_sw_in_f%from_file ) + & 2709 MERGE( rad_sw_in_dif_f%lod, 0, rad_sw_in_dif_f%from_file ),& 2710 KIND = wp ) / & 2711 ( MERGE( 1.0_wp, 0.0_wp, rad_sw_in_f%from_file ) + & 2712 MERGE( 1.0_wp, 0.0_wp, rad_sw_in_f%from_file ) + & 2713 MERGE( 1.0_wp, 0.0_wp, rad_sw_in_dif_f%from_file ) ) & 2714 /= 1.0_wp .AND. & 2715 REAL( MERGE( rad_sw_in_f%lod, 0, rad_sw_in_f%from_file ) + & 2716 MERGE( rad_sw_in_f%lod, 0, rad_sw_in_f%from_file ) + & 2717 MERGE( rad_sw_in_dif_f%lod, 0, rad_sw_in_dif_f%from_file ),& 2718 KIND = wp ) / & 2719 ( MERGE( 1.0_wp, 0.0_wp, rad_sw_in_f%from_file ) + & 2720 MERGE( 1.0_wp, 0.0_wp, rad_sw_in_f%from_file ) + & 2721 MERGE( 1.0_wp, 0.0_wp, rad_sw_in_dif_f%from_file ) ) & 2722 /= 2.0_wp ) THEN 2723 message_string = 'External radiation input should have the same '//& 2724 'lod.' 2725 CALL message( 'radiation_init', 'PA0673', 1, 2, 0, 6, 0 ) 2726 ENDIF 2727 2608 2728 ENDIF 2609 2729 ! … … 2690 2810 INTEGER(iwp) :: tm !< index of previous timestep 2691 2811 2692 REAL(wp) :: fac_dt !< interpolation factor 2693 REAL(wp) :: lw_in !< downwelling longwave radiation, interpolated value 2694 REAL(wp) :: sw_in !< downwelling shortwave radiation, interpolated value 2695 REAL(wp) :: sw_in_dif !< downwelling diffuse shortwave radiation, interpolated value 2812 REAL(wp) :: fac_dt !< interpolation factor 2696 2813 2697 2814 TYPE(surf_type), POINTER :: surf !< pointer on respective surface type, used to generalize routine … … 2708 2825 ELSE 2709 2826 t = 0 2710 DO WHILE ( time_rad_f%var (t) <= time_since_reference_point )2827 DO WHILE ( time_rad_f%var1d(t) <= time_since_reference_point ) 2711 2828 t = t + 1 2712 2829 ENDDO … … 2714 2831 tm = MAX( t-1, 0 ) 2715 2832 2716 fac_dt = ( time_since_reference_point - time_rad_f%var (tm) + dt_3d ) &2717 / ( time_rad_f%var (t) - time_rad_f%var(tm) )2833 fac_dt = ( time_since_reference_point - time_rad_f%var1d(tm) + dt_3d ) & 2834 / ( time_rad_f%var1d(t) - time_rad_f%var1d(tm) ) 2718 2835 fac_dt = MIN( 1.0_wp, fac_dt ) 2719 2836 ENDIF 2720 2721 sw_in = ( 1.0_wp - fac_dt ) * rad_sw_in_f%var(tm) &2722 + fac_dt * rad_sw_in_f%var(t)2723 2724 lw_in = ( 1.0_wp - fac_dt ) * rad_lw_in_f%var(tm) &2725 + fac_dt * rad_lw_in_f%var(t)2726 !2727 !-- Limit shortwave incoming radiation to positive values, in order2728 !-- to overcome possible observation errors.2729 sw_in = MAX( 0.0_wp, sw_in )2730 sw_in = MERGE( sw_in, 0.0_wp, sun_up )2731 2837 ! 2732 2838 !-- Call clear-sky calculation for each surface orientation. … … 2744 2850 CALL radiation_external_surf 2745 2851 ENDDO 2746 !2747 !-- In case of average radiation check if external diffuse shortwave2748 !-- radiation is available and, if required, interpolate it onto the2749 !-- respective arrays.2750 IF ( average_radiation ) THEN2751 IF ( rad_sw_in_dif_f%from_file ) THEN2752 sw_in_dif= ( 1.0_wp - fac_dt ) * rad_sw_in_dif_f%var(tm) &2753 + fac_dt * rad_sw_in_dif_f%var(t)2754 2755 IF ( ALLOCATED( rad_sw_in_diff ) ) rad_sw_in_diff = sw_in_dif2756 IF ( ALLOCATED( rad_sw_in_dir ) ) rad_sw_in_dir = sw_in - sw_in_dif2757 !2758 !-- Diffuse longwave radiation equals the total downwelling longwaver2759 !-- radiation2760 IF ( ALLOCATED( rad_lw_in_diff ) ) rad_lw_in_diff = lw_in2761 ENDIF2762 ENDIF2763 2852 2764 2853 CONTAINS … … 2770 2859 IMPLICIT NONE 2771 2860 2772 INTEGER(iwp) :: i !< grid index along x-dimension 2773 INTEGER(iwp) :: j !< grid index along y-dimension 2774 INTEGER(iwp) :: k !< grid index along z-dimension 2775 INTEGER(iwp) :: m !< running index for surface elements 2861 INTEGER(iwp) :: i !< grid index along x-dimension 2862 INTEGER(iwp) :: j !< grid index along y-dimension 2863 INTEGER(iwp) :: k !< grid index along z-dimension 2864 INTEGER(iwp) :: m !< running index for surface elements 2865 2866 REAL(wp) :: lw_in !< downwelling longwave radiation, interpolated value 2867 REAL(wp) :: sw_in !< downwelling shortwave radiation, interpolated value 2868 REAL(wp) :: sw_in_dif !< downwelling diffuse shortwave radiation, interpolated value 2776 2869 2777 2870 IF ( surf%ns < 1 ) RETURN 2871 ! 2872 !-- level-of-detail = 1. Note, here it must be distinguished between 2873 !-- averaged radiation and non-averaged radiation for the upwelling 2874 !-- fluxes. 2875 IF ( rad_sw_in_f%lod == 1 ) THEN 2876 2877 sw_in = ( 1.0_wp - fac_dt ) * rad_sw_in_f%var1d(tm) & 2878 + fac_dt * rad_sw_in_f%var1d(t) 2879 2880 lw_in = ( 1.0_wp - fac_dt ) * rad_lw_in_f%var1d(tm) & 2881 + fac_dt * rad_lw_in_f%var1d(t) 2882 ! 2883 !-- Limit shortwave incoming radiation to positive values, in order 2884 !-- to overcome possible observation errors. 2885 sw_in = MAX( 0.0_wp, sw_in ) 2886 sw_in = MERGE( sw_in, 0.0_wp, sun_up ) 2778 2887 2779 surf%rad_sw_in = sw_in 2780 surf%rad_lw_in = lw_in 2781 2782 IF ( average_radiation ) THEN 2888 surf%rad_sw_in = sw_in 2889 surf%rad_lw_in = lw_in 2783 2890 2784 surf%rad_sw_out = albedo_urb * surf%rad_sw_in 2891 IF ( average_radiation ) THEN 2892 surf%rad_sw_out = albedo_urb * surf%rad_sw_in 2893 2894 surf%rad_lw_out = emissivity_urb * sigma_sb * t_rad_urb**4 & 2895 + ( 1.0_wp - emissivity_urb ) * surf%rad_lw_in 2896 2897 surf%rad_net = surf%rad_sw_in - surf%rad_sw_out & 2898 + surf%rad_lw_in - surf%rad_lw_out 2899 2900 surf%rad_lw_out_change_0 = 4.0_wp * emissivity_urb & 2901 * sigma_sb & 2902 * t_rad_urb**3 2903 ELSE 2904 DO m = 1, surf%ns 2905 k = surf%k(m) 2906 surf%rad_sw_out(m) = ( surf%frac(ind_veg_wall,m) * & 2907 surf%albedo(ind_veg_wall,m) & 2908 + surf%frac(ind_pav_green,m) * & 2909 surf%albedo(ind_pav_green,m) & 2910 + surf%frac(ind_wat_win,m) * & 2911 surf%albedo(ind_wat_win,m) ) & 2912 * surf%rad_sw_in(m) 2913 2914 surf%rad_lw_out(m) = ( surf%frac(ind_veg_wall,m) * & 2915 surf%emissivity(ind_veg_wall,m) & 2916 + surf%frac(ind_pav_green,m) * & 2917 surf%emissivity(ind_pav_green,m) & 2918 + surf%frac(ind_wat_win,m) * & 2919 surf%emissivity(ind_wat_win,m) & 2920 ) & 2921 * sigma_sb & 2922 * ( surf%pt_surface(m) * exner(k) )**4 2923 2924 surf%rad_lw_out_change_0(m) = & 2925 ( surf%frac(ind_veg_wall,m) * & 2926 surf%emissivity(ind_veg_wall,m) & 2927 + surf%frac(ind_pav_green,m) * & 2928 surf%emissivity(ind_pav_green,m) & 2929 + surf%frac(ind_wat_win,m) * & 2930 surf%emissivity(ind_wat_win,m) & 2931 ) * 4.0_wp * sigma_sb & 2932 * ( surf%pt_surface(m) * exner(k) )**3 2933 ENDDO 2785 2934 2786 surf%rad_lw_out = emissivity_urb * sigma_sb * t_rad_urb**4 & 2787 + ( 1.0_wp - emissivity_urb ) * surf%rad_lw_in 2788 2789 surf%rad_net = surf%rad_sw_in - surf%rad_sw_out & 2790 + surf%rad_lw_in - surf%rad_lw_out 2791 2792 surf%rad_lw_out_change_0 = 4.0_wp * emissivity_urb * sigma_sb & 2793 * t_rad_urb**3 2794 ! 2795 !-- Calculate upwelling radiation fluxes and net radiation for each 2796 !-- surface element depending on its properties. 2935 ENDIF 2936 ! 2937 !-- If diffuse shortwave radiation is available, store it on 2938 !-- the respective files. 2939 IF ( rad_sw_in_dif_f%from_file ) THEN 2940 sw_in_dif= ( 1.0_wp - fac_dt ) * rad_sw_in_dif_f%var1d(tm) & 2941 + fac_dt * rad_sw_in_dif_f%var1d(t) 2942 2943 IF ( ALLOCATED( rad_sw_in_diff ) ) rad_sw_in_diff = sw_in_dif 2944 IF ( ALLOCATED( rad_sw_in_dir ) ) rad_sw_in_dir = sw_in & 2945 - sw_in_dif 2946 ! 2947 !-- Diffuse longwave radiation equals the total downwelling 2948 !-- longwave radiation 2949 IF ( ALLOCATED( rad_lw_in_diff ) ) rad_lw_in_diff = lw_in 2950 ENDIF 2951 ! 2952 !-- level-of-detail = 2 2797 2953 ELSE 2798 2954 … … 2801 2957 j = surf%j(m) 2802 2958 k = surf%k(m) 2959 2960 surf%rad_sw_in(m) = ( 1.0_wp - fac_dt ) & 2961 * rad_sw_in_f%var3d(tm,j,i) & 2962 + fac_dt * rad_sw_in_f%var3d(t,j,i) 2963 ! 2964 !-- Limit shortwave incoming radiation to positive values, in 2965 !-- order to overcome possible observation errors. 2966 surf%rad_sw_in(m) = MAX( 0.0_wp, surf%rad_sw_in(m) ) 2967 surf%rad_sw_in(m) = MERGE( surf%rad_sw_in(m), 0.0_wp, sun_up ) 2968 2969 surf%rad_lw_in(m) = ( 1.0_wp - fac_dt ) & 2970 * rad_lw_in_f%var3d(tm,j,i) & 2971 + fac_dt * rad_lw_in_f%var3d(t,j,i) 2803 2972 ! 2804 2973 !-- Weighted average according to surface fraction. … … 2833 3002 surf%rad_net(m) = surf%rad_sw_in(m) - surf%rad_sw_out(m) & 2834 3003 + surf%rad_lw_in(m) - surf%rad_lw_out(m) 3004 ! 3005 !-- If diffuse shortwave radiation is available, store it on 3006 !-- the respective files. 3007 IF ( rad_sw_in_dif_f%from_file ) THEN 3008 IF ( ALLOCATED( rad_sw_in_diff ) ) & 3009 rad_sw_in_diff(j,i) = ( 1.0_wp - fac_dt ) & 3010 * rad_sw_in_dif_f%var3d(tm,j,i) & 3011 + fac_dt * rad_sw_in_dif_f%var3d(t,j,i) 3012 ! 3013 !-- dir = sw_in - sw_in_dif. 3014 IF ( ALLOCATED( rad_sw_in_dir ) ) & 3015 rad_sw_in_dir(j,i) = surf%rad_sw_in(m) - & 3016 rad_sw_in_diff(j,i) 3017 ! 3018 !-- Diffuse longwave radiation equals the total downwelling 3019 !-- longwave radiation 3020 IF ( ALLOCATED( rad_lw_in_diff ) ) & 3021 rad_lw_in_diff = surf%rad_lw_in(m) 3022 ENDIF 2835 3023 2836 3024 ENDDO 2837 3025 2838 3026 ENDIF 2839 2840 ! 2841 !-- Fill out values in radiation arrays3027 ! 3028 !-- Store radiation also on 2D arrays, which are still used for 3029 !-- direct-diffuse splitting. 2842 3030 DO m = 1, surf%ns 2843 3031 i = surf%i(m) 2844 3032 j = surf%j(m) 3033 3034 rad_sw_in(0,:,:) = surf%rad_sw_in(m) 3035 rad_lw_in(0,:,:) = surf%rad_lw_in(m) 2845 3036 rad_sw_out(0,j,i) = surf%rad_sw_out(m) 2846 3037 rad_lw_out(0,j,i) = surf%rad_lw_out(m) 2847 3038 ENDDO 2848 2849 rad_sw_in(0,:,:) = sw_in2850 rad_lw_in(0,:,:) = lw_in2851 3039 2852 3040 END SUBROUTINE radiation_external_surf
Note: See TracChangeset
for help on using the changeset viewer.
