Changeset 2609 for trunk/LMDZ.MARS/libf/phymars/albedocaps.F90
- Timestamp:
- Jan 18, 2022, 3:29:00 PM (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.MARS/libf/phymars/albedocaps.F90
r2508 r2609 9 9 emisice, albedice, watercaptag, albedo_h2o_cap, & 10 10 emissiv, albedodat 11 USE mod_phys_lmdz_transfert_para, ONLY: bcast 12 USE mod_phys_lmdz_para, ONLY: is_master 11 13 implicit none 12 14 … … 23 25 ! local variables: 24 26 logical,save :: firstcall=.true. 27 !$OMP THREADPRIVATE(firstcall) 25 28 integer :: ig,icap 29 30 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 31 32 ! local variables: 33 real,save :: zls_old ! value of zls from a previous call 34 real,save :: pi,radeg ! to convert radians to degrees 35 36 !$OMP THREADPRIVATE(zls_old,pi,radeg) 37 38 ! TES datasets: (hard coded fixed length/sizes; for now) 39 integer,parameter :: TESlonsize=72 40 ! longitudes, in TES files, in degrees, from TESlon(1)=-177.5 to TESlon(72)=177.5 41 real,save :: TESlon(TESlonsize) 42 integer,parameter :: TESlatsize=30 43 ! latitudes (north hemisphere file), in degrees, from TESlatn(1)=31, 44 ! to TESlatn(30)=89 ; TESlatn(8)=45 45 real,save :: TESlatn(TESlatsize) 46 ! latitudes (south hemisphere file), in degrees, from TESlats(1)=-89, 47 ! to TESlats(30)=-31 ; TESlats(23)=-45 48 real,save :: TESlats(TESlatsize) 49 integer,parameter :: TESlssize=72 50 ! Solar longitude in TES files, TESls(1)=2.5 to TESls(72)=357.5 51 real,save :: TESls(TESlssize) 52 ! TES North albedo (=-1 for missing values) 53 real,save :: TESalbn(TESlonsize,TESlatsize,TESlssize) 54 ! TES South albedo (=-1 for missing values) 55 real,save :: TESalbs(TESlonsize,TESlatsize,TESlssize) 56 57 !$OMP THREADPRIVATE(TESlon,TESlatn,TESlats,TESls,TESalbn,TESalbs) 58 59 60 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 26 61 27 62 ! 1. Initializations … … 47 82 endif 48 83 84 call read_TES_icecap_albedo(zls_old,pi,radeg,TESlon,TESlatn,TESlats,TESls,TESalbn,TESalbs) 85 49 86 firstcall=.false. 50 87 endif ! of if (firstcall) 51 88 52 89 do ig=1,ngrid 53 90 if (latitude(ig).lt.0.) then … … 56 93 icap=1 ! Northern hemisphere 57 94 endif 58 95 59 96 if (piceco2(ig).gt.0) then 60 97 ! set emissivity of surface to be the ice emissivity … … 62 99 ! set the surface albedo to be the ice albedo 63 100 if (TESicealbedo) then 64 call TES_icecap_albedo(zls,ig,psolaralb(ig,1),icap )101 call TES_icecap_albedo(zls,ig,psolaralb(ig,1),icap,zls_old,pi,radeg,TESlon,TESlatn,TESlats,TESls,TESalbn,TESalbs) 65 102 psolaralb(ig,2)=psolaralb(ig,1) 66 103 else … … 82 119 endif ! of if (piceco2(ig).gt.0) 83 120 enddo ! of ig=1,ngrid 121 84 122 end subroutine albedocaps 85 123 86 124 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 87 subroutine TES_icecap_albedo(zls,ig,alb,icap) 88 89 use geometry_mod, only: latitude, longitude ! in radians 90 use surfdat_h, only: albedice, TESice_Ncoef, TESice_Scoef 125 126 subroutine read_TES_icecap_albedo(zls_old,pi,radeg,TESlon,TESlatn,TESlats,TESls,TESalbn,TESalbs) 127 91 128 use datafile_mod, only: datadir 92 129 use netcdf, only: nf90_open, NF90_NOWRITE, NF90_NOERR, & 93 130 nf90_strerror, nf90_inq_varid, nf90_get_var, nf90_close 131 USE mod_phys_lmdz_para, ONLY: is_master 132 USE mod_phys_lmdz_transfert_para, ONLY: bcast 94 133 95 134 implicit none 96 135 97 136 ! arguments: 98 real,intent(in) :: zls ! solar longitude (rad)99 integer,intent(in) :: ig ! grid point index100 real,intent(out) :: alb ! (interpolated) TES ice albedo at that grid point101 integer :: icap ! =1: Northern hemisphere =2: Southern hemisphere102 137 103 138 ! local variables: 104 logical,save :: firstcall=.true. 105 real,save :: zls_old ! value of zls from a previous call 106 integer,save :: tinf,tsup ! encompassing time indexes of TES data 107 real,save :: reltime ! relative position in-between time indexes (in [0;1]) 108 integer :: latinf,latsup ! encompassing latitude indexes of TES data 109 real :: rellat ! relative position in-between latitude indexes (in[0;1]) 110 integer :: loninf,lonsup ! encompassing longitude indexes of TES data 111 real :: rellon !relative position in-between longitude indexes (in[0;1]) 112 real,save :: pi,radeg ! to convert radians to degrees 113 real :: zlsd ! solar longitude, in degrees 114 real :: latd ! latitude, in degrees 115 real :: lond ! longitude, in degrees 116 integer :: i 139 real:: zls_old ! value of zls from a previous call 140 real:: pi,radeg ! to convert radians to degrees 117 141 character(len=20),parameter :: modname="TES_icecap_albedo" 118 142 119 143 ! TES datasets: (hard coded fixed length/sizes; for now) 120 integer,parameter :: TESlonsize=72 121 real,parameter :: TESdeltalon=5.0 ! step in longitude in TES files 144 integer,parameter:: TESlonsize=72 122 145 ! longitudes, in TES files, in degrees, from TESlon(1)=-177.5 to TESlon(72)=177.5 123 real,save :: TESlon(TESlonsize) 124 integer,parameter :: TESlatsize=30 125 real,parameter :: TESdeltalat=2.0 ! step in latitude in TES files 146 real:: TESlon(TESlonsize) 147 integer,parameter:: TESlatsize=30 126 148 ! latitudes (north hemisphere file), in degrees, from TESlatn(1)=31, 127 149 ! to TESlatn(30)=89 ; TESlatn(8)=45 128 real,parameter :: TESlatnmin=45. ! minimum TES latitude (North hemisphere) 129 real,parameter :: TESlatsmax=-45. ! maximum TES latitude (South hemisphere) 130 real,save :: TESlatn(TESlatsize) 150 real:: TESlatn(TESlatsize) 131 151 ! latitudes (south hemisphere file), in degrees, from TESlats(1)=-89, 132 152 ! to TESlats(30)=-31 ; TESlats(23)=-45 133 real,save :: TESlats(TESlatsize) 134 integer,parameter :: TESlssize=72 135 real,parameter :: TESdeltals=5.0 ! step in solar longitude in TES files 153 real:: TESlats(TESlatsize) 154 integer,parameter:: TESlssize=72 136 155 ! Solar longitude in TES files, TESls(1)=2.5 to TESls(72)=357.5 137 real ,save:: TESls(TESlssize)156 real:: TESls(TESlssize) 138 157 ! TES North albedo (=-1 for missing values) 139 real ,save:: TESalbn(TESlonsize,TESlatsize,TESlssize)158 real:: TESalbn(TESlonsize,TESlatsize,TESlssize) 140 159 ! TES South albedo (=-1 for missing values) 141 real,save :: TESalbs(TESlonsize,TESlatsize,TESlssize) 142 ! encompassing nodes arranged as follow : 4 3 143 real :: val(4) ! 1 2 160 real:: TESalbs(TESlonsize,TESlatsize,TESlssize) 144 161 145 162 !NetCDF variables: … … 149 166 150 167 ! 0. Preliminary stuff 151 if (firstcall) then 168 169 if(is_master) then 170 152 171 ! Load TES albedoes for Northern Hemisphere 153 172 ierr=nf90_open(trim(datadir)//"/npsc_albedo.nc",NF90_NOWRITE,nid) … … 273 292 274 293 zls_old=-999 ! dummy initialization 275 276 firstcall=.false. 277 endif ! of if firstcall 294 295 endif !is_master 296 297 call bcast(TESlon) 298 call bcast(TESlatn) 299 call bcast(TESlats) 300 call bcast(TESls) 301 call bcast(TESalbn) 302 call bcast(TESalbs) 303 call bcast(pi) 304 call bcast(zls_old) 305 call bcast(radeg) 306 307 end subroutine read_TES_icecap_albedo 308 309 subroutine TES_icecap_albedo(zls,ig,alb,icap,zls_old,pi,radeg,TESlon,TESlatn,TESlats,TESls,TESalbn,TESalbs) 310 311 use geometry_mod, only: latitude, longitude ! in radians 312 use surfdat_h, only: albedice, TESice_Ncoef, TESice_Scoef 313 use datafile_mod, only: datadir 314 use netcdf, only: nf90_open, NF90_NOWRITE, NF90_NOERR, & 315 nf90_strerror, nf90_inq_varid, nf90_get_var, nf90_close 316 USE mod_phys_lmdz_transfert_para, ONLY: bcast 317 318 implicit none 319 320 ! arguments: 321 real,intent(in) :: zls ! solar longitude (rad) 322 integer,intent(in) :: ig ! grid point index 323 real,intent(out) :: alb ! (interpolated) TES ice albedo at that grid point 324 integer :: icap ! =1: Northern hemisphere =2: Southern hemisphere 325 326 ! local variables: 327 integer,save :: tinf,tsup ! encompassing time indexes of TES data 328 real,save :: reltime ! relative position in-between time indexes (in [0;1]) 329 integer :: latinf,latsup ! encompassing latitude indexes of TES data 330 real :: rellat ! relative position in-between latitude indexes (in[0;1]) 331 integer :: loninf,lonsup ! encompassing longitude indexes of TES data 332 real :: rellon !relative position in-between longitude indexes (in[0;1]) 333 real :: zlsd ! solar longitude, in degrees 334 real :: latd ! latitude, in degrees 335 real :: lond ! longitude, in degrees 336 integer :: i 337 338 !$OMP THREADPRIVATE(tinf,tsup,reltime) 339 340 ! TES datasets: (hard coded fixed length/sizes; for now) 341 real,parameter :: TESdeltalon=5.0 ! step in longitude in TES files 342 real,parameter :: TESdeltalat=2.0 ! step in latitude in TES files 343 real,parameter :: TESlatnmin=45. ! minimum TES latitude (North hemisphere) 344 real,parameter :: TESlatsmax=-45. ! maximum TES latitude (South hemisphere) 345 ! to TESlats(30)=-31 ; TESlats(23)=-45 346 real,parameter :: TESdeltals=5.0 ! step in solar longitude in TES files 347 348 real:: zls_old ! value of zls from a previous call 349 real:: pi,radeg ! to convert radians to degrees 350 351 ! TES datasets: (hard coded fixed length/sizes; for now) 352 integer,parameter:: TESlonsize=72 353 ! longitudes, in TES files, in degrees, from TESlon(1)=-177.5 to TESlon(72)=177.5 354 real:: TESlon(TESlonsize) 355 integer,parameter:: TESlatsize=30 356 ! latitudes (north hemisphere file), in degrees, from TESlatn(1)=31, 357 ! to TESlatn(30)=89 ; TESlatn(8)=45 358 real:: TESlatn(TESlatsize) 359 ! latitudes (south hemisphere file), in degrees, from TESlats(1)=-89, 360 ! to TESlats(30)=-31 ; TESlats(23)=-45 361 real:: TESlats(TESlatsize) 362 integer,parameter:: TESlssize=72 363 ! Solar longitude in TES files, TESls(1)=2.5 to TESls(72)=357.5 364 real:: TESls(TESlssize) 365 ! TES North albedo (=-1 for missing values) 366 real:: TESalbn(TESlonsize,TESlatsize,TESlssize) 367 ! TES South albedo (=-1 for missing values) 368 real:: TESalbs(TESlonsize,TESlatsize,TESlssize) 369 ! encompassing nodes arranged as follow : 4 3 370 real :: val(4) ! 1 2 278 371 279 372 ! 1. Identify encompassing latitudes … … 281 374 ! Check that latitude is such that there is TES data to use 282 375 ! (ie: latitude 45 deg and poleward) otherwise use 'default' albedoes 376 283 377 latd=latitude(ig)*radeg ! latitude, in degrees 284 378 if (icap.eq.1) then
Note: See TracChangeset
for help on using the changeset viewer.