[2] | 1 | ;+ |
---|
| 2 | ; |
---|
[41] | 3 | ; @history |
---|
| 4 | ; Sebastien Masson (smasson\@lodyc.jussieu.fr) |
---|
[2] | 5 | ; |
---|
| 6 | ;- |
---|
| 7 | FUNCTION grad, field, direc |
---|
[41] | 8 | ; |
---|
[2] | 9 | @common |
---|
[41] | 10 | ; |
---|
[2] | 11 | res = litchamp(field) |
---|
| 12 | taille=size(res) |
---|
| 13 | grille, mask, glam, gphi, gdep, nx, ny,nz,premierx,premiery,premierz,dernierx, derniery, dernierz |
---|
| 14 | if n_elements(valmask) EQ 0 then valmask = 1e20 |
---|
| 15 | case strupcase(vargrid) of |
---|
| 16 | 'T':BEGIN |
---|
| 17 | case direc of |
---|
| 18 | 'x':BEGIN |
---|
| 19 | divi = e1u[premierx:dernierx, premiery:derniery] |
---|
| 20 | newmask = (umask())[premierx:dernierx, premiery:derniery, premierz:dernierz] |
---|
| 21 | vargrid = 'U' |
---|
| 22 | domdef, glamt[premierx, 0], glamu[dernierx, 0] $ |
---|
| 23 | , gphit[0, premiery], gphiu[0, derniery], grille = ['T','U'] |
---|
| 24 | END |
---|
| 25 | 'y':BEGIN |
---|
| 26 | divi = e2v[premierx:dernierx, premiery:derniery] |
---|
| 27 | newmask = (vmask())[premierx:dernierx, premiery:derniery, premierz:dernierz] |
---|
| 28 | vargrid = 'V' |
---|
| 29 | domdef, glamt[premierx, 0], glamv[dernierx, 0] $ |
---|
| 30 | , gphit[0, premiery], gphiv[0, derniery], grille = ['T','V'] |
---|
| 31 | END |
---|
| 32 | 'z':BEGIN |
---|
| 33 | divi = e3w[premierz:dernierz] |
---|
| 34 | newmask = mask |
---|
| 35 | vargrid = 'W' |
---|
| 36 | END |
---|
| 37 | ELSE:return, report('Bad definition of direction argument') |
---|
| 38 | ENDCASE |
---|
| 39 | END |
---|
| 40 | 'W':BEGIN |
---|
| 41 | case direc of |
---|
| 42 | 'x':divi = e1u[premierx:dernierx, premiery:derniery] |
---|
| 43 | 'y':divi = e2v[premierx:dernierx, premiery:derniery] |
---|
| 44 | 'z':BEGIN |
---|
| 45 | divi = e3t[premierz:dernierz] |
---|
| 46 | newmask = mask |
---|
| 47 | vargrid = 'T' |
---|
| 48 | END |
---|
| 49 | ELSE:return, report('Bad definition of direction argument') |
---|
| 50 | endcase |
---|
| 51 | END |
---|
| 52 | 'U':BEGIN |
---|
| 53 | case direc of |
---|
| 54 | 'x':BEGIN |
---|
| 55 | divi = (shift(e1t, -1, 0))[premierx:dernierx, premiery:derniery] |
---|
| 56 | newmask = tmask[premierx:dernierx, premiery:derniery, premierz:dernierz] |
---|
| 57 | vargrid = 'T' |
---|
| 58 | domdef, glamt[premierx, 0], glamu[dernierx] $ |
---|
| 59 | , gphit[0, premiery], gphiu[0, derniery], grille = ['T','U'] |
---|
| 60 | END |
---|
| 61 | 'y':BEGIN |
---|
| 62 | divi = e2f[premierx:dernierx, premiery:derniery] |
---|
| 63 | newmask = (fmask())[premierx:dernierx, premiery:derniery, premierz:dernierz] |
---|
| 64 | vargrid = 'F' |
---|
| 65 | domdef, glamu[premierx, 0], glamf[dernierx, 0] $ |
---|
| 66 | , gphiu[0, premiery], gphif[0, derniery], grille = ['U','F'] |
---|
| 67 | END |
---|
| 68 | 'z':BEGIN |
---|
| 69 | divi = e3w[premierz:dernierz] |
---|
| 70 | newmask = mask |
---|
| 71 | vargrid = 'W' |
---|
| 72 | END |
---|
| 73 | ELSE:return, report('Bad definition of direction argument') |
---|
| 74 | endcase |
---|
| 75 | END |
---|
| 76 | 'V':BEGIN |
---|
| 77 | case direc of |
---|
| 78 | 'x':BEGIN |
---|
| 79 | divi = e1f[premierx:dernierx, premiery:derniery] |
---|
| 80 | newmask = (fmask())[premierx:dernierx, premiery:derniery, premierz:dernierz] |
---|
| 81 | vargrid = 'F' |
---|
| 82 | domdef, glamv[premierx, 0], glamf[dernierx, 0] $ |
---|
| 83 | , gphiv[0, premiery], gphif[0, derniery], grille = ['V','F'] |
---|
| 84 | END |
---|
| 85 | 'y':BEGIN |
---|
| 86 | divi = (shift(e2t, 0, -1))[premierx:dernierx, premiery:derniery] |
---|
| 87 | newmask = tmask[premierx:dernierx, premiery:derniery, premierz:dernierz] |
---|
| 88 | vargrid = 'T' |
---|
| 89 | domdef, glamt[premierx, 0], glamv[dernierx, 0] $ |
---|
| 90 | , gphit[0, premiery], gphiv[0, derniery], grille = ['T','V'] |
---|
| 91 | END |
---|
| 92 | 'z':BEGIN |
---|
| 93 | divi = e3w[premierz:dernierz] |
---|
| 94 | newmask = mask |
---|
| 95 | vargrid = 'W' |
---|
| 96 | END |
---|
| 97 | ELSE:return, report('Bad definition of direction argument') |
---|
| 98 | endcase |
---|
| 99 | END |
---|
| 100 | ; 'F':BEGIN |
---|
| 101 | ; case direc of |
---|
| 102 | ; 'x':divi = (shift(e1v, -1, 0))[premierx:dernierx, premiery:derniery] |
---|
| 103 | ; 'y':divi = (shift(e2u, 0, -1))[premierx:dernierx, premiery:derniery] |
---|
| 104 | ; 'z':divi = e3w[premierz:dernierz] |
---|
| 105 | ; ELSE:return, report('Bad definition of direction argument') |
---|
| 106 | ; endcase |
---|
| 107 | ; END |
---|
| 108 | ELSE:return, report('Bad definition of vargrid') |
---|
| 109 | ENDCASE |
---|
| 110 | res = fitintobox(res) |
---|
| 111 | case 1 of |
---|
| 112 | ;---------------------------------------------------------------------------- |
---|
| 113 | ;xy |
---|
| 114 | ;---------------------------------------------------------------------------- |
---|
| 115 | taille[0] EQ 2:BEGIN |
---|
| 116 | earth = where(mask[*, *, premierz] EQ 0) |
---|
| 117 | if earth[0] NE -1 then res[earth] = !values.f_nan |
---|
| 118 | case direc of |
---|
| 119 | 'x':BEGIN |
---|
| 120 | res = (shift(res, -1, 0)-res)/divi |
---|
| 121 | if key_periodique EQ 0 OR nx NE jpi THEN res[nx-1, *] = !values.f_nan |
---|
| 122 | if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(res, 1, 0) |
---|
| 123 | END |
---|
| 124 | 'y':BEGIN |
---|
| 125 | res = (shift(res, 0, -1)-res)/divi |
---|
| 126 | res[*, ny-1] = !values.f_nan |
---|
| 127 | if vargrid EQ 'T' OR vargrid EQ 'U' then res = shift(res, 0, 1) |
---|
| 128 | END |
---|
| 129 | ELSE:return, report('Bad definition of direction argument for the type of array') |
---|
| 130 | ENDCASE |
---|
| 131 | earth = where(newmask[*, *, premierz] EQ 0) |
---|
| 132 | if earth[0] NE -1 then res[earth] = valmask |
---|
| 133 | END |
---|
| 134 | ;---------------------------------------------------------------------------- |
---|
| 135 | ;xyt |
---|
| 136 | ;---------------------------------------------------------------------------- |
---|
| 137 | taille[0] EQ 3 AND jpt NE 1:BEGIN |
---|
| 138 | earth = where(mask[*, *, premierz] EQ 0) |
---|
| 139 | if earth[0] NE -1 then BEGIN |
---|
| 140 | earth = earth#replicate(1, jpt)+replicate(1, n_elements(earth))#(nx*ny*lindgen(jpt)) |
---|
| 141 | res[earth[*]] = !values.f_nan |
---|
| 142 | ENDIF |
---|
| 143 | divi = divi[*]#replicate(1, jpt) |
---|
| 144 | case direc of |
---|
| 145 | 'x':BEGIN |
---|
| 146 | res = (shift(res, -1, 0, 0)-res)/divi |
---|
| 147 | if key_periodique EQ 0 OR nx NE jpi THEN res[nx-1, *, *] = !values.f_nan |
---|
| 148 | if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(res, 1, 0, 0) |
---|
| 149 | END |
---|
| 150 | 'y':BEGIN |
---|
| 151 | res = (shift(res, 0, -1, 0)-res)/divi |
---|
| 152 | res[*, ny-1, *] = !values.f_nan |
---|
| 153 | if vargrid EQ 'T' OR vargrid EQ 'U' then res = shift(res, 0, 1, 0) |
---|
| 154 | END |
---|
| 155 | ELSE:return, report('Bad definition of direction argument for the type of array') |
---|
| 156 | ENDCASE |
---|
| 157 | earth = where(newmask[*, *, premierz] EQ 0) |
---|
| 158 | if earth[0] NE -1 then BEGIN |
---|
| 159 | earth = earth#replicate(1, jpt)+replicate(1, n_elements(earth))#(nx*ny*lindgen(jpt)) |
---|
| 160 | res[earth] = valmask |
---|
| 161 | ENDIF |
---|
| 162 | END |
---|
| 163 | ;---------------------------------------------------------------------------- |
---|
| 164 | ;xyz |
---|
| 165 | ;---------------------------------------------------------------------------- |
---|
| 166 | taille[0] EQ 3 AND jpt EQ 1:BEGIN |
---|
| 167 | earth = where(mask EQ 0) |
---|
| 168 | if earth[0] NE -1 then res[earth] = !values.f_nan |
---|
| 169 | case direc OF |
---|
| 170 | 'x':BEGIN |
---|
| 171 | divi = divi[*]#replicate(1, nz) |
---|
| 172 | res = (shift(res, -1, 0, 0)-res)/divi |
---|
| 173 | if key_periodique EQ 0 OR nx NE jpi THEN res[nx-1, *, *] = !values.f_nan |
---|
| 174 | if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(res, 1, 0, 0) |
---|
| 175 | END |
---|
| 176 | 'y':BEGIN |
---|
| 177 | divi = divi[*]#replicate(1, nz) |
---|
| 178 | res = (shift(res, 0, -1, 0)-res)/divi |
---|
| 179 | res[*, ny-1, *] = !values.f_nan |
---|
| 180 | if vargrid EQ 'T' OR vargrid EQ 'U' then res = shift(res, 0, 1, 0) |
---|
| 181 | END |
---|
| 182 | 'z':BEGIN |
---|
| 183 | divi = reform(replicate(1, nx*ny)#divi, nx, ny, nz) |
---|
| 184 | if nx EQ 1 OR ny EQ 1 then res = reform(res, nx, ny, nz) |
---|
| 185 | if vargrid EQ 'W' THEN BEGIN |
---|
| 186 | res = (shift(res, 0, 0, 1)-res)/divi |
---|
| 187 | res[*, *, 0] = !values.f_nan |
---|
| 188 | ENDIF ELSE BEGIN |
---|
| 189 | res = (res-shift(res, 0, 0, -1))/divi |
---|
| 190 | res[*, *, nz-1] = !values.f_nan |
---|
| 191 | ENDELSE |
---|
| 192 | if earth[0] NE -1 then res[earth] = valmask |
---|
| 193 | END |
---|
| 194 | ENDCASE |
---|
| 195 | END |
---|
| 196 | ;------------------------------------------------------------ |
---|
| 197 | ;---------------------------------------------------------------------------- |
---|
| 198 | ;xyzt |
---|
| 199 | ;---------------------------------------------------------------------------- |
---|
| 200 | taille[0] EQ 4:BEGIN |
---|
| 201 | earth = where((mask[*]#replicate(1, jpt)) EQ 0) |
---|
| 202 | if earth[0] NE -1 then res[earth] = !values.f_nan |
---|
| 203 | case direc OF |
---|
| 204 | 'z':BEGIN |
---|
| 205 | divi = replicate(1, nx*ny)#divi |
---|
| 206 | divi = reform(divi[*]#replicate(1, jpt), nx, ny, nz, jpt, /over) |
---|
| 207 | if nx EQ 1 OR ny EQ 1 then res = reform(res, nx, ny, nz, jpt) |
---|
| 208 | if vargrid EQ 'W' THEN BEGIN |
---|
| 209 | res = (shift(res, 0, 0, 1, 0)-res)/divi |
---|
| 210 | res[*, *, 0, *] = !values.f_nan |
---|
| 211 | ENDIF ELSE BEGIN |
---|
| 212 | res = (res-shift(res, 0, 0, -1, 0))/divi |
---|
| 213 | res[*, *, nz-1, *] = !values.f_nan |
---|
| 214 | ENDELSE |
---|
| 215 | if earth[0] NE -1 then res[earth] = valmask |
---|
| 216 | END |
---|
| 217 | ENDCASE |
---|
| 218 | END |
---|
| 219 | ;------------------------------------------------------------ |
---|
| 220 | ;------------------------------------------------------------ |
---|
| 221 | endcase |
---|
| 222 | varname = 'grad of '+varname |
---|
| 223 | varunit = varunit+'/m' |
---|
| 224 | |
---|
| 225 | |
---|
| 226 | ;------------------------------------------------------------ |
---|
| 227 | return, res |
---|
| 228 | end |
---|
| 229 | |
---|
| 230 | |
---|
| 231 | |
---|
| 232 | |
---|
| 233 | |
---|
| 234 | |
---|