| [3435] | 1 | MODULE create_limit_unstruct_mod |
|---|
| 2 | PRIVATE |
|---|
| 3 | INTEGER, PARAMETER :: lmdep=12 |
|---|
| 4 | |
|---|
| 5 | PUBLIC create_limit_unstruct |
|---|
| 6 | |
|---|
| 7 | CONTAINS |
|---|
| 8 | |
|---|
| [3471] | 9 | |
|---|
| [3435] | 10 | SUBROUTINE create_limit_unstruct |
|---|
| 11 | USE dimphy |
|---|
| [4619] | 12 | USE lmdz_xios |
|---|
| [3435] | 13 | USE ioipsl, ONLY : ioget_year_len |
|---|
| 14 | USE time_phylmdz_mod, ONLY : annee_ref |
|---|
| 15 | USE indice_sol_mod |
|---|
| 16 | USE phys_state_var_mod |
|---|
| 17 | USE mod_phys_lmdz_para |
|---|
| 18 | IMPLICIT NONE |
|---|
| 19 | INCLUDE "iniprint.h" |
|---|
| [3471] | 20 | REAL, DIMENSION(:,:),ALLOCATABLE :: sic |
|---|
| 21 | REAL, DIMENSION(:,:),ALLOCATABLE :: sst |
|---|
| [3435] | 22 | REAL, DIMENSION(klon,lmdep) :: rugos |
|---|
| 23 | REAL, DIMENSION(klon,lmdep) :: albedo |
|---|
| [3471] | 24 | REAL, DIMENSION(:,:),ALLOCATABLE :: sic_mpi |
|---|
| 25 | REAL, DIMENSION(:,:),ALLOCATABLE :: sst_mpi |
|---|
| [3435] | 26 | REAL, DIMENSION(klon_mpi,lmdep) :: rugos_mpi |
|---|
| 27 | REAL, DIMENSION(klon_mpi,lmdep) :: albedo_mpi |
|---|
| 28 | INTEGER :: ndays |
|---|
| 29 | REAL :: fi_ice(klon) |
|---|
| 30 | REAL, ALLOCATABLE :: sic_year(:,:) |
|---|
| 31 | REAL, ALLOCATABLE :: sst_year(:,:) |
|---|
| 32 | REAL, ALLOCATABLE :: rugos_year(:,:) |
|---|
| 33 | REAL, ALLOCATABLE :: albedo_year(:,:) |
|---|
| 34 | REAL, ALLOCATABLE :: pctsrf_t(:,:,:) |
|---|
| 35 | REAL, ALLOCATABLE :: phy_bil(:,:) |
|---|
| 36 | REAL, ALLOCATABLE :: sst_year_mpi(:,:) |
|---|
| 37 | REAL, ALLOCATABLE :: rugos_year_mpi(:,:) |
|---|
| 38 | REAL, ALLOCATABLE :: albedo_year_mpi(:,:) |
|---|
| 39 | REAL, ALLOCATABLE :: pctsrf_t_mpi(:,:,:) |
|---|
| 40 | REAL, ALLOCATABLE :: phy_bil_mpi(:,:) |
|---|
| 41 | INTEGER :: l,k |
|---|
| 42 | INTEGER :: nbad |
|---|
| [3471] | 43 | INTEGER :: sic_time_axis_size |
|---|
| 44 | INTEGER :: sst_time_axis_size |
|---|
| 45 | CHARACTER(LEN=99) :: mess ! error message |
|---|
| [3435] | 46 | |
|---|
| [3471] | 47 | |
|---|
| [3435] | 48 | ndays=ioget_year_len(annee_ref) |
|---|
| 49 | |
|---|
| [3471] | 50 | IF (is_omp_master) CALL xios_get_axis_attr("time_sic",n_glo=sic_time_axis_size) |
|---|
| 51 | CALL bcast_omp(sic_time_axis_size) |
|---|
| 52 | ALLOCATE(sic_mpi(klon_mpi,sic_time_axis_size)) |
|---|
| 53 | ALLOCATE(sic(klon,sic_time_axis_size)) |
|---|
| 54 | |
|---|
| 55 | |
|---|
| 56 | IF (is_omp_master) CALL xios_get_axis_attr("time_sst",n_glo=sst_time_axis_size) |
|---|
| 57 | CALL bcast_omp(sst_time_axis_size) |
|---|
| 58 | ALLOCATE(sst_mpi(klon_mpi,sst_time_axis_size)) |
|---|
| 59 | ALLOCATE(sst(klon,sst_time_axis_size)) |
|---|
| 60 | |
|---|
| [3435] | 61 | IF (is_omp_master) THEN |
|---|
| 62 | CALL xios_recv_field("sic_limit",sic_mpi) |
|---|
| 63 | CALL xios_recv_field("sst_limit",sst_mpi) |
|---|
| 64 | CALL xios_recv_field("rugos_limit",rugos_mpi) |
|---|
| 65 | CALL xios_recv_field("albedo_limit",albedo_mpi) |
|---|
| 66 | ENDIF |
|---|
| 67 | CALL scatter_omp(sic_mpi,sic) |
|---|
| 68 | CALL scatter_omp(sst_mpi,sst) |
|---|
| 69 | CALL scatter_omp(rugos_mpi,rugos) |
|---|
| 70 | CALL scatter_omp(albedo_mpi,albedo) |
|---|
| 71 | |
|---|
| 72 | ALLOCATE(sic_year(klon,ndays)) |
|---|
| 73 | ALLOCATE(sst_year(klon,ndays)) |
|---|
| 74 | ALLOCATE(rugos_year(klon,ndays)) |
|---|
| 75 | ALLOCATE(albedo_year(klon,ndays)) |
|---|
| 76 | ALLOCATE(pctsrf_t(klon,nbsrf,ndays)) |
|---|
| 77 | ALLOCATE(phy_bil(klon,ndays)); phy_bil=0.0 |
|---|
| 78 | |
|---|
| [3471] | 79 | |
|---|
| [3435] | 80 | ! sic |
|---|
| [3471] | 81 | IF (sic_time_axis_size==lmdep) THEN |
|---|
| 82 | CALL time_interpolation(ndays,sic,'gregorian',sic_year) |
|---|
| 83 | ELSE IF (sic_time_axis_size==ndays) THEN |
|---|
| 84 | sic_year=sic |
|---|
| 85 | ELSE |
|---|
| 86 | WRITE(mess,*) 'sic time axis is nor montly, nor daily. sic time interpolation ',& |
|---|
| 87 | 'is requiered but is not currently managed' |
|---|
| 88 | CALL abort_physic('create_limit_unstruct',TRIM(mess),1) |
|---|
| 89 | ENDIF |
|---|
| 90 | |
|---|
| [3435] | 91 | sic_year(:,:)=sic_year(:,:)/100. ! convert percent to fraction |
|---|
| 92 | WHERE(sic_year(:,:)>1.0) sic_year(:,:)=1.0 ! Some fractions have some time large negative values |
|---|
| 93 | WHERE(sic_year(:,:)<0.0) sic_year(:,:)=0.0 ! probably better to apply alse this filter before horizontal interpolation |
|---|
| 94 | |
|---|
| 95 | ! sst |
|---|
| [3471] | 96 | IF (sst_time_axis_size==lmdep) THEN |
|---|
| 97 | CALL time_interpolation(ndays,sst,'gregorian',sst_year) |
|---|
| 98 | ELSE IF (sst_time_axis_size==ndays) THEN |
|---|
| 99 | sst_year=sst |
|---|
| 100 | ELSE |
|---|
| 101 | WRITE(mess,*)'sic time axis is nor montly, nor daily. sic time interpolation ',& |
|---|
| 102 | 'is requiered but is not currently managed' |
|---|
| 103 | CALL abort_physic('create_limit_unstruct',TRIM(mess),1) |
|---|
| 104 | ENDIF |
|---|
| [3435] | 105 | WHERE(sst_year(:,:)<271.38) sst_year(:,:)=271.38 |
|---|
| 106 | |
|---|
| [3471] | 107 | |
|---|
| [3435] | 108 | ! rugos |
|---|
| 109 | DO l=1, lmdep |
|---|
| 110 | WHERE(NINT(zmasq(:))/=1) rugos(:,l)=0.001 |
|---|
| 111 | ENDDO |
|---|
| [4361] | 112 | CALL time_interpolation(ndays,rugos,'360_day',rugos_year) |
|---|
| [3435] | 113 | |
|---|
| 114 | ! albedo |
|---|
| [4361] | 115 | CALL time_interpolation(ndays,albedo,'360_day',albedo_year) |
|---|
| [3435] | 116 | |
|---|
| 117 | |
|---|
| 118 | DO k=1,ndays |
|---|
| 119 | fi_ice=sic_year(:,k) |
|---|
| 120 | WHERE(fi_ice>=1.0 ) fi_ice=1.0 |
|---|
| 121 | WHERE(fi_ice<EPSFRA) fi_ice=0.0 |
|---|
| 122 | pctsrf_t(:,is_ter,k)=pctsrf(:,is_ter) ! land soil |
|---|
| 123 | pctsrf_t(:,is_lic,k)=pctsrf(:,is_lic) ! land ice |
|---|
| 124 | |
|---|
| 125 | !! IF (icefile==trim(fcpldsic)) THEN ! SIC=pICE*(1-LIC-TER) |
|---|
| 126 | !! pctsrf_t(:,is_sic,k)=fi_ice(:)*(1.-pctsrf(:,is_lic)-pctsrf(:,is_ter)) |
|---|
| 127 | !! ELSE IF (icefile==trim(fhistsic)) THEN ! SIC=pICE |
|---|
| 128 | !! pctsrf_t(:,is_sic,k)=fi_ice(:) |
|---|
| 129 | !! ELSE ! icefile==famipsic ! SIC=pICE-LIC |
|---|
| 130 | pctsrf_t(:,is_sic,k)=fi_ice-pctsrf_t(:,is_lic,k) |
|---|
| 131 | ! END IF |
|---|
| 132 | WHERE(pctsrf_t(:,is_sic,k)<=0) pctsrf_t(:,is_sic,k)=0. |
|---|
| 133 | WHERE(1.0-zmasq<EPSFRA) |
|---|
| 134 | pctsrf_t(:,is_sic,k)=0.0 |
|---|
| 135 | pctsrf_t(:,is_oce,k)=0.0 |
|---|
| 136 | ELSEWHERE |
|---|
| 137 | WHERE(pctsrf_t(:,is_sic,k)>=1.0-zmasq) |
|---|
| 138 | pctsrf_t(:,is_sic,k)=1.0-zmasq |
|---|
| 139 | pctsrf_t(:,is_oce,k)=0.0 |
|---|
| 140 | ELSEWHERE |
|---|
| 141 | pctsrf_t(:,is_oce,k)=1.0-zmasq-pctsrf_t(:,is_sic,k) |
|---|
| 142 | WHERE(pctsrf_t(:,is_oce,k)<EPSFRA) |
|---|
| 143 | pctsrf_t(:,is_oce,k)=0.0 |
|---|
| 144 | pctsrf_t(:,is_sic,k)=1.0-zmasq |
|---|
| 145 | END WHERE |
|---|
| 146 | END WHERE |
|---|
| 147 | END WHERE |
|---|
| 148 | nbad=COUNT(pctsrf_t(:,is_oce,k)<0.0) |
|---|
| 149 | IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb point = ',nbad |
|---|
| 150 | nbad=COUNT(abs(sum(pctsrf_t(:,:,k),dim=2)-1.0)>EPSFRA) |
|---|
| 151 | IF(nbad>0) WRITE(lunout,*) 'pb sous surface pour nb points = ',nbad |
|---|
| 152 | END DO |
|---|
| 153 | |
|---|
| 154 | ALLOCATE(sst_year_mpi(klon_mpi,ndays)) |
|---|
| 155 | ALLOCATE(rugos_year_mpi(klon_mpi,ndays)) |
|---|
| 156 | ALLOCATE(albedo_year_mpi(klon_mpi,ndays)) |
|---|
| 157 | ALLOCATE(pctsrf_t_mpi(klon_mpi,nbsrf,ndays)) |
|---|
| 158 | ALLOCATE(phy_bil_mpi(klon_mpi,ndays)) |
|---|
| 159 | |
|---|
| 160 | CALL gather_omp(pctsrf_t , pctsrf_t_mpi) |
|---|
| 161 | CALL gather_omp(sst_year , sst_year_mpi) |
|---|
| 162 | CALL gather_omp(phy_bil , phy_bil_mpi) |
|---|
| 163 | CALL gather_omp(albedo_year, albedo_year_mpi) |
|---|
| 164 | CALL gather_omp(rugos_year , rugos_year_mpi) |
|---|
| 165 | |
|---|
| 166 | IF (is_omp_master) THEN |
|---|
| 167 | CALL xios_send_field("foce_limout",pctsrf_t_mpi(:,is_oce,:)) |
|---|
| 168 | CALL xios_send_field("fsic_limout",pctsrf_t_mpi(:,is_sic,:)) |
|---|
| 169 | CALL xios_send_field("fter_limout",pctsrf_t_mpi(:,is_ter,:)) |
|---|
| 170 | CALL xios_send_field("flic_limout",pctsrf_t_mpi(:,is_lic,:)) |
|---|
| 171 | CALL xios_send_field("sst_limout", sst_year_mpi) |
|---|
| 172 | CALL xios_send_field("bils_limout",phy_bil_mpi) |
|---|
| 173 | CALL xios_send_field("alb_limout", albedo_year_mpi) |
|---|
| 174 | CALL xios_send_field("rug_limout", rugos_year_mpi) |
|---|
| 175 | ENDIF |
|---|
| 176 | END SUBROUTINE create_limit_unstruct |
|---|
| 177 | |
|---|
| 178 | |
|---|
| 179 | SUBROUTINE time_interpolation(ndays,field_in,calendar,field_out) |
|---|
| 180 | USE pchsp_95_m, only: pchsp_95 |
|---|
| 181 | USE pchfe_95_m, only: pchfe_95 |
|---|
| 182 | USE arth_m, only: arth |
|---|
| 183 | USE dimphy, ONLY : klon |
|---|
| 184 | USE ioipsl, ONLY : ioget_year_len |
|---|
| 185 | USE time_phylmdz_mod, ONLY : annee_ref |
|---|
| 186 | USE mod_phys_lmdz_para |
|---|
| 187 | IMPLICIT NONE |
|---|
| 188 | INCLUDE "iniprint.h" |
|---|
| 189 | |
|---|
| 190 | INTEGER, INTENT(IN) :: ndays |
|---|
| 191 | REAL, INTENT(IN) :: field_in(klon,lmdep) |
|---|
| 192 | CHARACTER(LEN=*),INTENT(IN) :: calendar |
|---|
| 193 | REAL, INTENT(OUT) :: field_out(klon,ndays) |
|---|
| 194 | |
|---|
| 195 | INTEGER :: ndays_in |
|---|
| 196 | REAL :: timeyear(lmdep) |
|---|
| 197 | REAL :: yder(lmdep) |
|---|
| 198 | INTEGER :: ij,ierr, n_extrap |
|---|
| 199 | LOGICAL :: skip |
|---|
| [3531] | 200 | |
|---|
| 201 | CHARACTER (len = 50) :: modname = 'create_limit_unstruct.time_interpolation' |
|---|
| 202 | CHARACTER (len = 80) :: abort_message |
|---|
| 203 | |
|---|
| [3435] | 204 | |
|---|
| 205 | IF (is_omp_master) ndays_in=year_len(annee_ref, calendar) |
|---|
| 206 | CALL bcast_omp(ndays_in) |
|---|
| [3469] | 207 | IF (is_omp_master) timeyear=mid_months(annee_ref, calendar, lmdep) |
|---|
| 208 | CALL bcast_omp(timeyear) |
|---|
| [3435] | 209 | |
|---|
| 210 | n_extrap = 0 |
|---|
| 211 | skip=.FALSE. |
|---|
| 212 | DO ij=1,klon |
|---|
| 213 | yder = pchsp_95(timeyear, field_in(ij, :), ibeg=2, iend=2, vc_beg=0., vc_end=0.) |
|---|
| 214 | CALL pchfe_95(timeyear, field_in(ij, :), yder, skip, arth(0., real(ndays_in) / ndays, ndays), field_out(ij, :), ierr) |
|---|
| [3531] | 215 | if (ierr < 0) then |
|---|
| 216 | abort_message='error in pchfe_95' |
|---|
| 217 | CALL abort_physic(modname,abort_message,1) |
|---|
| 218 | endif |
|---|
| [3435] | 219 | n_extrap = n_extrap + ierr |
|---|
| 220 | END DO |
|---|
| 221 | |
|---|
| 222 | IF (n_extrap /= 0) then |
|---|
| 223 | WRITE(lunout,*) "get_2Dfield pchfe_95: n_extrap = ", n_extrap |
|---|
| 224 | ENDIF |
|---|
| 225 | |
|---|
| 226 | |
|---|
| 227 | END SUBROUTINE time_interpolation |
|---|
| 228 | !------------------------------------------------------------------------------- |
|---|
| 229 | ! |
|---|
| 230 | FUNCTION year_len(y,cal_in) |
|---|
| 231 | ! |
|---|
| 232 | !------------------------------------------------------------------------------- |
|---|
| 233 | USE ioipsl, ONLY : ioget_calendar,ioconf_calendar,lock_calendar,ioget_year_len |
|---|
| 234 | IMPLICIT NONE |
|---|
| 235 | !------------------------------------------------------------------------------- |
|---|
| 236 | ! Arguments: |
|---|
| 237 | INTEGER :: year_len |
|---|
| 238 | INTEGER, INTENT(IN) :: y |
|---|
| 239 | CHARACTER(LEN=*), INTENT(IN) :: cal_in |
|---|
| 240 | !------------------------------------------------------------------------------- |
|---|
| 241 | ! Local variables: |
|---|
| 242 | CHARACTER(LEN=20) :: cal_out ! calendar (for outputs) |
|---|
| 243 | !------------------------------------------------------------------------------- |
|---|
| 244 | !--- Getting the input calendar to reset at the end of the function |
|---|
| 245 | CALL ioget_calendar(cal_out) |
|---|
| 246 | |
|---|
| 247 | !--- Unlocking calendar and setting it to wanted one |
|---|
| 248 | CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in)) |
|---|
| 249 | |
|---|
| 250 | !--- Getting the number of days in this year |
|---|
| 251 | year_len=ioget_year_len(y) |
|---|
| 252 | |
|---|
| 253 | !--- Back to original calendar |
|---|
| 254 | CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out)) |
|---|
| 255 | |
|---|
| 256 | END FUNCTION year_len |
|---|
| 257 | ! |
|---|
| 258 | !------------------------------------------------------------------------------- |
|---|
| 259 | |
|---|
| 260 | |
|---|
| 261 | !------------------------------------------------------------------------------- |
|---|
| 262 | ! |
|---|
| 263 | FUNCTION mid_months(y,cal_in,nm) |
|---|
| 264 | ! |
|---|
| 265 | !------------------------------------------------------------------------------- |
|---|
| 266 | USE ioipsl, ONLY : ioget_calendar,ioconf_calendar,lock_calendar,ioget_mon_len |
|---|
| 267 | IMPLICIT NONE |
|---|
| 268 | !------------------------------------------------------------------------------- |
|---|
| 269 | ! Arguments: |
|---|
| 270 | INTEGER, INTENT(IN) :: y ! year |
|---|
| 271 | CHARACTER(LEN=*), INTENT(IN) :: cal_in ! calendar |
|---|
| 272 | INTEGER, INTENT(IN) :: nm ! months/year number |
|---|
| 273 | REAL, DIMENSION(nm) :: mid_months ! mid-month times |
|---|
| 274 | !------------------------------------------------------------------------------- |
|---|
| 275 | ! Local variables: |
|---|
| 276 | CHARACTER(LEN=99) :: mess ! error message |
|---|
| 277 | CHARACTER(LEN=20) :: cal_out ! calendar (for outputs) |
|---|
| 278 | INTEGER, DIMENSION(nm) :: mnth ! months lengths (days) |
|---|
| 279 | INTEGER :: m ! months counter |
|---|
| 280 | INTEGER :: nd ! number of days |
|---|
| 281 | INTEGER :: k |
|---|
| 282 | !------------------------------------------------------------------------------- |
|---|
| 283 | nd=year_len(y,cal_in) |
|---|
| 284 | |
|---|
| 285 | IF(nm==12) THEN |
|---|
| 286 | |
|---|
| 287 | !--- Getting the input calendar to reset at the end of the function |
|---|
| 288 | CALL ioget_calendar(cal_out) |
|---|
| 289 | |
|---|
| 290 | !--- Unlocking calendar and setting it to wanted one |
|---|
| 291 | CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in)) |
|---|
| 292 | |
|---|
| 293 | !--- Getting the length of each month |
|---|
| 294 | DO m=1,nm; mnth(m)=ioget_mon_len(y,m); END DO |
|---|
| 295 | |
|---|
| 296 | !--- Back to original calendar |
|---|
| 297 | CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out)) |
|---|
| 298 | |
|---|
| 299 | ELSE IF(MODULO(nd,nm)/=0) THEN |
|---|
| 300 | WRITE(mess,'(a,i3,a,i3,a)')'Unconsistent calendar: ',nd,' days/year, but ',& |
|---|
| 301 | nm,' months/year. Months number should divide days number.' |
|---|
| 302 | CALL abort_physic('mid_months',TRIM(mess),1) |
|---|
| 303 | |
|---|
| 304 | ELSE |
|---|
| 305 | mnth=(/(m,m=1,nm,nd/nm)/) |
|---|
| 306 | END IF |
|---|
| 307 | |
|---|
| 308 | !--- Mid-months times |
|---|
| 309 | mid_months(1)=0.5*REAL(mnth(1)) |
|---|
| 310 | DO k=2,nm |
|---|
| 311 | mid_months(k)=mid_months(k-1)+0.5*REAL(mnth(k-1)+mnth(k)) |
|---|
| 312 | END DO |
|---|
| 313 | |
|---|
| 314 | END FUNCTION mid_months |
|---|
| 315 | |
|---|
| 316 | |
|---|
| 317 | END MODULE create_limit_unstruct_mod |
|---|
