[2] | 1 | ;+ |
---|
[231] | 2 | ; |
---|
[150] | 3 | ; @file_comments |
---|
[163] | 4 | ; Build the triangulation for a E-grid type |
---|
[2] | 5 | ; |
---|
[150] | 6 | ; @categories |
---|
[157] | 7 | ; Graphics |
---|
[231] | 8 | ; |
---|
[325] | 9 | ; @param maskentree {in}{optional}{type=2d array} |
---|
| 10 | ; It is a 2d array which will serve to mask the field we will trace after |
---|
| 11 | ; with CONTOUR, ...TRIANGULATION=triangule(mask) |
---|
[150] | 12 | ; If this argument is not specified, the function use tmask |
---|
[231] | 13 | ; |
---|
[150] | 14 | ; @keyword BASIC |
---|
[325] | 15 | ; Specify that the mask is on a basic grid (use the triangulation for vertical |
---|
| 16 | ; cuts and hovmoellers) |
---|
[2] | 17 | ; |
---|
[231] | 18 | ; @keyword COINMONTE {type=array} |
---|
| 19 | ; To obtain the array of "ascending land corner" to be treated with |
---|
| 20 | ; <pro>completecointerre</pro> in the variable array instead of make it pass |
---|
| 21 | ; by the global variable twin_corners_up. |
---|
| 22 | ; |
---|
[163] | 23 | ; @keyword COINDESCEND {type=array} |
---|
| 24 | ; See COINMONTE |
---|
[2] | 25 | ; |
---|
[150] | 26 | ; @keyword SHIFTED |
---|
[231] | 27 | ; |
---|
[150] | 28 | ; @uses |
---|
[370] | 29 | ; <pro>common</pro> |
---|
[231] | 30 | ; |
---|
[150] | 31 | ; @history |
---|
[157] | 32 | ; Sebastien Masson (smasson\@lodyc.jussieu.fr) |
---|
[2] | 33 | ; june 2001 |
---|
[231] | 34 | ; |
---|
| 35 | ; @version |
---|
[150] | 36 | ; $Id$ |
---|
[231] | 37 | ; |
---|
[150] | 38 | ; @todo |
---|
[231] | 39 | ; seb L.152->153 je ne pense pas que ce soit ce que tu voulais dire mais |
---|
[150] | 40 | ; c'est la traduction de ce qu'il y avait écrit. Correction si besoin. |
---|
[325] | 41 | ; |
---|
[2] | 42 | ;- |
---|
[325] | 43 | FUNCTION triangule_e, maskentree $ |
---|
[327] | 44 | , COINMONTE=coinmonte, COINDESCEND=coindescend $ |
---|
| 45 | , SHIFTED=shifted, BASIC=basic |
---|
[231] | 46 | ; |
---|
[114] | 47 | compile_opt idl2, strictarrsubs |
---|
| 48 | ; |
---|
[29] | 49 | @cm_4mesh |
---|
| 50 | IF NOT keyword_set(key_forgetold) THEN BEGIN |
---|
| 51 | @updatenew |
---|
| 52 | ENDIF |
---|
| 53 | ;--------------------------------------------------------- |
---|
[150] | 54 | tempsun = systime(1) ; For key_performance |
---|
[2] | 55 | ;------------------------------------------------------------ |
---|
[150] | 56 | ; Is the mask given or do we have to take tmask? |
---|
[2] | 57 | ;------------------------------------------------------------ |
---|
| 58 | ; |
---|
| 59 | msk = maskentree |
---|
| 60 | sizem = size(msk) |
---|
| 61 | nx = sizem[1] |
---|
| 62 | ny = sizem[2] |
---|
| 63 | ;------------------------------------------------------------ |
---|
[29] | 64 | if keyword_set(key_periodic)*(nx EQ jpi) $ |
---|
[231] | 65 | AND NOT keyword_set(basic) then BEGIN |
---|
[2] | 66 | msk = [msk, msk[0, *]] |
---|
| 67 | nx = nx+1 |
---|
| 68 | ENDIF |
---|
| 69 | ; |
---|
| 70 | ; we will find the diamond that must be cut in two triangle using the |
---|
| 71 | ; horizontal diagonal. |
---|
| 72 | ; |
---|
| 73 | index = lindgen(nx, ny) |
---|
| 74 | index = index[0:nx-2, 1:ny-2] |
---|
| 75 | if n_elements(shifted) EQ 0 then shifted = 1 |
---|
| 76 | oddeven = (index/nx+1-shifted) MOD 2 |
---|
| 77 | msk1 = msk[index] |
---|
| 78 | msk2 = msk[index+1] |
---|
| 79 | sum = msk[index-nx+oddeven]+msk[index+nx+oddeven] |
---|
| 80 | sum1 = msk2+sum |
---|
| 81 | sum2 = msk1+sum |
---|
| 82 | ; |
---|
| 83 | ; horizontal |
---|
| 84 | ; |
---|
| 85 | singularpoint = where((msk1 EQ 0 AND sum1 EQ 3) OR (msk1 EQ 1 AND sum1 EQ 0) $ |
---|
| 86 | OR (msk2 EQ 0 AND sum2 EQ 3) OR (msk2 EQ 1 AND sum2 EQ 0) $ |
---|
| 87 | OR (sum EQ 0 AND (msk1+msk2) EQ 2) ) |
---|
| 88 | if singularpoint[0] NE -1 then begin |
---|
| 89 | horizontal = index[singularpoint] |
---|
| 90 | triang = definetri_e(nx, ny, horizontal, SHIFTED = shifted) |
---|
| 91 | ENDIF ELSE triang = definetri_e(nx, ny, SHIFTED = shifted) |
---|
| 92 | ; coinmont = index[where(sum EQ 2 AND (msk1+msk2) EQ 0)] |
---|
| 93 | ; coindesc = index[where(sum EQ 0 AND (msk1+msk2) EQ 2)] |
---|
| 94 | ; |
---|
| 95 | ; we keep only the triangles which are outside the land |
---|
| 96 | ; but for some reasons we will in fact delete the land diamond |
---|
| 97 | ; |
---|
| 98 | ; allrecinland = where(sum1+msk1 EQ 0) |
---|
| 99 | ; indexallinland = index[allrecinland] |
---|
| 100 | ; otherrec = (lindgen(nx, ny))[0:nx-2, 1:ny-2] |
---|
| 101 | ; otherrec = different(otherrec, indexallinland) |
---|
| 102 | ; ; |
---|
| 103 | ; index = lindgen(nx, ny) |
---|
| 104 | ; index = index[0:nx-3, 2:ny-3] |
---|
| 105 | ; out = inter(index, indexallinland) |
---|
| 106 | ; IF out[0] NE -1 THEN begin |
---|
| 107 | ; out = inter(out+1, indexallinland) |
---|
| 108 | ; IF out[0] NE -1 THEN begin |
---|
| 109 | ; out = out-1 |
---|
| 110 | ; oddeven = (out/nx+1-shifted) MOD 2 |
---|
| 111 | ; out = inter(out-nx+oddeven, otherrec) |
---|
| 112 | ; IF out[0] NE -1 THEN begin |
---|
| 113 | ; out = inter(out+2*nx, otherrec) |
---|
| 114 | ; IF out[0] NE -1 THEN begin |
---|
| 115 | ; out = out-(nx+((out/nx+shifted) MOD 2)) |
---|
| 116 | ; endif |
---|
| 117 | ; endif |
---|
| 118 | ; endif |
---|
| 119 | ; ENDIF |
---|
| 120 | ; help, out |
---|
| 121 | ; ; |
---|
| 122 | ; index = lindgen(nx, ny) |
---|
| 123 | ; index = index[0:nx-3, 2:ny-3] |
---|
| 124 | ; out = inter(index, otherrec) |
---|
| 125 | ; IF out[0] NE -1 THEN begin |
---|
| 126 | ; out = inter(out+1, otherrec) |
---|
| 127 | ; IF out[0] NE -1 THEN begin |
---|
| 128 | ; out = out-1 |
---|
| 129 | ; oddeven = (out/nx+1-shifted) MOD 2 |
---|
| 130 | ; out = inter(out-nx+oddeven, indexallinland) |
---|
| 131 | ; IF out[0] NE -1 THEN begin |
---|
| 132 | ; out = inter(out+2*nx, indexallinland) |
---|
| 133 | ; IF out[0] NE -1 THEN begin |
---|
| 134 | ; out = out-(nx+((out/nx+shifted) MOD 2)) |
---|
| 135 | ; endif |
---|
| 136 | ; endif |
---|
| 137 | ; endif |
---|
| 138 | ; endif |
---|
| 139 | ; help, out |
---|
| 140 | ; ; |
---|
| 141 | ; IF out[0] EQ -1 THEN out = different(indexallinland, out) ELSE out = indexallinland |
---|
| 142 | ; triout = numtri(out, nx, ny) |
---|
| 143 | ; triout = [triout, triout+1] |
---|
| 144 | ; goodtri = lindgen(2*(nx-1)*(ny-1)) |
---|
| 145 | ; goodtri = different(goodtri, triout) |
---|
| 146 | ; triang = triang[*, temporary(goodtri)] |
---|
| 147 | |
---|
| 148 | ; ; |
---|
| 149 | ;------------------------------------------------------------ |
---|
[231] | 150 | ; When key_periodic equal 1, triang is a list of indexes's array which |
---|
[150] | 151 | ; have a surplus column. |
---|
[231] | 152 | ; We have to put it back to the initial matrix by putting indexes of |
---|
[150] | 153 | ; the last column equal to these of the last column... |
---|
[2] | 154 | ;------------------------------------------------------------ |
---|
[150] | 155 | tempdeux = systime(1) ; For key_performance =2 |
---|
[29] | 156 | if keyword_set(key_periodic)*(nx-1 EQ jpi) $ |
---|
[231] | 157 | AND NOT keyword_set(basic) then BEGIN |
---|
[2] | 158 | indicey = triang/nx |
---|
| 159 | indicex = triang-indicey*nx |
---|
| 160 | nx = nx-1 |
---|
| 161 | liste = where(indicex EQ nx) |
---|
| 162 | if liste[0] NE -1 then indicex[liste] = 0 |
---|
| 163 | triang = indicex+nx*indicey |
---|
| 164 | nx = nx+1 |
---|
| 165 | ; if coinmont[0] NE -1 then begin |
---|
| 166 | ; indicey = coinmont/nx |
---|
| 167 | ; indicex = coinmont-indicey*nx |
---|
| 168 | ; nx = nx-1 |
---|
| 169 | ; liste = where(indicex EQ nx) |
---|
| 170 | ; if liste[0] NE -1 THEN indicex[liste] = 0 |
---|
| 171 | ; coinmont = indicex+nx*indicey |
---|
| 172 | ; nx = nx+1 |
---|
| 173 | ; endif |
---|
| 174 | ; if coindesc[0] NE -1 then begin |
---|
| 175 | ; indicey = coindesc/nx |
---|
| 176 | ; indicex = coindesc-indicey*nx |
---|
| 177 | ; nx = nx-1 |
---|
| 178 | ; liste = where(indicex EQ nx) |
---|
| 179 | ; if liste[0] NE -1 THEN indicex[liste] = 0 |
---|
| 180 | ; coindesc = indicex+nx*indicey |
---|
| 181 | ; nx = nx+1 |
---|
| 182 | ; endif |
---|
| 183 | endif |
---|
| 184 | IF testvar(var = key_performance) EQ 2 THEN $ |
---|
| 185 | print, 'temps triangule: finitions', systime(1)-tempdeux |
---|
| 186 | |
---|
| 187 | ;------------------------------------------------------------ |
---|
[150] | 188 | if arg_present(coinmonte) THEN coinmonte = coinmont ELSE twin_corners_up = coinmont |
---|
| 189 | if arg_present(coindescend) THEN coindescend = coindesc ELSE twin_corners_dn = coindesc |
---|
| 190 | |
---|
| 191 | IF NOT keyword_set(key_forgetold) THEN BEGIN |
---|
| 192 | @updateold |
---|
| 193 | ENDIF |
---|
[2] | 194 | ;------------------------------------------------------------ |
---|
| 195 | |
---|
| 196 | |
---|
[231] | 197 | IF keyword_set(key_performance) THEN print, 'temps triangule', systime(1)-tempsun |
---|
[2] | 198 | |
---|
| 199 | return, triang |
---|
| 200 | |
---|
[231] | 201 | END |
---|