[2] | 1 | ;+ |
---|
| 2 | ; |
---|
[142] | 3 | ; @file_comments |
---|
| 4 | ; Overprint vectors in a field traced by plt. |
---|
[2] | 5 | ; |
---|
[142] | 6 | ; @categories |
---|
[157] | 7 | ; Graphics |
---|
[2] | 8 | ; |
---|
[163] | 9 | ; @param VECTEUR {in}{required}{type=vector} |
---|
[231] | 10 | ; It is a structure with 2 elements containing we 2 matrices U and V of |
---|
| 11 | ; values of the zonal and meridian component of the field of vectors to |
---|
[142] | 12 | ; be traced. |
---|
| 13 | ; For ex: |
---|
[2] | 14 | ; vecteur={matriceu:lec('unsurface'),matricev:lec('vnsurface')} |
---|
[142] | 15 | ; rq:the name of elements of vector does not have any importance. |
---|
| 16 | ; vecteur={u:lec('unsurface'),v:lec('vnsurface')} goes well too. |
---|
[2] | 17 | ; |
---|
[163] | 18 | ; @keyword UNVECTSUR {type=scalar or array} |
---|
[142] | 19 | ; It is a scalar n or an array with 2 elements [n1,n2]. |
---|
| 20 | ; In the first case, we will trace a vector on n following x and y. |
---|
[231] | 21 | ; In the second case, we will trace a vector on n1 following x and a |
---|
[142] | 22 | ; vector n2 following n2 |
---|
[231] | 23 | ; Comments: To trace all vectors following y and one vector on two |
---|
[142] | 24 | ; following x, put unvectsur=[2,1] |
---|
[2] | 25 | ; |
---|
[142] | 26 | ; @keyword VECTMIN {in}{required} |
---|
[231] | 27 | ; Minimum norme of vectors to be traced |
---|
[2] | 28 | ; |
---|
[142] | 29 | ; @keyword VECTMAX {in}{required} |
---|
[231] | 30 | ; Maximum norme of vectors to be traced |
---|
[2] | 31 | ; |
---|
[142] | 32 | ; @keyword _EXTRA |
---|
[231] | 33 | ; Used to pass keywords |
---|
| 34 | ; |
---|
[142] | 35 | ; @uses |
---|
| 36 | ; common.pro |
---|
[2] | 37 | ; |
---|
[142] | 38 | ; @history |
---|
[157] | 39 | ; Sebastien Masson (smasson\@lodyc.jussieu.fr) |
---|
[2] | 40 | ; 10/3/1999 |
---|
| 41 | ; 11/6/1999 compatibilite avec NAN et la lecture |
---|
| 42 | ; des structures. |
---|
[142] | 43 | ; |
---|
| 44 | ; @version |
---|
[231] | 45 | ; $Id$ |
---|
| 46 | ; |
---|
[2] | 47 | ;- |
---|
[114] | 48 | ; |
---|
[231] | 49 | PRO ajoutvect,vecteur, vectlegende, UNVECTSUR=unvectsur,VECTMIN=vectmin, VECTMAX=vectmax, _EXTRA = ex |
---|
| 50 | ; |
---|
[114] | 51 | compile_opt idl2, strictarrsubs |
---|
| 52 | ; |
---|
[2] | 53 | @common |
---|
[142] | 54 | tempsun = systime(1) ; For key_performance |
---|
[2] | 55 | ;---------------------------------------------------------------------------- |
---|
| 56 | ; |
---|
| 57 | u = litchamp(vecteur.(0)) |
---|
| 58 | u = checkfield(u, 'plt', TYPE = 'xy', /NOQUESTION) |
---|
| 59 | v = litchamp(vecteur.(1)) |
---|
| 60 | v = checkfield(v, 'plt', TYPE = 'xy', /NOQUESTION) |
---|
| 61 | ;----------------------------------------------------------- |
---|
[142] | 62 | ;----------------------------------------------------------- |
---|
| 63 | ; We recuperate possible informations on fields |
---|
| 64 | ;----------------------------------------------------------- |
---|
[2] | 65 | grilleu = litchamp(vecteur.(0), /grid) |
---|
| 66 | if grilleu EQ '' then grilleu = 'U' |
---|
| 67 | grillev = litchamp(vecteur.(1), /grid) |
---|
| 68 | if grillev EQ '' then grillev = 'V' |
---|
[41] | 69 | |
---|
[2] | 70 | IF grilleu EQ 'V' AND grillev EQ 'U' THEN inverse = 1 |
---|
| 71 | IF grilleu EQ grillev THEN interpolle = 0 ELSE interpolle = 1 |
---|
| 72 | if keyword_set(inverse) then begin |
---|
| 73 | rien = u |
---|
| 74 | u = v |
---|
| 75 | v = rien |
---|
| 76 | endif |
---|
| 77 | ;------------------------------------------------------------ |
---|
[142] | 78 | ; We find common points between u and v |
---|
[2] | 79 | ;------------------------------------------------------------ |
---|
| 80 | if interpolle then begin |
---|
[41] | 81 | indicexu = (lindgen(jpi))[firstxu:firstxu+nxu-1] |
---|
| 82 | indicexv = (lindgen(jpi))[firstxv:firstxv+nxv-1] |
---|
[2] | 83 | indicex = inter(indicexu, indicexv) |
---|
[41] | 84 | indiceyu = (lindgen(jpj))[firstyu:firstyu+nyu-1] |
---|
| 85 | indiceyv = (lindgen(jpj))[firstyv:firstyv+nyv-1] |
---|
[2] | 86 | indicey = inter(indiceyu, indiceyv) |
---|
[231] | 87 | nx = n_elements(indicex) |
---|
[2] | 88 | ny = n_elements(indicey) |
---|
| 89 | indice2d = lindgen(jpi, jpj) |
---|
| 90 | indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1] |
---|
| 91 | ;------------------------------------------------------------ |
---|
[142] | 92 | ; extraction of u and v on the appropriated domain |
---|
[2] | 93 | ;------------------------------------------------------------ |
---|
| 94 | case 1 of |
---|
| 95 | (size(u))[0] NE 2 OR (size(v))[0] NE 2: return |
---|
| 96 | (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $ |
---|
[231] | 97 | (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN |
---|
[2] | 98 | if nxu NE nx then $ |
---|
[231] | 99 | if indicex[0] EQ firstxu then u = u[0:nx-1, *] ELSE u = u[1: nx, *] |
---|
[2] | 100 | IF nxv NE nx THEN $ |
---|
[41] | 101 | if indicex[0] EQ firstxv then v = v[0:nx-1, *] ELSE v = v[1: nx, *] |
---|
[2] | 102 | IF nyu NE ny THEN $ |
---|
[231] | 103 | if indicey[0] EQ firstyu then u = u[*, 0:ny-1] ELSE u = u[*, 1: ny] |
---|
[2] | 104 | IF nyv NE ny THEN $ |
---|
[41] | 105 | if indicey[0] EQ firstyv then v = v[*, 0:ny-1] ELSE v = v[*, 1: ny] |
---|
[2] | 106 | END |
---|
| 107 | (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $ |
---|
[231] | 108 | (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN |
---|
[2] | 109 | u = u[indice2d] |
---|
| 110 | v = v[indice2d] |
---|
| 111 | END |
---|
[231] | 112 | ELSE:BEGIN |
---|
[2] | 113 | ras = report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs') |
---|
| 114 | return |
---|
| 115 | end |
---|
| 116 | endcase |
---|
| 117 | ;------------------------------------------------------------------ |
---|
[142] | 118 | ; We reshape u and v to make sure that none dimension has been erased. |
---|
[2] | 119 | ;------------------------------------------------------------------ |
---|
| 120 | if ny EQ 1 then begin |
---|
| 121 | u = reform(u, nx, ny) |
---|
| 122 | v = reform(v, nx, ny) |
---|
| 123 | endif |
---|
| 124 | ;------------------------------------------------------------------ |
---|
[142] | 125 | ; construction of u and v at points T |
---|
[2] | 126 | ;----------------------------------------------------------- |
---|
[114] | 127 | a=u[0,*] |
---|
[2] | 128 | u=(u+shift(u,1,0))/2. |
---|
[114] | 129 | if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*]=a |
---|
| 130 | a=v[*,0] |
---|
[2] | 131 | v=(v+shift(v,0,1))/2. |
---|
[114] | 132 | if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0]=a |
---|
[2] | 133 | ;---------------------------------------------------------------------------- |
---|
[142] | 134 | ; attribution of the mask and of longitude and latitude arrays. |
---|
[231] | 135 | ; We recuperate the complete grid to establish a big mask extensive |
---|
| 136 | ; in the four directions to cover points for which a land point has |
---|
[142] | 137 | ; been considerated (do a small drawing) |
---|
[2] | 138 | ;---------------------------------------------------------------------------- |
---|
| 139 | vargrid='T' |
---|
[41] | 140 | msku = (umask())[indice2d+jpi*jpj*firstzt] |
---|
| 141 | mskv = (vmask())[indice2d+jpi*jpj*firstzt] |
---|
[2] | 142 | glam = glamt[indice2d] |
---|
| 143 | gphi = gphit[indice2d] |
---|
| 144 | if ny EQ 1 then begin |
---|
| 145 | msku = reform(msku, nx, ny) |
---|
| 146 | mskv = reform(mskv, nx, ny) |
---|
| 147 | ; glam = reform(glam, nx, ny) |
---|
| 148 | ; gphi = reform(gphi, nx, ny) |
---|
| 149 | endif |
---|
| 150 | ;----------------------------------------------------------- |
---|
[296] | 151 | ; We mask u and v the long of coasts (the place where we |
---|
[142] | 152 | ; can not calculate the average) |
---|
[231] | 153 | ;----------------------------------------------------------- |
---|
[297] | 154 | ; extension of the mask |
---|
[2] | 155 | u = u*msku*shift(msku,1,0) |
---|
| 156 | v = v*mskv*shift(mskv,0,1) |
---|
[231] | 157 | ENDIF ELSE BEGIN |
---|
[41] | 158 | u = u*tmask[firstxt:lastxt,firstyt:lastyt,firstzt] |
---|
| 159 | v = v*tmask[firstxt:lastxt,firstyt:lastyt,firstzt] |
---|
[2] | 160 | indice2d = lindgen(jpi, jpj) |
---|
[41] | 161 | indice2d = indice2d[firstxt:lastxt, firstyt:lastyt] |
---|
[2] | 162 | nx = nxt |
---|
| 163 | ny = nyt |
---|
| 164 | endelse |
---|
| 165 | tabnorme=sqrt(u^2+v^2) |
---|
| 166 | nan = where(finite(u, /nan) EQ 1) |
---|
| 167 | if nan[0] NE -1 then u[nan] = 1e5 |
---|
| 168 | nan = where(finite(v, /nan) EQ 1) |
---|
| 169 | if nan[0] NE -1 then v[nan] = 1e5 |
---|
[231] | 170 | if keyword_set(vectmin) then BEGIN |
---|
[41] | 171 | toosmall=where(tabnorme lt vectmin) |
---|
| 172 | if toosmall[0] NE -1 then begin |
---|
| 173 | u[toosmall] = 1e5 |
---|
| 174 | v[toosmall] = 1e5 |
---|
[2] | 175 | ENDIF |
---|
| 176 | endif |
---|
| 177 | if keyword_set(vectmax) then BEGIN |
---|
[41] | 178 | toobig=where(tabnorme gt vectmax) |
---|
| 179 | if toobig[0] NE -1 then begin |
---|
| 180 | u[toobig] = 1e5 |
---|
| 181 | v[toobig] = 1e5 |
---|
[2] | 182 | ENDIF |
---|
| 183 | ENDIF |
---|
| 184 | ;----------------------------------------------------------- |
---|
[142] | 185 | ; Put back of a big value on all points for which we can do the calculation. |
---|
[2] | 186 | ;----------------------------------------------------------- |
---|
| 187 | if interpolle then t2 = msku*shift(msku,1,0)*mskv*shift(mskv,0,1) $ |
---|
[41] | 188 | ELSE t2 = tmask[firstxt:lastxt,firstyt:lastyt,firstzt] |
---|
[114] | 189 | if NOT keyword_set(key_periodic) OR nx NE jpi then t2[0, *]=0. |
---|
| 190 | t2[*,0]=0. |
---|
[2] | 191 | terre=where(t2 eq 0) |
---|
| 192 | if terre[0] ne -1 then begin |
---|
[114] | 193 | u[terre]=1e5 |
---|
| 194 | v[terre]=1e5 |
---|
[2] | 195 | ENDIF |
---|
| 196 | ;----------------------------------------------------------- |
---|
[142] | 197 | ; trace only one vector one two |
---|
[2] | 198 | ;----------------------------------------------------------- |
---|
| 199 | if keyword_set(unvectsur) then BEGIN ; |
---|
[142] | 200 | ; indx is a vector containing number of columns to be selected. |
---|
| 201 | ; indy is a vector containing number of lines to be selected. |
---|
[2] | 202 | if n_elements(unvectsur) EQ 1 then begin |
---|
| 203 | indx = where((lindgen(nx) MOD unvectsur[0]) eq 0) |
---|
| 204 | indy = where((lindgen(ny) MOD unvectsur[0]) eq 0) |
---|
[231] | 205 | ENDIF ELSE BEGIN |
---|
[2] | 206 | indx = where((lindgen(nx) MOD unvectsur[0]) eq 0) |
---|
| 207 | indy = where((lindgen(ny) MOD unvectsur[1]) eq 0) |
---|
[142] | 208 | ENDELSE |
---|
| 209 | ; From indx and indy, we will construct an array which will give indexes |
---|
| 210 | ; of intersections points of columns specified by indx. |
---|
[2] | 211 | indicereduit = indx#replicate(1,n_elements(indy))+nx*replicate(1,n_elements(indx))#indy |
---|
[142] | 212 | ; We reduce arrays which will be passed to vecteur. |
---|
[2] | 213 | u = u[indicereduit] |
---|
| 214 | v = v[indicereduit] |
---|
| 215 | tabnorme = tabnorme[indicereduit] |
---|
[241] | 216 | t2 = t2[indicereduit] |
---|
[2] | 217 | ; |
---|
| 218 | endif |
---|
| 219 | ;----------------------------------------------------------- |
---|
[231] | 220 | ; |
---|
[2] | 221 | ;----------------------------------------------------------- |
---|
| 222 | if keyword_set(inverse) then begin |
---|
| 223 | rien = u |
---|
| 224 | u = v |
---|
| 225 | v = rien |
---|
| 226 | endif |
---|
| 227 | ;----------------------------------------------------------- |
---|
[142] | 228 | ; Drawing of vectors. |
---|
[2] | 229 | ;---------------------------------------------------------- |
---|
| 230 | vecteur, u, v, tabnorme, indice2d, indicereduit, missing=1e5, _extra = ex |
---|
| 231 | ;----------------------------------------------------------- |
---|
[142] | 232 | ; We complete the caption. |
---|
[2] | 233 | ;----------------------------------------------------------- |
---|
| 234 | if terre[0] ne -1 then mini = min(tabnorme[where(t2 eq 1)], max = maxi, /nan) $ |
---|
| 235 | ELSE mini = min(tabnorme, max = maxi, /nan) |
---|
[231] | 236 | |
---|
[2] | 237 | if litchamp(vecteur.(0), /u) NE '' then $ |
---|
| 238 | vectlegende = {minmax:[mini, maxi], unite:litchamp(vecteur.(0), /u)} $ |
---|
| 239 | ELSE vectlegende = {minmax:[mini, maxi], unite:varunit} |
---|
| 240 | |
---|
| 241 | |
---|
| 242 | sortie: |
---|
[231] | 243 | if keyword_set(key_performance) NE 0 THEN print, 'temps ajoutvect', systime(1)-tempsun |
---|
[2] | 244 | return |
---|
| 245 | end |
---|