Oasis3 4.0.2
|
00001 SUBROUTINE fiasco (poutp, ptlon, ptlat, ptsurf, 00002 $ ktmsk, ktlon, ktlat, cdtper, ktper, 00003 $ pinpt, pslon, pslat, pssurf, 00004 $ ksmsk, kslon, kslat, cdsper, ksper, 00005 $ plonz, platz, psgrb, psgra, psfrb, psfra, 00006 $ kfield, cdint, cdtyp, cddim) 00007 C**** 00008 C ***************************** 00009 C * OASIS ROUTINE - LEVEL 3 * 00010 C * ------------- ------- * 00011 C ***************************** 00012 C 00013 C**** *fiasco* - Interpolation driver 00014 C 00015 C Purpose: 00016 C ------- 00017 C Interpolation from source grid to target grid 00018 C Use Fast Scalar Interpolator (FSCINT) package from Yves Chartier 00019 C or ANAIS software from Olivier Thual et al. 00020 C 00021 C** Interface: 00022 C --------- 00023 C *CALL* *fiasco* (poutp, ptlon, ptlat, ptsurf, 00024 C ktmsk, ktlon, ktlat, cdtper, ktper, 00025 C pinpt, pslon, pslat, pssurf, 00026 C ksmsk, kslon, kslat, cdsper, ksper, 00027 C plonz, platz, psgrb, psgra, psfrb, psfra, 00028 C kfield, cdint, cdtyp, cddim)* 00029 C 00030 C Input: 00031 C ----- 00032 C ptlon : longitudes of target grid (real 2D) 00033 C ptlat : latitudes of target grid (real 2D) 00034 C ptsurf : target grid square surfaces (real 2D) 00035 C ktmsk : mask of target grid (integer 2D) 00036 C ktlon : number of longitudes for target grid 00037 C ktlat : number of latitudes for target grid 00038 C cdtper : target grid periodicity 00039 C ktper : number of overlapped points for target grid 00040 C pinpt : input field on source grid (real 2D) 00041 C pslon : longitudes of source grid (real 2D) 00042 C pslat : latitudes of source grid (real 2D) 00043 C pssurf : source grid square surfaces (real 2D) 00044 C ksmsk : mask of source grid (integer 2D) 00045 C kslon : number of longitudes for source grid 00046 C kslat : number of latitudes for source grid 00047 C cdsper : source grid periodicity 00048 C ksper : number of overlapped points for source grid 00049 C kfield : current field number 00050 C cdint : type of interpolation to be performed 00051 C cdtyp : type of source grid 00052 C cddim : type of field (scalar or vector) 00053 C 00054 C Output: 00055 C ------ 00056 C poutp : output field on target grid (real 2D) 00057 C 00058 C Workspace: 00059 C --------- 00060 C These are work arrays passed as arguments : 00061 C plonz, platz, psgrb, psgra, psfrb, psfra 00062 C 00063 C Externals: 00064 C --------- 00065 C rgoptc, namset, nagset, cxgaig, igscint, namsst, 00066 C nagsst, naflux. 00067 C 00068 C Reference: 00069 C --------- 00070 C See OASIS manual (1995) 00071 C 00072 C History: 00073 C ------- 00074 C Version Programmer Date Description 00075 C ------- ---------- ---- ----------- 00076 C 1.0 L. Terray 94/01/01 created 00077 C 1.1 L. Terray 94/08/01 modified: interpolation differ 00078 C for scalar or vector in fscint 00079 C 2.0 L. Terray 95/12/15 modified: new structure 00080 C 2.2 L. Terray 97/11/28 added: grid type A,B,L in fscint 00081 C 2.3 L. Terray 99/03/01 corrected: calls to na(g-m)sst with 00082 C inclusion of pointer ipdeb 00083 C 2.3 S. Valcke 99/04/30 added: printing levels 00084 C 2.3 L. Terray 99/08/11 modified: fscint extrapolation 00085 C for periodic Z grids 00086 C 2.3 L. Terray 99/09/15 changed periodicity variables 00087 C 00088 C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 00089 C 00090 C* ---------------- Include files and USE of modules--------------------------- 00091 C 00092 USE mod_unit 00093 USE mod_anais 00094 USE mod_parameter 00095 USE mod_timestep 00096 USE mod_printing 00097 C 00098 C* ---------------------------- Argument declarations ------------------- 00099 C 00100 REAL (kind=ip_realwp_p) poutp(ktlon,ktlat), 00101 $ ptlon(ktlon,ktlat), ptlat(ktlon,ktlat) 00102 REAL (kind=ip_realwp_p) pinpt(kslon,kslat), 00103 $ pslon(kslon,kslat), pslat(kslon,kslat) 00104 REAL (kind=ip_realwp_p) pssurf(kslon,kslat), 00105 $ ptsurf(ktlon,ktlat) 00106 REAL (kind=ip_realwp_p) plonz(ktlon), platz(ktlat), 00107 $ psgrb(kslon,kslat), 00108 $ psgra(ktlon,ktlat), psfrb(kslon,kslat), psfra(ktlon,ktlat) 00109 INTEGER (kind=ip_intwp_p) ktmsk(ktlon,ktlat), ksmsk(kslon,kslat) 00110 CHARACTER*8 cdsper, cdtper 00111 C 00112 C* ---------------------------- Local declarations ---------------------- 00113 C 00114 CHARACTER*1 cdtyp, clref, cltyp 00115 CHARACTER*6 cddim, cldim 00116 CHARACTER*8 cdint, clord, clflg, clind 00117 LOGICAL llinit, llsym 00118 C 00119 C* ---------------------------- Poema verses ---------------------------- 00120 C 00121 C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 00122 C 00123 C* 1. Initialization 00124 C -------------- 00125 C 00126 IF (nlogprt .GE. 2) THEN 00127 WRITE (UNIT = nulou,FMT = *) ' ' 00128 WRITE (UNIT = nulou,FMT = *) 00129 $ ' ROUTINE fiasco - Level 3' 00130 WRITE (UNIT = nulou,FMT = *) 00131 $ ' ************** *******' 00132 WRITE (UNIT = nulou,FMT = *) ' ' 00133 WRITE (UNIT = nulou,FMT = *) 00134 $ ' Interpolation package driving routine' 00135 WRITE (UNIT = nulou,FMT = *) ' ' 00136 CALL FLUSH(nulou) 00137 ENDIF 00138 C 00139 C* Assign local value to grid and field type variables to deal 00140 C with FSCINT set-up 00141 C 00142 cltyp = cdtyp 00143 cldim = cddim 00144 C 00145 C* set up interpolation method 00146 C 00147 IF (cdint .EQ. 'BILINEAR') THEN 00148 clflg = 'FSCINT' 00149 CALL rgoptc('INTERP', 'LINEAIR', .TRUE.) 00150 ELSEIF (cdint .EQ. 'BICUBIC') THEN 00151 clflg = 'FSCINT' 00152 CALL rgoptc('INTERP', 'CUBIQUE', .TRUE.) 00153 ELSEIF (cdint .eq. 'NNEIBOR') THEN 00154 clflg = 'FSCINT' 00155 CALL rgoptc('INTERP', 'VOISIN', .TRUE.) 00156 ELSEIF (cdint .EQ. 'SURFMESH') THEN 00157 clflg = 'ANAISM' 00158 llinit = linit(kfield) 00159 ELSEIF (cdint .EQ. 'GAUSSIAN') THEN 00160 clflg = 'ANAISG' 00161 llinit = linit(kfield) 00162 ELSE 00163 WRITE (UNIT = nulou,FMT = *) ' ***WARNING***' 00164 WRITE (UNIT = nulou,FMT = *) 00165 $ ' ===>>>> : unknown kind of interpolator' 00166 WRITE (UNIT = nulou,FMT = *) 00167 $ ' ======= ============' 00168 WRITE (UNIT = nulou,FMT = *) 00169 $ ' --->>> CDINT = ',cdint 00170 WRITE (UNIT = nulou,FMT = *) ' ' 00171 WRITE (UNIT = nulou,FMT = *) 00172 $ ' We STOP!!! check variable syntax ' 00173 WRITE (UNIT = nulou,FMT = *) ' ' 00174 CALL HALTE ('STOP in fiasco') 00175 ENDIF 00176 IF (clflg .eq. 'FSCINT') THEN 00177 CALL rgoptc('INTERP', clord, .FALSE.) 00178 ELSEIF (clflg .EQ. 'ANAISM') THEN 00179 clord = 'SURFMESH' 00180 ELSEIF (clflg .EQ. 'ANAISG') THEN 00181 clord = 'GAUSSIAN' 00182 ELSE 00183 clord ='UNKNOWN' 00184 ENDIF 00185 IF (nlogprt .GE. 2) THEN 00186 CALL prcout 00187 $ ('Current interpolator is clord = ', clord, 1) 00188 ENDIF 00189 C 00190 C* set up value to be given in case of extrapolation in fscint 00191 C 00192 IF (clflg .eq. 'FSCINT') THEN 00193 IF (cltyp .EQ. 'Z') THEN 00194 C 00195 C* No extrapolation for Z grids (global grid) 00196 C 00197 CALL rgoptc ('EXTRAP', 'OUI', .TRUE.) 00198 ELSE 00199 CALL rgoptc ('EXTRAP', 'VOISIN', .TRUE.) 00200 ENDIF 00201 CALL rgoptc ('EXTRAP', clind , .false.) 00202 IF (nlogprt .GE. 2) THEN 00203 CALL prcout 00204 $ ('Current extrapolation in FSCINT is clind =', clind, 1) 00205 ENDIF 00206 ENDIF 00207 C 00208 C* set up grid descriptors for fscint interpolator 00209 C 00210 IF (clflg .eq. 'FSCINT') THEN 00211 IF (cltyp .EQ. 'Z' .OR. cltyp .EQ. 'Y') THEN 00212 C 00213 C* here ig* descriptors are arbitrary (date of louis XVI's death) 00214 C 00215 ig1 = 21 00216 ig2 = 0 00217 ig3 = 1 00218 ig4 = 1793 00219 zdllg = 1.0 00220 zdllt = 1.0 00221 zlgin = 0.0 00222 zltin = 0.0 00223 clref = 'L' 00224 CALL cxgaig (clref, ig1z, ig2z, ig3z, ig4z, 00225 $ zltin, zlgin, zdllt, zdllg) 00226 DO 110 ji = 1, kslon 00227 plonz(ji) = pslon(ji,1) + 1.0 00228 110 CONTINUE 00229 DO 120 jj = 1, kslat 00230 platz(jj) = pslat(1,jj) + 1.0 00231 120 CONTINUE 00232 ELSE IF (cltyp .EQ. 'G'.or. cltyp .EQ. 'A' 00233 $ .OR. cltyp .EQ. 'B') THEN 00234 ig1 = 0 00235 ig2 = 0 00236 ig3 = 0 00237 ig4 = 0 00238 llsym = .TRUE. 00239 ELSE IF (cltyp .EQ. 'L') THEN 00240 zdllt = pslat(1,2) - pslat(1,1) 00241 zdllg = pslon(2,1) - pslon(1,1) 00242 zlgin = pslon(1,1) 00243 zltin = pslat(1,1) 00244 CALL cxgaig (cltyp, ig1, ig2, ig3, ig4, 00245 $ zltin, zlgin, zdllt, zdllg) 00246 llsym = .TRUE. 00247 ELSE 00248 WRITE (UNIT = nulou,FMT = *) ' ***WARNING***' 00249 WRITE (UNIT = nulou,FMT = *) 00250 $ ' ===>>> : grid type not implemented yet' 00251 WRITE (UNIT = nulou,FMT = *) 00252 $ ' ====== ==== =========== ' 00253 WRITE (UNIT = nulou,FMT = *) 00254 $ ' --->>> cltyp = ', 00255 $ cltyp 00256 WRITE (UNIT = nulou,FMT = *) ' ' 00257 WRITE (UNIT = nulou,FMT = *) 00258 $ ' We STOP !!! Check source grid type' 00259 WRITE (UNIT = nulou,FMT = *) ' ' 00260 CALL HALTE ('STOP in fiasco') 00261 ENDIF 00262 ELSEIF (clflg .eq. 'ANAISM') THEN 00263 IF (llinit) THEN 00264 C 00265 C* Mask values for both models 00266 C 00267 ismsq = 1 00268 itmsq = 1 00269 C 00270 C* If no intersection, give value imesh to nmesh array element 00271 C 00272 imesh = 0 00273 iwlun = nulcc 00274 C 00275 C* Get pointers for ANAISM interpolator parameters 00276 C 00277 ipdeb = (naismfl(kfield)-1)*ig_maxwoa*ig_maxgrd+1 00278 isdeb = (naismfl(kfield)-1)*ig_maxgrd+1 00279 C 00280 C* Set up surfmesh interpolation 00281 C 00282 CALL namset (pslon, pslat, ksmsk, ismsq, kslon, kslat, 00283 $ cdsper, ksper, 00284 $ ptlon, ptlat, ktmsk, itmsq, ktlon, ktlat, 00285 $ cdtper, ktper, 00286 $ cltyp, 00287 C 00288 C* Specific ANAISM parameters 00289 C 00290 $ amint(ipdeb), 00291 $ nmint(ipdeb), 00292 $ naismvoi(kfield), psgrb, psgra, 00293 $ niwtm(kfield), iwlun, imesh, 00294 $ nmesh(isdeb), 00295 $ naismfl(kfield)) 00296 C 00297 C* Switch initialization flag 00298 C 00299 linit(kfield) = .FALSE. 00300 ENDIF 00301 ELSEIF (clflg .eq. 'ANAISG') THEN 00302 IF (llinit) THEN 00303 ismsq = 1 00304 itmsq = 1 00305 iwlun = nulgg 00306 C 00307 C* Get pointers for Anaisg interpolator parameters 00308 C 00309 ipdeb = (naisgfl(kfield)-1)*ig_maxnoa*ig_maxgrd+1 00310 CALL nagset (pslon, pslat, pssurf, 00311 $ ksmsk, ismsq, kslon, kslat, 00312 $ ptlon, ptlat, ptsurf, 00313 $ ktmsk, itmsq, ktlon, ktlat, 00314 $ cltyp, 00315 C 00316 C* Specific Anaisg parameters 00317 C 00318 $ agint(ipdeb), 00319 $ ngint(ipdeb), 00320 $ naisgvoi(kfield), psfrb, psfra, 00321 $ varmul(kfield), niwtg(kfield), iwlun, 00322 $ naisgfl(kfield)) 00323 C 00324 C* Switch off initialization flag 00325 C 00326 linit(kfield) = .FALSE. 00327 ENDIF 00328 ELSE 00329 WRITE (UNIT = nulou,FMT = *) ' ***WARNING***' 00330 WRITE (UNIT = nulou,FMT = *) 00331 $ ' ===>>> : unknown interpolation package' 00332 WRITE (UNIT = nulou,FMT = *) 00333 $ ' ====== ======= =======' 00334 WRITE (UNIT = nulou,FMT = *) 00335 $ ' --->>> clflg = ', 00336 $ clflg 00337 WRITE (UNIT = nulou,FMT = *) ' ' 00338 WRITE (UNIT = nulou,FMT = *) 00339 $ ' We STOP !!! Check interpolation package' 00340 WRITE (UNIT = nulou,FMT = *) ' ' 00341 CALL HALTE ('STOP in fiasco') 00342 ENDIF 00343 C 00344 C 00345 C* 2. Interpolation 00346 C ------------- 00347 C 00348 IF (clflg .eq. 'FSCINT') THEN 00349 IF (cltyp .eq. 'Z' .OR. cltyp .EQ. 'Y') THEN 00350 CALL igscint (poutp, ktlon, ktlat, ptlat, ptlon, 00351 $ pinpt, kslon, kslat, cltyp, clref, 00352 $ ig1z, ig2z, ig3z, ig4z, .true., 00353 $ plonz, platz, cldim) 00354 ELSEIF (cltyp .eq. 'G'.or. cltyp .eq. 'A' 00355 $ .OR. cltyp .eq. 'B'.or. cltyp .eq. 'L') THEN 00356 CALL rgscint (poutp, ktlon, ktlat, ptlat, ptlon, 00357 $ pinpt, kslon, kslat, cltyp, 00358 $ ig1, ig2, ig3, ig4, llsym, cldim) 00359 ELSE 00360 WRITE (UNIT = nulou,FMT = *) ' ***WARNING***' 00361 WRITE (UNIT = nulou,FMT = *) 00362 $ ' ===>>> grid type not implemented yet' 00363 WRITE (UNIT = nulou,FMT = *) 00364 $ ' ====== ==== =========== ' 00365 WRITE (UNIT = nulou,FMT = *) 00366 $ ' --->>> cltyp = ', 00367 $ cltyp 00368 WRITE (UNIT = nulou,FMT = *) ' ' 00369 WRITE (UNIT = nulou,FMT = *) 00370 $ ' WE STOP !!! check source grid type' 00371 WRITE (UNIT = nulou,FMT = *) ' ' 00372 CALL HALTE ('STOP in fiasco') 00373 ENDIF 00374 ELSEIF (clflg .eq. 'ANAISM') THEN 00375 itmsq = 1 00376 ipdeb = (naismfl(kfield)-1)*ig_maxwoa*ig_maxgrd+1 00377 CALL namsst (poutp, ktmsk, itmsq, ktlon, ktlat, 00378 $ amint(ipdeb), nmint(ipdeb), naismvoi(kfield), 00379 $ pinpt, kslon, kslat) 00380 ELSEIF (clflg .eq. 'ANAISG') THEN 00381 itmsq = 1 00382 ipdeb = (naisgfl(kfield)-1)*ig_maxnoa*ig_maxgrd+1 00383 CALL nagsst (poutp, ktmsk, itmsq, ktlon, ktlat, 00384 $ agint(ipdeb), ngint(ipdeb), naisgvoi(kfield), 00385 $ pinpt, kslon, kslat) 00386 ELSE 00387 WRITE (UNIT = nulou,FMT = *) ' ***WARNING***' 00388 WRITE (UNIT = nulou,FMT = *) 00389 $ ' ===>>> : unknown interpolation package' 00390 WRITE (UNIT = nulou,FMT = *) 00391 $ ' ====== ======= =======' 00392 WRITE (UNIT = nulou,FMT = *) 00393 $ ' --->>> clflg = ', 00394 $ clflg 00395 WRITE (UNIT = nulou,FMT = *) ' ' 00396 WRITE (UNIT = nulou,FMT = *) 00397 $ ' We STOP !!! Check interpolation package' 00398 WRITE (UNIT = nulou,FMT = *) ' ' 00399 CALL HALTE ('STOP in fiasco') 00400 ENDIF 00401 C 00402 C 00403 C* 3. End of routine 00404 C -------------- 00405 C 00406 IF (nlogprt .GE. 2) THEN 00407 WRITE (UNIT = nulou,FMT = *) ' ' 00408 WRITE (UNIT = nulou,FMT = *) 00409 $ ' --------- End of ROUTINE fiasco ---------' 00410 CALL FLUSH (nulou) 00411 ENDIF 00412 RETURN 00413 END