- Timestamp:
- Aug 4, 2025, 3:03:07 PM (11 days ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/contrails/libf/phylmd/lmdz_cloud_optics_prop.f90
r5790 r5796 4 4 CONTAINS 5 5 6 SUBROUTINE cloud_optics_prop(klon, klev, paprs, pplay, temp, radocond, picefra, pcl c, &6 SUBROUTINE cloud_optics_prop(klon, klev, paprs, pplay, temp, radocond, picefra, pclf, pclc, & 7 7 pcltau, pclemi, pch, pcl, pcm, pct, radocondwp, xflwp, xfiwp, xflwc, xfiwc, & 8 8 mass_solu_aero, mass_solu_aero_pi, pcldtaupi, distcltop, temp_cltop, re, fl, reliq, reice, & … … 11 11 icefrac_optics, dNovrN, ptconv, rnebcon, ccwcon, & 12 12 !--AB contrails 13 contfra , qice_cont, Nice_cont, pclc_nocont, pcltau_nocont, &13 contfravol, contfra, qice_cont, Nice_cont, pclc_cont, pcltau_nocont, & 14 14 pclemi_nocont, pcltau_cont, pclemi_cont, pch_nocont, pct_cont, & 15 xfiwp_nocont, xfiwc_nocont, reice_nocont) 15 xfiwp_cont, xfiwc_cont, reice_cont, & 16 missing_val) 16 17 17 18 USE lmdz_cloud_optics_prop_ini , ONLY : flag_aerosol, ok_cdnc … … 33 34 USE lmdz_cloud_optics_prop_ini , ONLY : rei_coef, rei_min_temp 34 35 USE lmdz_cloud_optics_prop_ini , ONLY : zepsec, novlp, iflag_ice_thermo, ok_new_lscp 35 USE lmdz_cloud_optics_prop_ini , ONLY : ok_plane_contrail, eff2vol_radius_contrails, rho_ice 36 USE lmdz_cloud_optics_prop_ini , ONLY : ok_plane_contrail, rho_ice 37 USE lmdz_cloud_optics_prop_ini , ONLY : eff2vol_radius_contrails, rei_min_contrails 36 38 37 39 … … 59 61 ! input: 60 62 INTEGER, INTENT(IN) :: klon, klev ! number of horizontal and vertical grid points 63 REAL, INTENT(IN) :: missing_val 61 64 REAL, INTENT(IN) :: paprs(klon, klev+1)! pressure at bottom interfaces [Pa] 62 65 REAL, INTENT(IN) :: pplay(klon, klev) ! pressure at the middle of layers [Pa] … … 73 76 REAL, INTENT(OUT) :: distcltop(klon,klev) ! distance from large scale cloud top [m] 74 77 REAL, INTENT(OUT) :: temp_cltop(klon,klev)!temperature at large scale cloud top [K] 78 REAL, INTENT(IN) :: pclf(klon, klev) ! cloud fraction for radiation [-] 75 79 76 80 LOGICAL, INTENT(IN) :: ptconv(klon, klev) ! flag for grid points affected by deep convection 77 81 78 82 ! inout: 79 REAL, INTENT(INOUT) :: pclc(klon, klev) ! cloud fractionfor radiation [-]83 REAL, INTENT(INOUT) :: pclc(klon, klev) ! cloud cover for radiation [-] 80 84 81 85 ! out: … … 119 123 120 124 !--AB for contrails. All these are used / outputed only if ok_plane_contrail=y 121 REAL, INTENT(IN) :: contfra(klon, klev) ! contrails fraction [-] 125 REAL, INTENT(IN) :: contfravol(klon, klev) ! contrails volumic fraction [-] 126 REAL, INTENT(IN) :: contfra(klon, klev) ! contrails fraction for radiation [-] 122 127 REAL, INTENT(IN) :: qice_cont(klon, klev) ! contrails condensed water [kg/kg] 123 128 REAL, INTENT(IN) :: Nice_cont(klon, klev) ! contrails ice crystals concentration [#/kg] 124 129 REAL, INTENT(OUT) :: pch_nocont(klon) ! 2D high cloud cover without contrails[-] 125 130 REAL, INTENT(OUT) :: pct_cont(klon) ! 2D total contrails cover[-] 126 REAL, INTENT(OUT) :: xfiwp_ nocont(klon) ! ice water path (seen by radiation) withoutcontrails [kg/m2]127 REAL, INTENT(OUT) :: xfiwc_ nocont(klon, klev) ! ice water content seen by radiation withoutcontrails [kg/kg]128 REAL, INTENT(OUT) :: pclc_ nocont(klon, klev) ! cloud fraction for radiation withoutcontrails [-]131 REAL, INTENT(OUT) :: xfiwp_cont(klon) ! ice water path (seen by radiation) of contrails [kg/m2] 132 REAL, INTENT(OUT) :: xfiwc_cont(klon, klev) ! ice water content seen by radiation of contrails [kg/kg] 133 REAL, INTENT(OUT) :: pclc_cont(klon, klev) ! cloud fraction for radiation of contrails [-] 129 134 REAL, INTENT(OUT) :: pcltau_nocont(klon, klev) ! cloud optical depth without contrails [-] 130 135 REAL, INTENT(OUT) :: pclemi_nocont(klon, klev) ! cloud emissivity without contrails [-] 131 136 REAL, INTENT(OUT) :: pcltau_cont(klon, klev) ! contrails optical depth [-] 132 137 REAL, INTENT(OUT) :: pclemi_cont(klon, klev) ! contrails emissivity [-] 133 REAL, INTENT(OUT) :: reice_ nocont(klon, klev) ! ice effective radius withoutcontrails [microns]138 REAL, INTENT(OUT) :: reice_cont(klon, klev) ! ice effective radius of contrails [microns] 134 139 !--AB 135 140 … … 172 177 REAL zflwp_var, zfiwp_var 173 178 REAL d_rei_dt 174 REAL pclc_cont(klon,klev) 179 REAL pclc_nocont(klon,klev) 180 REAL pclf_nocont(klon,klev) 181 REAL xfiwc_nocont(klon,klev) 175 182 REAL mice_cont, rei_cont 176 183 … … 194 201 xfiwc = 0.D0 195 202 !--AB 196 IF ( ok_plane_contrail ) THEN 197 xfiwp_nocont = 0.D0 198 xfiwc_nocont = 0.D0 199 reice_nocont = 0. 200 ENDIF 203 xfiwp_cont = 0.D0 204 xfiwc_cont = 0.D0 205 reice_cont = 1. 201 206 202 207 reliq = 0. … … 254 259 DO k = 1, klev 255 260 DO i = 1, klon 261 pclf_nocont(i,k) = MAX(0., pclf(i, k) - contfravol(i, k)) 256 262 pclc_nocont(i,k) = MAX(0., pclc(i, k) - contfra(i, k)) 257 263 xfiwc_nocont(i, k) = MAX(0., xfiwc(i, k) - qice_cont(i, k)) 264 xfiwc_cont(i,k) = qice_cont(i,k) 258 265 ENDDO 259 266 ENDDO … … 343 350 IF (iflag_rei .EQ. 2) THEN 344 351 ! in-cloud ice water content in g/m3 345 iwc = icefrac_optics(i,k) * radocond(i,k) / pclc(i,k) * zrho(i,k) * 1000. 352 IF ( ptconv(i,k) ) THEN 353 !--Needed because pclf does not contain convective clouds (should be fixed...) 354 iwc = icefrac_optics(i,k) * radocond(i,k) / pclc(i,k) * zrho(i,k) * 1000. 355 ELSE 356 iwc = icefrac_optics(i,k) * radocond(i,k) / pclf(i,k) * zrho(i,k) * 1000. 357 ENDIF 346 358 ! this formula is a simplified version of the Sun 2001 one (as in the IFS model, 347 359 ! and which is activated for iflag_rei = 1). … … 355 367 dei = rei_coef * (iwc**0.2445) * ( temp(i,k) - rei_min_temp ) 356 368 ! we clip the results 357 deimin = 20. 369 deimin = 20. + 40. * COS(ABS(latitude_deg(i)) / 180. * RPI) 358 370 deimax = 155. 359 371 dei = MIN(MAX(dei, deimin), deimax) … … 364 376 ! It is recommended to use the rei formula from Sun and Rikkus 1999 with a revision 365 377 ! from Sun 2001 (as in the IFS model) 366 iwc=icefrac_optics(i, k)*radocond(i, k)/pclc(i,k)*zrho(i,k)*1000. !in cloud ice water content in g/m3 378 ! in cloud ice water content in g/m3 379 IF ( ptconv(i,k) ) THEN 380 !--Needed because pclf does not contain convective clouds (should be fixed...) 381 iwc = icefrac_optics(i,k) * radocond(i,k) / pclc(i,k) * zrho(i,k) * 1000. 382 ELSE 383 iwc = icefrac_optics(i,k) * radocond(i,k) / pclf(i,k) * zrho(i,k) * 1000. 384 ENDIF 367 385 dei=(1.2351+0.0105*(temp(i,k)-273.15))*(45.8966*(iwc**0.2214) + & 368 386 & 0.7957*(iwc**0.2535)*(temp(i,k)-83.15)) … … 371 389 !Etienne: deimax and deimin controled by rei_max and rei_min in physiq.def 372 390 deimax=rei_max*2.0 373 deimin=2.0*rei_min+40 *cos(abs(latitude_deg(i))/180.*RPI)391 deimin=2.0*rei_min+40.*cos(abs(latitude_deg(i))/180.*RPI) 374 392 dei=min(dei,deimax) 375 393 dei=max(dei,deimin) … … 472 490 IF (iflag_rei .EQ. 2) THEN 473 491 ! in-cloud ice water content in g/m3 474 iwc = icefrac_optics(i,k) * radocond(i,k) / pclc(i,k) * zrho(i,k) * 1000. 492 IF ( ptconv(i,k) ) THEN 493 !--Needed because pclf does not contain convective clouds (should be fixed...) 494 iwc = icefrac_optics(i,k) * radocond(i,k) / pclc(i,k) * zrho(i,k) * 1000. 495 ELSE 496 iwc = icefrac_optics(i,k) * radocond(i,k) / pclf(i,k) * zrho(i,k) * 1000. 497 ENDIF 475 498 ! this formula is a simplified version of the Sun 2001 one (as in the IFS model, 476 499 ! and which is activated for iflag_rei = 1). … … 484 507 dei = rei_coef * (iwc**0.2445) * ( temp(i,k) - rei_min_temp ) 485 508 ! we clip the results 486 deimin = 20. 509 deimin = 20. + 40. * COS(ABS(latitude_deg(i)) / 180. * RPI) 487 510 deimax = 155. 488 511 dei = MIN(MAX(dei, deimin), deimax) … … 494 517 ! we use the rei formula from Sun and Rikkus 1999 with a revision 495 518 ! from Sun 2001 (as in the IFS model) 496 iwc=icefrac_optics(i, k)*radocond(i, k)/pclc(i,k)*zrho(i,k)*1000. !in cloud ice water content in g/m3 519 ! in cloud ice water content in g/m3 520 IF ( ptconv(i,k) ) THEN 521 !--Needed because pclf does not contain convective clouds (should be fixed...) 522 iwc = icefrac_optics(i,k) * radocond(i,k) / pclc(i,k) * zrho(i,k) * 1000. 523 ELSE 524 iwc = icefrac_optics(i,k) * radocond(i,k) / pclf(i,k) * zrho(i,k) * 1000. 525 ENDIF 497 526 dei=(1.2351+0.0105*(temp(i,k)-273.15))*(45.8966*(iwc**0.2214) + & 498 527 &0.7957*(iwc**0.2535)*(temp(i,k)-83.15)) … … 501 530 !Etienne: deimax and deimin controled by rei_max and rei_min in physiq.def 502 531 deimax=rei_max*2.0 503 deimin=2.0*rei_min+40 *cos(abs(latitude_deg(i))/180.*RPI)532 deimin=2.0*rei_min+40.*cos(abs(latitude_deg(i))/180.*RPI) 504 533 dei=min(dei,deimax) 505 534 dei=max(dei,deimin) … … 565 594 re(i, k) = rad_chaud(i, k)*fl(i, k) 566 595 rel = 0. 567 rei = 0.596 rei = 1. 568 597 pclc(i, k) = 0.0 569 598 pcltau(i, k) = 0.0 … … 571 600 572 601 !--AB contrails 573 rei ce_nocont(i,k) = 0.602 rei_cont = 1. 574 603 pclc_nocont(i,k) = 0. 575 604 pclc_cont(i,k) = 0. 576 pcltau_cont(i,k) = 0.577 pclemi_cont(i,k) = 0.578 pcltau_nocont(i,k) = 0.579 pclemi_nocont(i,k) = 0.605 pcltau_cont(i,k) = missing_val 606 pclemi_cont(i,k) = missing_val 607 pcltau_nocont(i,k) = missing_val 608 pclemi_nocont(i,k) = missing_val 580 609 581 610 ELSE … … 604 633 IF (iflag_rei .EQ. 2) THEN 605 634 ! in-cloud ice water content in g/m3 606 iwc = xfiwc_nocont(i, k) / pclc_nocont(i,k) * zrho(i,k) * 1000. 635 IF ( ptconv(i,k) ) THEN 636 !--Needed because pclf does not contain convective clouds (should be fixed...) 637 iwc = xfiwc_nocont(i, k) / pclc_nocont(i,k) * zrho(i,k) * 1000. 638 ELSE 639 iwc = xfiwc_nocont(i, k) / pclf_nocont(i,k) * zrho(i,k) * 1000. 640 ENDIF 607 641 ! this formula is a simplified version of the Sun 2001 one (as in the IFS model, 608 642 ! and which is activated for iflag_rei = 1). … … 616 650 dei = rei_coef * (iwc**0.2445) * ( temp(i,k) - rei_min_temp ) 617 651 ! we clip the results 618 !deimin = 20.652 deimin = 20. + 40. * COS(ABS(latitude_deg(i)) / 180. * RPI) 619 653 deimax = 155. 620 !dei = MIN(MAX(dei, deimin), deimax) 621 dei = MIN(dei, deimax) 654 dei = MIN(MAX(dei, deimin), deimax) 622 655 ! formula to convert effective diameter to effective radius 623 656 rei = 3. * SQRT(3.) / 8. * dei … … 628 661 ! we use the rei formula from Sun and Rikkus 1999 with a revision 629 662 ! from Sun 2001 (as in the IFS model) 630 iwc = xfiwc_nocont(i, k) / pclc_nocont(i,k) * zrho(i,k) * 1000. !in cloud ice water content in g/m3 663 ! in cloud ice water content in g/m3 664 IF ( ptconv(i,k) ) THEN 665 !--Needed because pclf does not contain convective clouds (should be fixed...) 666 iwc = xfiwc_nocont(i, k) / pclc_nocont(i,k) * zrho(i,k) * 1000. 667 ELSE 668 iwc = xfiwc_nocont(i, k) / pclf_nocont(i,k) * zrho(i,k) * 1000. 669 ENDIF 631 670 dei=(1.2351+0.0105*(temp(i,k)-273.15))*(45.8966*(iwc**0.2214) + & 632 671 &0.7957*(iwc**0.2535)*(temp(i,k)-83.15)) … … 635 674 !Etienne: deimax and deimin controled by rei_max and rei_min in physiq.def 636 675 deimax=rei_max*2.0 637 deimin=2.0*rei_min+40 *cos(abs(latitude_deg(i))/180.*RPI)676 deimin=2.0*rei_min+40.*cos(abs(latitude_deg(i))/180.*RPI) 638 677 dei=min(dei,deimax) 639 678 dei=max(dei,deimin) … … 662 701 !--Diagnostics of clouds emissivity, optical depth and ice crystal radius 663 702 !--without contrails 664 reice_nocont(i,k) = rei665 703 pcltau_nocont(i,k) = 3.0/2.0*(zflwp_var/rel) + zfiwp_var*(3.448E-03+2.431/rei) 666 704 ! -- cloud infrared emissivity: … … 678 716 !--Diagnostics of clouds emissivity, optical depth and ice crystal radius 679 717 !--without contrails 680 reice_nocont(i,k) = 1.681 718 pclc_nocont(i,k) = 0. 682 719 pclc(i,k) = contfra(i,k) … … 693 730 rei_cont = (mice_cont / rho_ice * 3. / 4. / RPI)**(1./3.) 694 731 !--Effective radius (in microns) 695 rei_cont = MIN(100., MAX(1., rei_cont / eff2vol_radius_contrails * 1e6)) 696 zfiwp_var = 1000.*(xfiwc(i, k)-xfiwc_nocont(i, k))& 697 / (pclc(i, k)-pclc_nocont(i, k))*rhodz(i, k) 732 rei_cont = MIN(100., MAX(rei_min_contrails, rei_cont / eff2vol_radius_contrails * 1e6)) 733 zfiwp_var = 1000.*xfiwc_cont(i, k)/pclc_cont(i, k)*rhodz(i, k) 698 734 699 735 pcltau_cont(i, k) = zfiwp_var*(3.448E-03+2.431/rei_cont) … … 711 747 ENDIF 712 748 713 rei = ( rei_cont * pclc_cont(i,k) + reice_nocont(i, k) * pclc_nocont(i, k) ) & 714 / ( pclc_cont(i,k) + pclc_nocont(i,k) ) 715 716 zflwp_var = 1000.*xflwc(i, k)/pclc(i, k)*rhodz(i, k) 717 zfiwp_var = 1000.*xfiwc(i, k)/pclc(i, k)*rhodz(i, k) 718 719 pcltau(i,k) = 3.0/2.0*(zflwp_var/rel) + zfiwp_var*(3.448E-03+2.431/rei) 720 ! -- cloud infrared emissivity: 721 ! [the broadband infrared absorption coefficient is PARAMETERized 722 ! as a function of the effective cld droplet radius] 723 ! Ebert and Curry (1992) formula as used by Kiehl & Zender (1995): 724 k_ice = k_ice0 + 1.0/rei 725 pclemi(i, k) = 1.0 - exp(-coef_chau*zflwp_var-df*k_ice*zfiwp_var) 749 pcltau(i,k) = (pclc_nocont(i,k)*pcltau_nocont(i,k) & 750 & + pclc_cont(i,k)*pcltau_cont(i,k)) & 751 & / pclc(i,k) 752 pclemi(i,k) = (pclc_nocont(i,k)*pclemi_nocont(i,k) & 753 & + pclc_cont(i,k)*pclemi_cont(i,k)) & 754 & / pclc(i,k) 755 756 IF ( pclc_nocont(i,k) .EQ. 0. ) THEN 757 pcltau_nocont(i,k) = missing_val 758 pclemi_nocont(i,k) = missing_val 759 ENDIF 760 761 IF ( pclc_cont(i,k) .EQ. 0. ) THEN 762 pcltau_cont(i,k) = missing_val 763 pclemi_cont(i,k) = missing_val 764 ENDIF 726 765 727 766 ENDIF 728 767 729 768 reice(i, k) = rei 769 reice_cont(i,k) = rei_cont 730 770 731 771 xflwp(i) = xflwp(i) + xflwc(i, k)*rhodz(i, k) 732 772 xfiwp(i) = xfiwp(i) + xfiwc(i, k)*rhodz(i, k) 733 xfiwp_nocont(i) = xfiwp_nocont(i) + xfiwc_nocont(i, k)*rhodz(i, k) 734 735 !--We weight the optical properties with the cloud fractions 736 !--This is only used for the diagnostics 737 pcltau_nocont(i, k) = pcltau_nocont(i, k) * pclc_nocont(i,k) 738 pclemi_nocont(i, k) = pclemi_nocont(i, k) * pclc_nocont(i,k) 739 pcltau_cont(i, k) = pcltau_cont(i, k) * pclc_cont(i,k) 740 pclemi_cont(i, k) = pclemi_cont(i, k) * pclc_cont(i,k) 741 pcltau(i, k) = pcltau(i, k) * pclc(i,k) 742 pclemi(i, k) = pclemi(i, k) * pclc(i,k) 773 xfiwp_cont(i) = xfiwp_cont(i) + xfiwc_cont(i, k)*rhodz(i, k) 743 774 744 775 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.
