|
Oasis3 4.0.2
|
|
Oasis3 4.0.2
|
00001 subroutine fracnnei (src_size, dst_size, 00002 $ ld_srcmask, ld_dstmask, 00003 $ src_lon, src_lat, dst_lon, dst_lat, 00004 $ num_links, num_wgts, 00005 $ weights, src_addr, dst_addr) 00006 C**** 00007 C ***************************** 00008 C * OASIS ROUTINE - LEVEL 4 * 00009 C * ------------- ------- * 00010 C ***************************** 00011 C 00012 C**** *fracnnei* - SCRIP remapping 00013 C 00014 C Purpose: 00015 C ------- 00016 C Treatment of the tricky points in an interpolation 00017 C 00018 C Interface: 00019 C --------- 00020 C *CALL* * 00021 C 00022 C Called from: 00023 C ----------- 00024 C scriprmp 00025 C 00026 C Input: 00027 C ----- 00028 C src_size : source grid size (integer) 00029 C dst_size : target grid size (integer) 00030 C ld_srcmask : mask of the source grid 00031 C ld_dstmask : mask of the target grid 00032 C src_lon : longitudes of the source grid 00033 C src_lat : latitudes of the source grid 00034 C dst_lon : longitudes of the target grid 00035 C dst_lat : latitudes of the target grid 00036 C num_links : total number of links 00037 C num_wgts : number of weights for each link 00038 C InOut 00039 C ----- 00040 C weights : remapping weights 00041 C src_addr : remapping source addresses 00042 C dst_addr : remapping target addresses 00043 C 00044 C History: 00045 C ------- 00046 C Version Programmer Date Description 00047 C ------- ---------- ---- ----------- 00048 C 2.5 D.Declat 2002/08/20 adapted from S. Valcke ptmsq 00049 C 3.0 S. Valcke 2002/10/30 test and corrections 00050 C 00051 C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 00052 C* ---------------------------- Modules used ---------------------------- 00053 C 00054 use kinds_mod ! defines common data types 00055 use constants ! defines common constants 00056 use grids ! module containing grid information 00057 use remap_vars ! module containing remap information 00058 C 00059 C* ---------------------------- Implicit -------------------------------- 00060 C 00061 implicit none 00062 C 00063 C* ---------------------------- Include files --------------------------- 00064 C 00065 00066 C INCLUDE 'netcdf.inc' 00067 C 00068 C* ---------------------------- Intent In ------------------------------- 00069 C 00070 INTEGER (kind=int_kind) :: 00071 $ src_size, ! size of the source grid 00072 $ dst_size ! size of the destination grid 00073 C 00074 REAL (kind=dbl_kind) :: 00075 $ src_lat(src_size), src_lon(src_size), 00076 $ dst_lat(dst_size), dst_lon(dst_size) 00077 00078 C 00079 LOGICAL :: 00080 $ ld_srcmask(src_size), ! source grid mask 00081 $ ld_dstmask(dst_size) ! target grid mask 00082 C 00083 INTEGER (kind=int_kind) :: 00084 $ num_links, ! number of links between src and tgt 00085 $ num_wgts ! number of weights 00086 00087 C 00088 C* ---------------------------- Intent InOut ------------------------------ 00089 C 00090 REAL (kind=dbl_kind) :: 00091 $ weights(num_wgts, num_links) ! remapping weights 00092 C 00093 INTEGER (kind=int_kind) :: 00094 $ src_addr(num_links), ! remapping source addresses 00095 $ dst_addr(num_links) ! remapping target addresses 00096 C 00097 C* ---------------------------- Local declarations ---------------------- 00098 C 00099 C 00100 INTEGER (kind=int_kind) :: 00101 $ ila_nneiadd ! Nearest-neighbor address 00102 C 00103 INTEGER (kind=int_kind) :: 00104 $ ib_dst, ! INDEX loop for the distance grid 00105 $ ib_src, ! INDEX loop for the source grid 00106 $ ib_neigh, ! INDEX loop on the corresponding src pts 00107 $ ib_links ! INDEX loop for the links 00108 C 00109 INTEGER (kind=int_kind) :: 00110 $ nb_Vmm, ! number of Mtt points 00111 $ ntottotland, ! number of land points 00112 $ ntotland, ! number of land points 00113 $ ntotoce, ! number of oceanic points 00114 $ ntotngh ! number of corresponding source points 00115 C 00116 INTEGER (kind=int_kind) :: 00117 $ beg_links, ! begining of the serie of links 00118 $ nb_links ! number of links 00119 C 00120 REAL (kind=dbl_kind) :: 00121 $ coslat, ! cosinus of the latitude 00122 $ sinlat, ! sinus of the latitude 00123 $ coslon, ! cosinus of the longitude 00124 $ sinlon, ! sinus of the longitude 00125 $ distance, 00126 $ dist_min, 00127 $ arg 00128 C 00129 INTEGER (kind=int_kind) :: n, il 00130 C 00131 C* ---------------------------- Poema verses ---------------------------- 00132 C 00133 C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 00134 C 00135 C* 1. Initialization 00136 C -------------- 00137 C 00138 IF (nlogprt .GE. 2) THEN 00139 WRITE (UNIT = nulou,FMT = *) ' ' 00140 WRITE (UNIT = nulou,FMT = *) ' ' 00141 WRITE (UNIT = nulou,FMT = *) 00142 $ ' Entering ROUTINE fracnnei - Level 4' 00143 WRITE (UNIT = nulou,FMT = *) 00144 $ ' **************** *******' 00145 WRITE (UNIT = nulou,FMT = *) ' ' 00146 WRITE (UNIT = nulou,FMT = *) 00147 $ ' Treating the tricky points of the remapping' 00148 WRITE (UNIT = nulou,FMT = *) ' ' 00149 CALL FLUSH(nulou) 00150 ENDIF 00151 C 00152 C *---------------------------------------------------------------------- 00153 C 00154 C* 2. Treating Vmm points V 00155 C ------------------- m m 00156 C The target point is a non-masked Valid point while the source points 00157 C are all masked points. Use of the non-masked nearest neighbours. 00158 C 00159 nb_Vmm = 0 00160 ntottotland = 0 00161 beg_links = 1 00162 00163 C* -- Loop all other the target points 00164 DO ib_dst = 1, dst_size 00165 C* -- If the point is a sea point 00166 IF (ld_dstmask(ib_dst)) THEN 00167 00168 beg_links = 0 00169 00170 DO ib_links = 1, num_links 00171 IF (dst_addr(ib_links) .eq. ib_dst) THEN 00172 beg_links = ib_links 00173 exit 00174 ENDIF 00175 END DO 00176 00177 IF (beg_links .ne. 0) THEN 00178 00179 ntotland = 0 00180 ntotoce = 0 00181 ntotngh = 0 00182 00183 C* -- Find the number of associated src points to the non-masked tgt point 00184 nb_links = 1 00185 DO il = beg_links+1, num_links 00186 IF (dst_addr(il) .eq. dst_addr(beg_links)) 00187 $ nb_links = nb_links + 1 00188 END DO 00189 00190 C* -- For each point on the src grid associated to the non-masked tgt point 00191 DO ib_neigh = 1, nb_links 00192 C 00193 ntotngh = ntotngh + 1 00194 00195 C* -- Check IF the point is a masked or non-masked point and treat it 00196 IF (.not. ld_srcmask(src_addr(beg_links+ib_neigh-1))) 00197 $ THEN 00198 ntotland = ntotland + 1 00199 ntottotland = ntottotland + 1 00200 ELSE IF (ld_srcmask(src_addr(beg_links+ib_neigh-1))) 00201 $ THEN 00202 ntotoce = ntotoce + 1 00203 ELSE 00204 WRITE (nulou, *) 'Pb with ocean mask with Mtt 1' 00205 END IF 00206 00207 END DO 00208 00209 C* -- If all the src points are land, treat it ! 00210 IF (ntotland .EQ. ntotngh) THEN 00211 nb_Vmm = nb_Vmm + 1 00212 WRITE(nulou,*) 00213 $ '************ Doing FRACNNEI for point', ib_dst 00214 00215 C* -- Find the nearest neighbours and change weights and address 00216 00217 coslat = cos(dst_lat(ib_dst)) 00218 sinlat = sin(dst_lat(ib_dst)) 00219 coslon = cos(dst_lon(ib_dst)) 00220 sinlon = sin(dst_lon(ib_dst)) 00221 00222 dist_min = bignum 00223 ila_nneiadd = 0 00224 DO ib_src = 1, src_size 00225 IF (ld_srcmask(ib_src)) THEN 00226 arg = 00227 & coslat*cos(src_lat(ib_src))* 00228 & (coslon*cos(src_lon(ib_src)) + 00229 & sinlon*sin(src_lon(ib_src)))+ 00230 & sinlat*sin(src_lat(ib_src)) 00231 IF (arg < -1.0d0) THEN 00232 arg = -1.0d0 00233 ELSE IF (arg > 1.0d0) THEN 00234 arg = 1.0d0 00235 END IF 00236 distance = acos(arg) 00237 IF (distance < dist_min) THEN 00238 ila_nneiadd = ib_src 00239 dist_min = distance 00240 ENDIF 00241 ENDIF 00242 END DO 00243 src_addr(beg_links) = ila_nneiadd 00244 weights(1,beg_links) = 1.0 00245 WRITE(nulou,*) 00246 $ '*************** Nearest source neighbour is ', 00247 $ ila_nneiadd 00248 IF (nb_links > 1) then 00249 DO n=2, nb_links 00250 src_addr(beg_links+n-1)=0 00251 weights(1,beg_links+n-1)=0.0 00252 END DO 00253 END IF 00254 END IF 00255 END IF 00256 ENDIF 00257 END DO 00258 C 00259 C 00260 C *---------------------------------------------------------------------- 00261 C 00262 IF (nlogprt .GE. 2) THEN 00263 WRITE (UNIT = nulou,FMT = *) ' ' 00264 WRITE (UNIT = nulou,FMT = *) 00265 $ ' Leaving ROUTINE fracnnei - Level 4' 00266 WRITE (UNIT = nulou,FMT = *) ' ' 00267 CALL FLUSH(nulou) 00268 ENDIF 00269 00270 END SUBROUTINE fracnnei 00271 00272 !*********************************************************************** 00273 00274
1.7.4 
