Ignore:
Timestamp:
08/09/06 12:12:54 (18 years ago)
Author:
navarro
Message:

english and nicer header (3a)

Location:
trunk/SRC/ToBeReviewed/TRIANGULATION
Files:
15 edited

Legend:

Unmodified
Added
Removed
  • trunk/SRC/ToBeReviewed/TRIANGULATION/ciseauxtri.pro

    r134 r150  
    33;------------------------------------------------------------ 
    44;+ 
    5 ; NAME: 
    6 ; 
    7 ; PURPOSE:vire les tableaux qui ne doivent pas etre dessines grace a 2 
    8 ; tests: 
    9 ;     1) les coins du tableau doivent etre ds la fenetre 
    10 ;     2) les clongeurs des cotes des triangfles exprimes en 
    11 ;     coordonnees normalisesne doivent pas depasser une certaine 
    12 ;     longueur seuil.  
    13 ; 
    14 ; CATEGORY: 
    15 ; 
    16 ; CALLING SEQUENCE: 
    17 ; 
    18 ; INPUTS: 
    19 ; 
    20 ; KEYWORD PARAMETERS: 
    21 ; 
    22 ; OUTPUTS: 
    23 ; 
    24 ; COMMON BLOCKS: 
    25 ;       common.pro 
    26 ; 
    27 ; SIDE EFFECTS: 
    28 ; 
    29 ; RESTRICTIONS: 
    30 ; 
    31 ; EXAMPLE: 
    32 ; 
    33 ; MODIFICATION HISTORY: Sebastien Masson (smasson@lodyc.jussieu.fr) 
     5; @file_comments  
     6; Delete arrays which do not have to be drawn thanks to 2 tests: 
     7;     1) Corners of the array must be in the window 
     8;     2) Lenghtes of side of triangles expressed in normalized  
     9;        coordinates must not surpass a sill lenght.  
     10; 
     11; @categories 
     12; 
     13; @param TRIANG 
     14;  
     15; @param GLAM  
     16; 
     17; @param GPHI  
     18; 
     19; @keyword ALL 
     20;  
     21; @keyword _EXTRA 
     22; Used to pass your keywords. 
     23; 
     24; @uses 
     25; common.pro 
     26; 
     27; @history 
     28; Sebastien Masson (smasson@lodyc.jussieu.fr) 
    3429;                       20/2/99 
     30; 
     31; @version 
     32; $Id$ 
     33;  
    3534;- 
    3635;------------------------------------------------------------ 
    3736;------------------------------------------------------------ 
    3837;------------------------------------------------------------ 
    39 FUNCTION ciseauxtri, triang, glam, gphi, TOUT = tout, _EXTRA = ex 
     38FUNCTION ciseauxtri, triang, glam, gphi, ALL = all, _EXTRA = ex 
    4039;--------------------------------------------------------- 
    4140; 
     
    5251             OR !map.projection EQ 18) THEN return, triang 
    5352; 
    54    tempsun = systime(1)         ; pour key_performance 
     53   tempsun = systime(1)         ; For key_performance 
    5554; 
    5655   taille = size(glam) 
     
    5857   ny = taille[2] 
    5958 
    60    tempdeux = systime(1)        ; pour key_performance =2 
     59   tempdeux = systime(1)        ; For key_performance =2 
    6160   z = convert_coord(glam[*],gphi[*],/data,/to_normal)  
    6261   x = z[0, *] 
     
    6665    print, 'temps ciseauxtri: convert_coord data to_normal', systime(1)-tempdeux 
    6766; 
    68 ; attention, suivant la projection certains points x ou y peuvent 
    69 ; devenir NaN (cf. points deriere la terre ds une projection 
    70 ; orthographique) il faut dans ce cas enlever tous les triangles qui 
    71 ; contiennent un de ces points 
     67; Beware, following the projection, some points x or y can become NaN  
     68; (see points behind the Earth in an orthographic projection).  
     69; In this case, we have to remove all triangle which contain one of these points. 
    7270; 
    7371   if (!map.projection LE 7 AND !map.projection NE 0) $ 
    7472    OR !map.projection EQ 14 OR !map.projection EQ 15 OR !map.projection EQ 18 then begin 
    75       tempdeux = systime(1)     ; pour key_performance =2 
     73      tempdeux = systime(1)     ; For key_performance =2 
    7674; 
    7775      test = (x*y)[triang] 
     
    8886   seuil = 5 < (min([nx, ny])-2) 
    8987; 
    90 ; maintenant on vire les triangles dont un des cotes a une taille 
    91 ; superieure a 1/seuil du domaine reserve au dessin. 
     88; Now we delete tringles whose one side has a size superior to 1/sill from the domain reserved for the drawing. 
    9289; 
    9390 
    9491   if keyword_set(key_periodic)  then begin 
    95       tempdeux = systime(1)     ; pour key_performance =2 
     92      tempdeux = systime(1)     ; For key_performance =2 
    9693; 
    9794      xtri = x[triang] 
     
    111108       print, 'temps ciseauxtri: trouver les triangles trop grands', systime(1)-tempdeux 
    112109; 
    113       tempdeux = systime(1)     ; pour key_performance =2 
     110      tempdeux = systime(1)     ; For key_performance =2 
    114111      if ind[0] NE -1 then triang = triang[*, temporary(ind)] ELSE return,  -1 
    115112      trichanged = 1b 
     
    118115   endif 
    119116; 
    120 ; on supprime les triangles qui sont un peu trop hors du cadre 
    121 ; valable qd /TOUT est active 
     117; We delete all triangle which are out of the valid frame when  on supprime les triangles qui sont un peu trop hors du cadre 
     118; valable qd ALL is activated 
    122119; 
    123120    if keyword_set(key_irregular) then begin 
    124121 
    125        tempdeux = systime(1)     ; pour key_performance =2 
     122       tempdeux = systime(1)     ; For key_performance =2 
    126123       xtri = x[triang] 
    127124       test1 = xtri GE !p.position[0] 
     
    148145    ENDIF 
    149146; 
    150 ;        dernier tri 
     147;        Last sort. 
    151148; 
    152149   if keyword_set(trichanged) then BEGIN 
    153150; 
    154 ; il faut verifier que les triangles que l''on garde ne 
    155 ; forment pas une triangulation dans laquelle 2 triangles ont un 
    156 ; sommet en commun sans avoir de cote de commun. 
    157 ; 
    158 ; on constitue la liste des rectangles que l''on souhaite garder (on 
    159 ; garde un rectangle des qu''il y a un triangle dedans) 
    160 ; 
    161 ; dans definetri, on a construit les triangles tels que le premier et  
    162 ; le dernier sommets soient ceux de la diagonale du rectangle definit 
    163 ; par le maillage. 
    164 ; pour retouver de quel rectangle provient un triangle, on cherche le 
    165 ; min de l''indice suivant x et suivant y de chaque triangle. Apres on 
    166 ; repasse cette indissage suivant x et y en indicage suivant nx*ny 
    167 ; 
    168       tempdeux = systime(1)     ; pour key_performance =2 
    169 ; indices en y des rectangles 
     151; We have to check that triangles we keep do not form a triangulation  
     152; in which 2 triangles have a common summit without have a common side.  
     153; 
     154; We constitute the list of rectangles we want to keep (we keep every  
     155; rectangle containing a triangle) 
     156; 
     157; In definetri, we have construct triangles just so the first and the 
     158; last summit are those of the diagonale of the rectangle defined by 
     159; the mesh size. To find from which rectangle a triangle comes, we look 
     160; for the min of the index folowing x and following y of each triangle. 
     161; Then we go by again this indexion following x and y in an indexion 
     162; following nx*ny/ 
     163; 
     164      tempdeux = systime(1)     ; For key_performance =2 
     165; y indexes of rectangles 
    170166      indytriang = (triang[0, *]/nx) < (triang[2, *]/nx) 
    171 ; indices en x des rectangles 
     167; x indexes of rectangles 
    172168      indxtriang0 = triang[0, *]-nx*(triang[0, *]/nx) 
    173169      indxtriang2 = triang[2, *]-nx*(triang[2, *]/nx) 
    174170      indxmin = indxtriang0 < indxtriang2 
    175 ; attention dans le cas ou la grille est periodique et ou on a tous 
    176 ; les points suivant x, les triangles qui assurent la periodicite en x 
    177 ; ont des indices suivant x qui sont 0 et nx-1. Ils appatienent au 
    178 ; rectangles de la colonne nx-1 et non de la colonne 0 
     171; Beware in the case where the grid is periodic and where we have all points  
     172; following x, triangles which assure the periodicity in x have indexes  
     173; following x which are 0 and nx-1. They belong to rectangles of the column  
     174; nx-1 and not to column 0. 
    179175      if keyword_set(key_periodic) AND nx EQ jpi then BEGIN 
    180176      indxmax = indxtriang0 > indxtriang2 
     
    185181       print, 'temps ciseauxtri: liste des rectangles', systime(1)-tempdeux 
    186182; 
    187 ; maintenant qu''on a cette liste on va s''assuter que l''on a pas de 
    188 ; triangles qui n''ont qu''on sommet en commun. 
     183; Now we have this list, we will make sure that we do not have triangles  
     184; with only a common summit. 
    189185; 
    190186      test = bytarr(nx, ny) 
     
    192188      dejavire = 1b-test 
    193189; 
    194       tempdeux = systime(1)     ; pour key_performance =2 
     190      tempdeux = systime(1)     ; For key_performance =2 
    195191      vire1 = 0 
    196192      vire2 = 0 
     
    198194         vire1 = where( (test*shift(test, -1, -1) $ 
    199195                         *(1-shift(test, 0, -1))*(1-shift(test, -1, 0))) EQ 1) 
    200          if vire1[0] NE -1 THEN test[vire1] = 0 ; on vire le rectangle  
     196         if vire1[0] NE -1 THEN test[vire1] = 0 ; We delete the rectangle 
    201197         vire2 = where( ((1-test)*(1-shift(test, -1, -1)) $ 
    202198                         *shift(test, 0, -1)*shift(test, -1, 0)) EQ 1) 
    203 ; on vire le rectangle du dessus (meme indice x mais egale a 1) 
     199; We delete the top rectangle (same x index but equal to 1) 
    204200         if vire2[0] NE -1 THEN test[vire2+nx] = 0  
    205201      ENDWHILE 
     
    212208; 
    213209      if avirer[0] NE -1 then begin 
    214       tempdeux = systime(1)     ; pour key_performance =2 
     210      tempdeux = systime(1)     ; For key_performance =2 
    215211      indnx = n_elements(listrect) 
    216212      indny = n_elements(avirer) 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/completecointerre.pro

    r134 r150  
    1 ;------------------------------------------------------------ 
    2 ;------------------------------------------------------------ 
    3 ;------------------------------------------------------------ 
    4 ;+ 
    5 ; NAME: COMPLETECOINTERRE 
    6 ; 
    7 ; PURPOSE: pour colorier proprement les continents! (c''est une longue 
    8 ; histoire...) 
    9 ; 
    10 ; CATEGORY: pour plt 
    11 ; 
    12 ; CALLING SEQUENCE: completecointerre 
    13 ;  
    14 ; INPUTS: non 
    15 ; 
    16 ; KEYWORD PARAMETERS:  _EXTRA  
    17 ; 
    18 ;        CONT_COLOR: the color of the continent. defaut value is 
    19 ;        (!d.n_colors - 1) < 255 => white 
    20 ; 
    21 ; OUTPUTS: non 
    22 ; 
    23 ; COMMON BLOCKS: common.pro 
    24 ; 
    25 ; SIDE EFFECTS: 
    26 ; 
    27 ; RESTRICTIONS: 
    28 ; 
    29 ; EXAMPLE: 
    30 ; 
    31 ; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr) 
    32 ;                      01/10/1999 
    33 ;- 
    34 ;------------------------------------------------------------ 
    35 ;------------------------------------------------------------ 
    36 ;------------------------------------------------------------ 
    371PRO draw_corner_triangle, lons, lats, seuil, CONT_COLOR = cont_color, _extra = ex 
    382; 
     
    5519;------------------------------------------------------------ 
    5620;------------------------------------------------------------ 
     21;------------------------------------------------------------ 
     22;+ 
     23; 
     24; @file_comments 
     25; To color cleanly continents 
     26; 
     27; @categories 
     28; graphic 
     29; 
     30; @keyword _EXTRA 
     31; Used to pass your keywords  
     32; 
     33; @keyword CONT_COLOR 
     34; The color of the continent. defaut value is 
     35; (!d.n_colors - 1) < 255 => white 
     36; 
     37; @uses 
     38; common.pro 
     39; 
     40; @history 
     41; Sebastien Masson (smasson@lodyc.jussieu.fr) 
     42;                      01/10/1999 
     43; 
     44; @version 
     45; $Id$ 
     46; 
     47;- 
     48;------------------------------------------------------------ 
     49;------------------------------------------------------------ 
     50;------------------------------------------------------------ 
    5751PRO completecointerre, COINMONTE = coinmonte, COINDESCEND = coindescend $ 
    5852                       , CONT_COLOR = cont_color, INDICEZOOM = indicezoom $ 
     
    6660;   if NOT keyword_set(coindescend) then return 
    6761;   if NOT keyword_set(indicezoom) then return 
    68   tempsun = systime(1)          ; pour key_performance 
     62  tempsun = systime(1)          ; For key_performance 
    6963;------------------------------------------------------------ 
    70 ; definitions des vecteurs coinmont et coindesc 
     64; definitions of vectors coinmont and coindesc 
    7165;------------------------------------------------------------ 
    7266  if keyword_set(coinmonte) then coinmont = coinmonte $ 
     
    7670  IF NOT keyword_set(cont_color) THEN cont_color = (!d.n_colors-1) <  255 
    7771;------------------------------------------------------------ 
    78 ; definition descoordonnees des points numerotes 1,2,3,4,5,6 cf. les 
    79 ; schemas en dessous! 
     72; definition of coordinates of points numbered 1,2,3,4,5,6 (see figures below) 
    8073;------------------------------------------------------------ 
    81   tempdeux = systime(1)         ; pour key_performance =2 
     74  tempdeux = systime(1)         ; For key_performance =2 
    8275  if coinmont[0] NE -1 OR coindesc[0] NE -1 then BEGIN 
    8376    if keyword_set(indicezoom) then BEGIN 
     
    122115; 
    123116; 
    124 ; cas coin terre en montee: 
    125 ;      2 points terre en diagonale montante avec 2 points mer sur 
    126 ;      la diagonale descendante. 
     117; Case land corner in ascent: 
     118;      2 land points in diagonal ascending with 2 ocean points on the descendant diagonal. 
    127119; 
    128120;                     4     
     
    139131; 
    140132  if coinmont[0] NE -1 then BEGIN 
    141     tempdeux = systime(1)       ; pour key_performance =2 
     133    tempdeux = systime(1)       ; For key_performance =2 
    142134    for id = 0, n_elements(coinmont)-1 do BEGIN 
    143135      i = coinmont[id] 
     
    159151  ENDIF 
    160152;------------------------------------------------------------ 
    161 ; cas coin terre en descendante.: 
    162 ;      2 points terre en diagonale descendante avec 2 points mer sur 
    163 ;      la diagonale montante 
    164 ; 
     153; Case land corner in descent: 
     154;      2 land points in diagonal descending with 2 ocean points on the ascendant diagonal. 
     155;  
    165156;                     4 
    166157;     t(i+nx)=1    u(i+nx)       t(i+nx+1)=0 
     
    175166; 
    176167  if coindesc[0] NE -1 then begin 
    177     tempdeux = systime(1)       ; pour key_performance =2 
     168    tempdeux = systime(1)       ; For key_performance =2 
    178169    for id = 0, n_elements(coindesc)-1 do BEGIN 
    179170      i = coindesc[id] 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/definetri.pro

    r134 r150  
    11;+ 
    2 ; NAME:definetri 
    32; 
    4 ; PURPOSE:Define a triangulation array like TRIANGULATE.  
     3; @file_comments 
     4; efine a triangulation array like TRIANGULATE.  
    55;         But in a VERY SIMPLE CASE: 
    66; the points are regulary-gridded on nx*ny array. 
     
    3535; 
    3636; 
    37 ; CATEGORY: to understand how TRIANGULATE and TRIANGULATION work! 
     37; @categories 
     38; utilities 
     39;  
     40; @param NX {in}{required} 
     41; The x dimension array 
    3842; 
    39 ; CALLING SEQUENCE:triangles=definetri(nx, ny [,downward]) 
    40 ;  
    41 ; INPUTS: nx and ny are the array dimensions  
     43; @param NY {in}{required} 
     44; The y dimension array 
    4245; 
    43 ; OPTIONAL INPUTS: 
     46; @param DOWNWARD {in}{optional} 
     47; When downward is undefine all rectangles are cut 
     48; in using the upward diagonal. Downward is a vector which 
     49; contains the rectangles numbers which are cut in using the 
     50; downward diagonal. 
     51; The rectangle number is define by the index (in a nx*ny 
     52; vector) of the lower-left corner of the rectangle. 
    4453; 
    45 ;        downward: When downward is undefine all rectangles are cut 
    46 ;        in using the upward diagonal. Downward is a vector which 
    47 ;        contains the rectangles numbers which are cut in using the 
    48 ;        downward diagonal. 
    49 ;        The rectangle number is define by the index (in a nx*ny 
    50 ;        vector) of the lower-left corner of the rectangle. 
     54; @returns 
     55; triangles is a 2d array and is dimensions are 3 and 
     56; 2*(nx-1)*(ny-1) 
     57; triangles is define like in the TRIANGULATE procedure. 
    5158; 
    52 ; KEYWORD PARAMETERS: 
    53 ; 
    54 ; OUTPUTS: 
    55 ;        triangles is a 2d array and is dimensions are 3 and 
    56 ;        2*(nx-1)*(ny-1) 
    57 ;        triangles is define like in the TRIANGULATE procedure. 
    58 ;         
    59 ; 
    60 ; OPTIONAL OUTPUTS: 
    61 ; 
    62 ; COMMON BLOCKS: 
    63 ; 
    64 ; SIDE EFFECTS: 
    65 ; 
    66 ; RESTRICTIONS: 
    67 ; 
    68 ; PROCEDURE: 
    69 ; 
    70 ; EXAMPLE: 
     59; @examples 
    7160; 
    7261; triangles=definetri(3,3,[1,3]) 
     
    8473; 
    8574; 
    86 ; MODIFICATION HISTORY: sebastien Masson (smlod@ipsl.jussieu.fr) 
     75; @history 
     76; sebastien Masson (smlod@ipsl.jussieu.fr) 
    8777;                       4/3/1999 
    8878; 
     79; @version 
     80; $Id$ 
    8981;- 
    9082FUNCTION definetri, nx, ny, downward 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/definetri_e.pro

    r134 r150  
    1313 
    1414;+ 
    15 ; NAME:definetri 
    1615; 
    17 ; PURPOSE:Define a triangulation array like TRIANGULATE but for a 
     16; @file_comments 
     17; Define a triangulation array like TRIANGULATE but for a 
    1818; E-grid type 
    1919; 
    20 ; CATEGORY: make contours with E-grid type 
     20; @categories 
     21; Make contours with E-grid type 
     22;  
     23; @param NX {in}{required} 
     24; The x dimension array 
    2125; 
    22 ; CALLING SEQUENCE:triangles=definetri(nx, ny [,vertical]) 
    23 ;  
    24 ; INPUTS: nx and ny are the array dimensions  
     26; @param NY {in}{required} 
     27; The y dimension array 
    2528; 
    26 ; OPTIONAL INPUTS: 
     29; @param SINGULAR {in}{optional} 
     30; When singular is undefined all rectangles are cut 
     31; in using the vertical diagonal. Singular is a vector which 
     32; contains the rectangles numbers which are cut in using the 
     33; horizontal diagonal. 
     34; The rectangle number is define by the index (in a nx*ny 
     35; vector) of the lower-left corner of the rectangle. 
    2736; 
    28 ;        vertical: When vertical is undefine all rectangles are cut 
    29 ;        in using the horizontal diagonal. Vertical is a vector which 
    30 ;        contains the rectangles numbers which are cut in using the 
    31 ;        vertical diagonal. 
    32 ;        The rectangle number is define by the index (in a nx*ny 
    33 ;        vector) of the lower-left corner of the rectangle. 
    34 ; 
    35 ; KEYWORD PARAMETERS: 
    36 ; 
    37 ; OUTPUTS: 
    38 ;        triangles is a 2d array and is dimensions are 3 and 
    39 ;        2*(nx-1)*(ny-1) 
    40 ;        triangles is define like in the TRIANGULATE procedure. 
    41 ;         
    42 ; 
    43 ; OPTIONAL OUTPUTS: 
    44 ; 
    45 ; COMMON BLOCKS: 
    46 ; 
    47 ; SIDE EFFECTS: 
    48 ; 
    49 ; RESTRICTIONS: 
    50 ; 
    51 ; PROCEDURE: 
    52 ; 
    53 ; EXAMPLE: 
     37; @keyword SHIFTED 
    5438; 
    5539; 
    56 ; MODIFICATION HISTORY: sebastien Masson (smlod@ipsl.jussieu.fr) 
     40; @returns 
     41; Triangles is a 2d array and is dimensions are 3 and 
     42; (nx-1)*(ny-1) 
     43; Triangles is defined like in the TRIANGULATE procedure. 
     44;         
     45; @history 
     46; Sebastien Masson (smlod@ipsl.jussieu.fr) 
    5747;                       June 2001 
    5848; 
     49; @version 
     50; $Id$ 
     51; 
     52; @todo  
     53; seb: documenter SHIFTED 
    5954;- 
    6055FUNCTION definetri_e, nx, ny, singular, SHIFTED = shifted 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/dessinetri.pro

    r134 r150  
    33;------------------------------------------------------------ 
    44;+ 
    5 ; NAME:dessinetri 
    65; 
    7 ; PURPOSE:dessine la triangulation 
     6; @file_comments 
     7; Draw the triangulation 
    88; 
    9 ; CATEGORY:pour comprendre comment ca marche 
     9; @categories 
     10; utilities 
     11;  
     12; @param TRI {in}{optional} 
     13; Array which define the triangulation (provided by triangule.pro or triangulate) 
     14;  
     15; @param X {in}{optional} 
     16; The x position of points to which the trangulation refer to  
     17; (see the x array provided in triangulate) 
    1018; 
    11 ; CALLING SEQUENCE:dessinetri [, tri, x, y] 
    12 ;  
    13 ; INPUTS:optionnels 
    14 ;        par defaut on choisit la triangulation qui est utilise pour 
    15 ;        les plots et on la trace aux points definites par vargrid 
     19; @param Y {in}{optional} 
     20; The y position of points to which the trangulation refer to  
     21; (see the y array provided in triangulate) 
    1622; 
    17 ;        sinon il faut fournir les tableaux  
    18 ;        tri definissant la triangulation (fournis par triangule.pro 
    19 ;        ou triangulate)  
    20 ;        x et y qui sont les positions de points a laquelle se raporte 
    21 ;        la triangulation (cf. les tableau x et y fournis ds 
    22 ;        triangulate) 
     23; @keyword WAIT 
     24; =x. to call wait x second between each triangle draw. 
    2325; 
    24 ; KEYWORD PARAMETERS:  
     26; @keyword ONEBYONE 
     27; To draw the triangles one by one 
    2528; 
    26 ;        All plots or polyfill keywords. 
     29; @keyword FILL 
     30; To fill the triangles (using polyfill) instead of plotting them 
    2731; 
    28 ;        WAIT=x. to call wait x second between each triangle draw. 
     32; @keyword CHANGECOLOR 
     33; =n. To change the color of each traingle. n colors 
     34; will be used and repeted if necessary. 
    2935; 
    30 ;        /ONEBYONE: to draw the triangles one by one 
     36; @uses 
     37; common.pro 
    3138; 
    32 ;        /FILL: to fill the triangles (using polyfill) instead of plotting them 
     39; @history 
     40; Sebastien Masson (smasson@lodyc.jussieu.fr) 
    3341; 
    34 ;        CHANGECOLOR=n. to change the color of each traingle. n colors 
    35 ;        will be used and repeted if necessary. 
    36 ; 
    37 ; 
    38 ; OUTPUTS: 
    39 ; 
    40 ; COMMON BLOCKS:common.pro 
    41 ; 
    42 ; SIDE EFFECTS: 
    43 ; 
    44 ; RESTRICTIONS: 
    45 ; 
    46 ; EXAMPLE: 
    47 ; 
    48 ; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr) 
     42; @version 
     43; $Id$ 
    4944; 
    5045;- 
     
    5853; 
    5954@common 
    60    tempsun = systime(1)         ; pour key_performance 
     55   tempsun = systime(1)         ;  For key_performance 
    6156   a = '' 
    6257   if n_params() EQ 3 then BEGIN 
     
    9994     color = color#replicate(1, n_elements(tri)/3/n_elements(color)+1)  
    10095; 
    101    tempdeux = systime(1)         ; pour key_performance =2 
     96   tempdeux = systime(1)         ; For key_performance =2 
    10297   for i = 0L, n_elements(tri)/3-1 do begin 
    10398      t = [tri[*, i], tri[0, i]] 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/drawcoast_c.pro

    r134 r150  
     1;+ 
     2; @file_comments 
     3; 
     4; 
     5; @categories 
     6; 
     7;  
     8; @param MASK 
     9; 
     10;  
     11; @param XF 
     12;  
     13; 
     14; @param YF 
     15;  
     16; 
     17; @param NX 
     18;  
     19; 
     20; @param NY 
     21;  
     22; 
     23; @keyword COAST_COLOR 
     24; The color of the coastline. 
     25; defaut value is 0 => black 
     26; 
     27; @keyword COAST_THICK 
     28; The thick of the trait to trace continents 
     29; By default, it is 1. 
     30; 
     31; @keyword XSEUIL 
     32; To eliminate segments of coasts which are to big (which link points which can  
     33; be close on the sphere but distant on the drawing). We delete all segments  
     34; whose the size surpass the size of the window following X/XSEUIL. By default,  
     35; XSEUIL=5. But it can be to big if we do a big zoom or a little one for some  
     36; projections... We specify it the keyword thanks to this keyword! 
     37; 
     38; @keyword YSEUIL 
     39; See XSEUIL 
     40;  
     41; @keyword _EXTRA 
     42; Used to pass our keywords 
     43 
     44; @returns 
     45;  
     46;  
     47; @uses 
     48;  
     49;  
     50; @restrictions 
     51;  
     52;  
     53; @examples 
     54;  
     55;  
     56; @history 
     57;  
     58;  
     59; @version  
     60;  
     61; @todo  
     62; Seb: remplir le header 
     63; 
     64;- 
     65;--------------------------------------------------------- 
    166PRO drawcoast_c, mask, xf, yf, nx, ny, COAST_COLOR = coast_color, COAST_THICK = coast_thick, YSEUIL = yseuil, XSEUIL = xseuil, _extra = ex 
    2 ;--------------------------------------------------------- 
    367; 
    468  compile_opt idl2, strictarrsubs 
     
    1074  ENDIF 
    1175;--------------------------------------------------------- 
    12    tempsun = systime(1)         ; pour key_performance 
     76   tempsun = systime(1)         ; For key_performance 
    1377;--------------------------------------------------------- 
    1478; 
    15 ; on trace les segments verticaux: 
     79; We trace vertical segments: 
    1680; 
    1781   if NOT keyword_set(yseuil) then yseuil = 5. < (min([nx, ny])-2) 
    1882   distanceseuil = (!p.position[3]-!p.position[1])/yseuil 
    19 ; liste: liste des points i pourlesquels on va tracer un segment entre 
    20 ; le point i,j-1 et i,j 
    21    tempdeux = systime(1)        ; pour key_performance =2 
     83; list: list of points i for which we will trace a segment between the point i,j-1 and i,j 
     84   tempdeux = systime(1)        ; For key_performance =2 
    2285   liste = where((mask+shift(mask, -1, 0)) EQ 1 $ 
    2386                 AND ((xf-shift(xf, 0, 1))^2+(yf-shift(yf, 0, 1))^2) LE distanceseuil^2) 
    2487   IF liste[0] NE -1 THEN BEGIN 
    25 ; on recupere lx et ly qui sont les indices ds un tableau 2d des 
    26 ; points donnes par liste 
     88; We recuperate lx an dly which are indexes in a 2d array of points given by list 
    2789      ly = liste/nx & lx = temporary(liste)-nx*ly 
    28       indice = where(ly NE 0)   ; on ne prend pas les points concernant 
     90      indice = where(ly NE 0)   ; We do not take points concerning 
    2991      if indice[0] NE -1 then begin 
    30 ; la premiere ligne car ds ce cas le pt j-1 n''est pas definit 
     92; the first line because in this case, the point j-1 is undefined. 
    3193         lx = lx[indice] & ly = ly[temporary(indice)] 
    32 ; boucle sur les points concernes et trace du segment 
    33 ; rq: on utilise plost au lieu de plot car plots est bcp plus rapide. 
     94; Loop on concerned points and drawing of the segment. 
     95; Comment: we use plots instead of plot because plots goes faster. 
    3496         IF testvar(var = key_performance) EQ 2 THEN $ 
    3597          print, 'temps tracecote: determiner liste des points concernes par un trait vertical', systime(1)-tempdeux 
    36          tempdeux = systime(1)  ; pour key_performance =2 
     98         tempdeux = systime(1)  ; For key_performance =2 
    3799         for pt = 0L, n_elements(lx)-1 do BEGIN  
    38100            i = lx[pt] & j = ly[pt] 
     
    45107   ENDIF 
    46108; 
    47 ; pour le trace des segments horizontaux, c''est la meme chose sauf 
    48 ; qu'il faut faire attention si on est periodique:  
     109; For the drawing of horizontal segments , it is the same thing but we have to be careful if it is periodic. 
    49110; 
    50 ; si on est periodique on duplique la premiere colonne et on la met a 
    51 ; la fin. (ceci est fait non pas pour le shift qui est par defaut 
    52 ; periodique mais pour le plots  
    53    tempdeux = systime(1)        ; pour key_performance =2 
     111; If it is periodic, we duplicate the first column and we put it at the end.  
     112; (This is made not for the shift which is periodic by default but for the plots) 
     113   tempdeux = systime(1)        ; For key_performance =2 
    54114   if keyword_set(key_periodic) AND nx EQ jpi then begin 
    55115      mask = [mask, mask[0, *]] 
     
    66126      indice = where(ly NE ny-1 AND lx NE 0) 
    67127      if indice[0] NE -1 then begin 
    68 ; on ne prend pas les points de la  
    69 ; premiere colonne et de la derniere ligne (car on l''a rajoute artificiellement!)) 
     128; We do not take points of the first column and of the last line (because we added artificially) 
    70129         lx = lx[indice] & ly = ly[temporary(indice)] 
    71130         IF testvar(var = key_performance) EQ 2 THEN $ 
    72131          print, 'temps tracecote: determiner liste des points concernes par un trait horizontal', systime(1)-tempdeux 
    73          tempdeux = systime(1)  ; pour key_performance =2 
     132         tempdeux = systime(1)  ; For key_performance =2 
    74133         for pt = 0L, n_elements(lx)-1 do BEGIN  
    75134            i = lx[pt] & j = ly[pt] 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/drawcoast_e.pro

    r134 r150  
     1;+ 
     2; @file_comments 
     3; 
     4; 
     5; @categories 
     6; 
     7;  
     8; @param MASK 
     9; 
     10;  
     11; @param XF 
     12;  
     13; 
     14; @param YF 
     15;  
     16; 
     17; @param NX 
     18;  
     19; 
     20; @param NY 
     21;  
     22; 
     23; @keyword COAST_COLOR 
     24; The color of the coastline. 
     25; defaut value is 0 => black 
     26; 
     27; @keyword COAST_THICK 
     28; The thick of the trait to trace continents 
     29; By default, it is 1. 
     30; 
     31; @keyword XSEUIL 
     32; To eliminate segments of coasts which are to big (which link points which can  
     33; be close on the sphere but distant on the drawing). We delete all segments  
     34; whose the size surpass the size of the window following X/XSEUIL. By default,  
     35; XSEUIL=5. But it can be to big if we do a big zoom or a little one for some  
     36; projections... We specify it the keyword thanks to this keyword! 
     37; 
     38; @keyword YSEUIL 
     39; See XSEUIL 
     40; 
     41; @keyword _EXTRA 
     42; Used to pass our keywords 
     43 
     44; @returns 
     45;  
     46;  
     47; @uses 
     48;  
     49;  
     50; @restrictions 
     51;  
     52;  
     53; @examples 
     54;  
     55;  
     56; @history 
     57;  
     58;  
     59; @version  
     60; 
     61; @todo 
     62; Seb: remplir le header 
     63;  
     64;- 
     65 
    166PRO drawcoast_e, mask, xf, yf, nx, ny, COAST_COLOR = coast_color, COAST_THICK = coast_thick, YSEUIL = yseuil, XSEUIL = xseuil, onemore = onemore, _extra = ex 
    267;--------------------------------------------------------- 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/drawsectionbottom.pro

    r134 r150  
    33;------------------------------------------------------------ 
    44;+ 
    5 ; NAME:drawsectionbottom 
    65; 
    7 ; PURPOSE:fill and draw the bottom continents for a real section. 
     6; @file_comments 
     7; Fill and draw the bottom continents for a real section. 
    88; 
    9 ; CATEGORY: 
     9; @categories 
     10;  
     11; @param MASKIN {in}{requierd} 
    1012; 
    11 ; CALLING SEQUENCE: 
    12 ;  
    13 ; INPUTS: 
     13; @param XXAXISIN {in}{requierd}  
    1414; 
    15 ; KEYWORD PARAMETERS: 
     15; @param DEPTHSIN {in}{requierd} 
    1616; 
    17 ;        COAST_COLOR: the color of the coastline. 
    18 ;                     defaut value is 0 => black 
     17; @keyword COAST_COLOR 
     18; The color of the coastline. 
     19; defaut value is 0 => black 
    1920; 
    20 ;        COAST_THICK: the thickness of the coastline. 
    21 ;                     defaut value is 1 
     21; @keyword COAST_THICK 
     22; The thickness of the coastline. 
     23; defaut value is 1 
    2224; 
    23 ;        CONT_COLOR: the color of the continent. defaut value is 
    24 ;        (!d.n_colors - 1) < 255 => white 
     25; @keyword CONT_COLOR 
     26; The color of the continent. defaut value is 
     27; (!d.n_colors - 1) < 255 => white 
    2528; 
    26 ; OUTPUTS: 
     29; @uses 
     30; common.pro 
    2731; 
    28 ; COMMON BLOCKS:common.pro 
    29 ; 
    30 ; SIDE EFFECTS: 
    31 ; 
    32 ; RESTRICTIONS:simple way to fill continents for a section (using the 
     32; @restrictions 
     33; Simple way to fill continents for a section (using the 
    3334; fact that continents are wider at the bottom than at the top).  
    3435; 
    35 ; EXAMPLE: 
     36; @history 
     37; Sebastien Masson (smasson@lodyc.jussieu.fr) 
     38;                      June 14, 2002 
    3639; 
    37 ; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr) 
    38 ;                      June 14, 2002 
     40; @version 
     41; $Id$ 
     42; 
     43; @todo 
     44; Seb: definir params 
    3945;- 
    4046;------------------------------------------------------------ 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/fillcornermask.pro

    r134 r150  
    33;------------------------------------------------------------ 
    44;+ 
    5 ; NAME: FILLCORNERMASK 
    65; 
    7 ; PURPOSE: pour colorier proprement les continents! (c''est une longue 
    8 ; histoire...) 
     6; @file_comments 
     7; To color cleanly continents 
    98; 
    10 ; CATEGORY: pour plt 
     9; @categories 
     10; graphic 
    1111; 
    12 ; CALLING SEQUENCE: completecointerre 
    13 ;  
    14 ; INPUTS: non 
     12; @keyword _EXTRA 
     13; Used to pass your keywords  
    1514; 
    16 ; KEYWORD PARAMETERS:  _EXTRA  
     15; @keyword CONT_COLOR 
     16; The color of the continent. defaut value is 
     17; (!d.n_colors - 1) < 255 => white 
    1718; 
    18 ;        CONT_COLOR: the color of the continent. defaut value is 
    19 ;        (!d.n_colors - 1) < 255 => white 
     19; @uses 
     20; common.pro 
    2021; 
    21 ; OUTPUTS: non 
     22; @history 
     23; Sebastien Masson (smasson@lodyc.jussieu.fr) 
     24;                      8/8/2002 
    2225; 
    23 ; COMMON BLOCKS: common.pro 
     26; @version 
     27; $Id$ 
    2428; 
    25 ; SIDE EFFECTS: 
    26 ; 
    27 ; RESTRICTIONS: 
    28 ; 
    29 ; EXAMPLE: 
    30 ; 
    31 ; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr) 
    32 ;                      8/8/2002 
    3329;- 
    3430;------------------------------------------------------------ 
     
    4541  if NOT keyword_set(coinmonte) AND NOT keyword_set(coindescend) then return 
    4642; 
    47   tempsun = systime(1)          ; pour key_performance 
     43  tempsun = systime(1)          ; For key_performance 
    4844; 
    4945  IF NOT keyword_set(cont_color) THEN cont_color = (!d.n_colors-1) <  255     
    5046;------------------------------------------------------------ 
    51 ; definition descoordonnees des points numerotes 1,2,3,4,5,6 cf. les 
    52 ; schemas en dessous! 
     47; definition of coordinates of points numbered 1,2,3,4,5,6 (see figures below) 
    5348;------------------------------------------------------------ 
    5449; 
     
    6459; 
    6560; 
    66 ; cas coin terre en montee: 
    67 ;      2 points terre en diagonale montante avec 2 points mer sur 
    68 ;      la diagonale descendante. 
     61; Case land corner in ascent: 
     62;      2 land points in diagonal ascending with 2 ocean points on the descendant diagonal. 
    6963; 
    7064;                     3     
     
    9690  endif 
    9791;------------------------------------------------------------ 
    98 ; cas coin terre en descendante.: 
    99 ;      2 points terre en diagonale descendante avec 2 points mer sur 
    100 ;      la diagonale montante 
     92; Case land corner in descent: 
     93;      2 land points in diagonal descending with 2 ocean points on the ascendant diagonal. 
    10194; 
    10295;                     4 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/section.pro

    r134 r150  
    1 ;------------------------------------------------------------ 
    2 ;------------------------------------------------------------ 
    3 ;------------------------------------------------------------ 
    41;+ 
    5 ; NAME: 
    6 ; 
    7 ; PURPOSE: 
    8 ; 
    9 ; CATEGORY: 
    10 ; 
    11 ; CALLING SEQUENCE: 
    12 ;  
    13 ; INPUTS: 
    14 ; 
    15 ; KEYWORD PARAMETERS: 
    16 ; 
    17 ; OUTPUTS: 
    18 ; 
    19 ; COMMON BLOCKS:common.pro 
    20 ; 
    21 ; SIDE EFFECTS: 
    22 ; 
    23 ; RESTRICTIONS: 
    24 ; 
    25 ; EXAMPLE: 
    26 ; 
    27 ; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr) 
    28 ; 
     2; @file_comments 
     3; 
     4; 
     5; @categories 
     6; 
     7;  
     8; @param FIELD 
     9; 
     10;  
     11; @param RES 
     12; 
     13;  
     14; @param GLAMAXE 
     15; 
     16;  
     17; @param GPHIAXE 
     18; 
     19;  
     20; @keyword ENDPOINTS  
     21;  
     22;  
     23; @keyword BOXZOOM 
     24; 
     25;  
     26; @keyword TYPE  
     27; 
     28; 
     29; @keyword WDEPTH 
     30; 
     31; 
     32; @keyword DIREC 
     33; 
     34; 
     35; @keyword SHOWBUILD 
     36; 
     37; 
     38; @keyword ONLYBOX 
     39; 
     40; 
     41; @keyword _EXTRA  
     42; Used to pass your keywords 
     43;  
     44; @returns 
     45;  
     46;  
     47; @uses 
     48; common.pro 
     49;  
     50; @restrictions 
     51;  
     52;  
     53; @examples 
     54;  
     55;  
     56; @history 
     57; Sebastien Masson (smasson@lodyc.jussieu.fr) 
     58;  
     59; @version  
     60;  
    2961;- 
    30 ;------------------------------------------------------------ 
    31 ;------------------------------------------------------------ 
    32 ;------------------------------------------------------------ 
    3362 
    3463PRO section, field, res, glamaxe, gphiaxe, ENDPOINTS = endpoints $ 
     
    5281;--------------------- 
    5382;--------------------- 
    54 ; definition de boxzoom en fonction de endpoints 
    55 ; puis redefinition du domaine 
     83; definition of boxzoom in function of endpoints, then redefinition of the domain 
    5684  boxzoom2d = [min([endpoints[0], endpoints[2]], max = ma02), ma02 $ 
    5785               , min([endpoints[1], endpoints[3]], max = ma13), ma13] 
     
    83111; how can we find a good test value? 
    84112  IF nx * ny LE 10 THEN domdef, GRIDTYPE = grillechoice, _extra = ex 
    85 ; on redefinit lon1, ... au cas ou findalways ait ete utilise ds domdef 
     113; We redefine lon1, ... in case findalways has been used in domdef 
    86114  lon1 = min([endpoints[0], endpoints[2]], max = lon2) 
    87115  lat1 = min([endpoints[1], endpoints[3]], max = lat2) 
     
    89117; until its bottom part. 
    90118  if strpos(type, 'z') NE -1 THEN BEGIN 
    91 ;on garde les yranges (axe z) avant de changer la boxzoom. 
     119; We keep yranges (axis z) before changing the boxzoom. 
    92120    !y.range = [localbox[nelbox-1], localbox[nelbox-2]]  
    93121    if vargrid EQ 'W' OR keyword_set(wdepth) then BEGIN 
     
    105133  ENDIF 
    106134  !y.range = [localbox[nelbox-1], localbox[nelbox-2]]  
    107 ; on recupere la grille sur la boxzoom 
     135; We recuperate the grid on the boxzoom. 
    108136; 
    109137  IF NOT keyword_set(onlybox) THEN saveboxparam, 'boxparam4section.dat' 
     
    126154 
    127155;--------------------- 
    128 ; on definit la triangulation qui va nous permetre de determiner la 
    129 ; section. on la recalcule car elle doit etre definie sur la terre 
    130 ; aussi bien que sur la mer. suivant le sens de la section -plutot 
    131 ; longitude ou plutot latitude- on definit la facon de trianguler 
     156; We define the triangulation which will allows us to determinate the section.  
     157; We recalculate it because it must be defined on the Earth and on oceans. 
     158; Following the direction of the section -rather longitude or rather latitude-  
     159; we define the way to triangulate. 
    132160  if strpos(type, 'x') NE -1 then BEGIN  
    133161    downward = (lindgen(nx, ny))[0:nx-2, 0:ny-2] 
     
    146174  ENDIF 
    147175;--------------------- 
    148 ; equation de la droite suivant laquelle on fait la section 
     176; Equation of the line on which we do the section. 
    149177;--------------------- 
    150178  abc = linearequation(endpoints[0:1], endpoints[2:3]) 
    151179  glamtri = glam[tri] 
    152180  gphitri = gphi[tri] 
    153 ; quels sont les points de la triangulation qui sont au dessus et au 
    154 ; dessous de la droite ? 
     181; Which points of the triangulation are above and below the line? 
    155182  if abc[1] NE 0 THEN $ 
    156183    test = temporary(gphitri) GE (-abc[0]/abc[1]*temporary(glamtri)-abc[2]/abc[1]) $ 
     
    158185 
    159186  zero123 = total(test, 1) 
    160 ; to keep: triangles de la triangulation qui sont a cheval sur la droite. 
     187; to keep: triangles of the triangulation which are over the line. 
    161188  tokeep1 = where(zero123 EQ 1) 
    162189  tokeep2 = where(temporary(zero123) EQ 2) 
     
    165192  test = test[*, tokeep] 
    166193  tri = tri[*, tokeep] 
    167 ; quel est le sommet du triangle qui est seul d''un cote de la droite? 
     194; Which summit of the triangle is alone in a side of the line? 
    168195  single1 = where(test[*, 0:n_elements(tokeep1)-1] EQ 1) 
    169196  single1 = single1-(single1/3)*3 
     
    177204 
    178205  single = [temporary(single1), temporary(single2)] 
    179 ; points1 le point du triangles qui est seul d''un cote de la droite. 
    180 ; point2 l''autre point du triangle de l''autre cote de la droite 
     206; points1 the point (of the triangle) alone in a side of the line. 
     207; point2 the other point of the triangle in the other side of the line. 
    181208  point1 = [single, single] 
    182209  point2 = [single EQ 0, 1 + (single LE 1)] 
     
    189216  points1 = tri[point1, index] 
    190217  points2 = tri[point2, temporary(index)] 
    191 ; points : complexe contenant les couples de points de part et 
    192 ; d''autre de la droite. Ils faut supprimer les doublons. 
     218; points : complex containing couples of points in a side and the other  
     219; side of the line. We have to delete duplicates. 
    193220  points = dcomplex(points1, points2) 
    194221  points = points[uniq(points, sort(points))] 
    195222  symetrique = dcomplex(imaginary(points), double(points)) 
    196223  points = points[where(points-shift(temporary(symetrique), 1) NE 0)] 
    197 ; points1 les coordnnees du point du triangles qui est seul d''un cote de la droite. 
    198 ; point2 les coordnnees de l''autre point du triangle de l''autre cote de la droite 
     224; points1 coordinates of the point of the triangle which is alone in a side of the line. 
     225; point2 coordinates of the other point of the triangle in the other side of the line. 
    199226  points1 = complex(glam[   double(points)], gphi[   double(points)]) 
    200227  points2 = complex(glam[imaginary(points)], gphi[imaginary(points)]) 
    201 ; droites les equations des droites dont on cherche l''intersection 
    202 ; avec la section. 
     228; droites equations of line whose we look for the intersection wit the section. 
    203229  droites = linearequation(points1, points2) 
    204230  inter = lineintersection(droites, abc#replicate(1, n_elements(points1) ) ) 
    205231 
    206 ; les ccordonnes geographiques des points que l''on cherche sur la section. 
     232; Geographic coordinates of points we look for on the section. 
    207233  glamaxe = float(inter) 
    208234  gphiaxe = imaginary(inter) 
    209 ; on les range ds l''ordre croissant entre les bornes de la section 
     235; We arrange them in the growing order between boundaries of the section. 
    210236  if strpos(type, 'x') NE -1 then BEGIN  
    211237    sort = sort(glamaxe) 
     
    272298      poids = poids#replicate(1, nz) 
    273299      res = poids*value1+(1-poids)*value2 
    274 ; moyenne suivant z ? 
     300; average following z ? 
    275301      if strpos(type, 'z') EQ -1 then begin 
    276302        nan = where(finite(res) EQ 0) 
     
    319345      poids = reform(poids, npoints, nz, jpt, /over) 
    320346      res = poids*value1+(1-poids)*value2 
    321 ; moyenne suivant z ? 
     347; average following z ? 
    322348      if strpos(type, 'z') EQ -1 then begin 
    323349        nan = where(finite(res) EQ 0) 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/tracecote.pro

    r134 r150  
    33;------------------------------------------------------------ 
    44;+ 
    5 ; NAME:tracecote 
    65; 
    7 ; PURPOSE: dessine les cotes ds plt 
     6; @file_comments 
     7; Draw coasts in plt. 
    88; 
    9 ; CATEGORY: pour faire un joli dessin 
     9; @categories 
     10; graphic 
    1011; 
    11 ; CALLING SEQUENCE:tracecote,mask 
    1212;  
    1313; INPUTS:mask le tableau mask sur la zone consideree pour le dessin 
    1414; 
    15 ; KEYWORD PARAMETERS:  
     15; @keyword SURFACE_COASTLINE 
     16; To draw the surface coast line instead of 
     17; the coast line at level firstz[tw]. Usefull only for deep 
     18; plots! 
    1619; 
    17 ;        COAST_COLOR: the color of the coastline. 
    18 ;                     defaut value is 0 => black 
     20; @keyword _EXTRA 
     21; used to pass your keywords 
    1922; 
    20 ;        COAST_THICK: l''epaisseur du trait pour tracer les 
    21 ;        continents. par defaut c''est 1. 
     23; @uses 
     24; common.pro 
    2225; 
    23 ;        /SURFACE_COASTLINE: to draw the furface coast line instead of 
    24 ;        the coast line at level firstz[tw]. Usefull only for deep 
    25 ;        plots! 
     26; @history 
     27; Sebastien Masson (smasson@lodyc.jussieu.fr) 
     28;                      30/9/1999 
    2629; 
    27 ;        XSEUIL: pour eliminer les segments de cote qui sont trop 
    28 ;        grand (qui relient des points qui peuvent etre tres proches 
    29 ;        sur la sphere mais tres eloignes sur le dessin) on supprime 
    30 ;        tous les egments dot la taille depasse: 
    31 ;              taille de la fenetre suivant X/ xseuil. 
    32 ;        Par defaut xseuil est egale a 5. masi peut etre trop grand si 
    33 ;        on fait un fort zoom ou trout petit pour certaines 
    34 ;        projections... le specifier alors a l''aide de ce mot cle! 
     30; @version 
     31; $Id$ 
    3532; 
    36 ;        YSEUIL: cf. xseuil 
    37 ; 
    38 ; OUTPUTS: rien 
    39 ; 
    40 ; COMMON BLOCKS:common.pro 
    41 ; 
    42 ; SIDE EFFECTS: 
    43 ; 
    44 ; RESTRICTIONS: 
    45 ; 
    46 ; EXAMPLE: 
    47 ; 
    48 ; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr) 
    49 ;                      30/9/1999 
    5033;- 
    5134;------------------------------------------------------------ 
     
    6447  ENDIF 
    6548;-------------------------------------------------------------- 
    66    tempsun = systime(1)         ; pour key_performance 
     49   tempsun = systime(1)         ; For key_performance 
    6750   if n_elements(key_gridtype) EQ 0 then key_gridtype = 'c' 
    6851; 
    69 ; on agrandi un peu le cadre definit par les premier..., dernier... de 
    70 ; facon a bien recuperer les bords de cote qui sont en bordure du 
    71 ; domaine a tracer 
     52; We enlarge a bit the frame defined by firsts..., lasts... in order to  
     53; recuperate edges of the coast which are in the edging of the domain. 
    7254;    
    73    tempdeux = systime(1)        ; pour key_performance =2 
     55   tempdeux = systime(1)        ; For key_performance =2 
    7456   firstx = 0 > (min([firstxt, firstxf])-1) 
    7557   lastx = (max([lastxt, lastxf])+1) < (jpi-1) 
     
    7860   nx = lastx-firstx+1 
    7961   ny = lasty-firsty+1 
    80 ; quel niveau vertical choisir ? 
     62; Which vertical level choose? 
    8163   IF keyword_set(surface_coastline) THEN firstz = 0 ELSE $ 
    8264     IF strupcase(vargrid) eq 'W' THEN firstz = firstzw ELSE firstz = firstzt 
    83 ; attribution du masque et des coordonnes delimitant les limites de la 
    84 ; terre (coordonnees f) 
     65; Attribution of the mask and of coordinates  delimiting limits of the land (coordinates f). 
    8566   mask = tmask[firstx:lastx, firsty:lasty, firstz] 
    8667   xf = glamf[firstx:lastx, firsty:lasty] 
     
    9071; 
    9172   if key_gridtype EQ 'e' then onemore = xf[0, 0] gT xf[0, 1] 
    92 ; on passe en coordonnee normaliser pour pouvoir s'affranchir du type 
    93 ; de projection choisie et du suport surlequel on fait le dessin 
    94 ; (ecran ou postscript) 
     73; We pass in normalized coordinates to be able to become independant from the projection's  
     74; type choosen and from the support on which we do the drawing (screen or postscript) 
    9575   z = convert_coord(xf[*],yf[*],/data,/to_normal)  
    9676   xf = reform(z[0, *], nx, ny) 
     
    9878   tempvar = SIZE(TEMPORARY(z)) 
    9979; 
    100 ; attention, suivant la projection certains points x ou y peuvent 
    101 ; devenir NaN (cf. points deriere la terre ds une projection 
    102 ; orthographique) 
     80; Beware, following the projection, some points x or y can become NaN (see point  
     81; behind the earth in an orthographic projection). 
    10382; 
    104 ; on met les points a eliminer a une tres gande valeur comme ca il ne 
    105 ; passerons pas le test avec distanceseuil (cf. plus bas) 
     83; We put points to be eliminated at a very big value so that they will not pass the  
     84; test with distanceseuil (see further). 
    10685; 
    10786   if (!map.projection LE 7 AND !map.projection NE 0) $ 
     
    11796   ind = where(yf LT !p.position[1] OR yf GT !p.position[3]) 
    11897   IF ind[0] NE -1 THEN yf[ind] = 1e5 
    119    tempvar = SIZE(TEMPORARY(ind)) ; on efface ind 
     98   tempvar = SIZE(TEMPORARY(ind)) ; we delete ind 
    12099; 
    121100   if n_elements(key_gridtype) EQ 0 then key_gridtype = 'c' 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/tracemask.pro

    r134 r150  
    33;------------------------------------------------------------ 
    44;+ 
    5 ; NAME:tracemask 
    65; 
    7 ; PURPOSE:dessiner des contour d''un mask 
     6; @file_comments 
     7; Draw contours of a mask 
    88; 
    9 ; CATEGORY:plus simple que tracecote, car ne s''occuppe pas du type de 
    10 ; projection et de la periodicite de la grille 
     9; @categories 
     10; utilities 
     11;  
     12; @param MASKENTREE {in}{required} 
     13; 2d array specifing the mask 
     14;  
     15; @param XIN {in}{required}, 
     16; 2d array specifing longitude coordinates. 
     17;  
     18; @param YIN {in}{required}, 
     19; 2d array specifing latitude coordinates. 
    1120; 
    12 ; CALLING SEQUENCE: tracemask, maskentree, xentree, yentree 
    13 ;  
    14 ; INPUTS:maskentree, xentree, yentree tableaux 2d specifiant le mask 
    15 ; et ses coordonees en longitude te latitude. 
     21; @keyword COAST_COLOR 
     22; The color of the coastline. 
     23; defaut value is 0 => black 
    1624; 
    17 ; KEYWORD PARAMETERS: 
     25; @keyword COAST_THICK 
     26; The thick of the trait to trace continents 
     27; It is 1 by default. 
    1828; 
    19 ;        COAST_COLOR: the color of the coastline. 
    20 ;                     defaut value is 0 => black 
     29; @keyword OVERPLOT 
     30; To do a plot over an other one. 
    2131; 
    22 ;        COAST_THICK: l''epaisseur du trait pour tracer les 
    23 ;        continents. par defaut c''est 1. 
     32; @keyword _EXTRA 
     33; used to pass your keywords 
    2434; 
    25 ; OUTPUTS: none 
     35; @uses 
     36; common.pro 
    2637; 
    27 ; COMMON BLOCKS:common.pro 
     38; @history 
     39; Sebastien Masson (smasson@lodyc.jussieu.fr) 
    2840; 
    29 ; SIDE EFFECTS: 
    30 ; 
    31 ; RESTRICTIONS: 
    32 ; 
    33 ; EXAMPLE: 
    34 ; 
    35 ; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr) 
     41; @version 
     42; $Id$ 
    3643; 
    3744;- 
     
    5158  ENDIF 
    5259;---------------------------------------------------------  
    53    tempsun = systime(1)         ; pour key_performance 
    54 ; on s''afranchit des problemes de bord: 
    55    tempdeux = systime(1)        ; pour key_performance =2 
     60   tempsun = systime(1)         ; For key_performance 
     61; We avoid edging problems: 
     62   tempdeux = systime(1)        ; For key_performance =2 
    5663   tailleentree = size(maskentree) 
    5764   nx = tailleentree[1]+1 
     
    6269  IF n_elements(yin) EQ 0 THEN yentree = findgen(ny-1) ELSE yentree = yin 
    6370  IF (size(yentree))[0] EQ 1 THEN yentree = replicate(1,nx-1)#yentree 
    64 ; on agrandi le mask de une colonne a gauche et de une colonne en bas 
     71; We enlarge the mask by 1 column to the left an d1 line to the bottom 
    6572   mask = intarr(tailleentree[1]+1, tailleentree[2]+1) 
    6673   mask[1:tailleentree[1], 1:tailleentree[2]] = maskentree 
    67 ; les 2 premieres colonnes sont identiques 
     74; The 2 first columns are identical. 
    6875   mask[0, 1:tailleentree[2]] = maskentree[0, *] 
    69 ; les 2 premieres lignes sont identiques 
     76; The 2 first lines are identical. 
    7077   mask[1:tailleentree[1], 0] = maskentree[*, 0] 
    71 ; on calcul la position suivant x des points qui seviront a tracer le 
    72 ; masque. ils sont situes entre chaque points du masque, sauf pour la 
    73 ; derniere colonne que l''on ne peut pas calculer et que l''on met 
    74 ; donc a max(!x.range) 
    75    xrange = !x.range[sort(!x.range)] ; si reverse_x est utilise! 
     78; We calculate the position following x of points which will serve to trace the mask. They are situated between each points of the mask, exept for the last column we can not calculate and so we put at max (!x.range). 
     79   xrange = !x.range[sort(!x.range)] ; if REVERSE_X is used 
    7680   xentree = .5*(xentree+shift(xentree, -1, 0)) 
    7781   IF not keyword_set(overplot) THEN xentree[nx-2, *] = xrange[1] $ 
    7882   ELSE xentree[nx-2, *] = xentree[nx-3, *] 
    79 ; on seuil 
     83; we sill 
    8084   xentree = xrange[0] > xentree < xrange[1] 
    81 ; on agrandit le tableau 
     85; we enlarge the array  
    8286   xf = fltarr(nx, ny) 
    8387   xf[1:nx-1, 1:ny-1] = xentree 
     
    103107    print, 'temps tracemask: determination du mask et des ses coordonnes', systime(1)-tempdeux 
    104108; 
    105 ; on trace les segments verticaux: 
     109; We trace vertical segments: 
    106110; 
    107    tempdeux = systime(1)        ; pour key_performance =2 
     111   tempdeux = systime(1)        ; For key_performance =2 
    108112   liste = where(mask+shift(mask, -1, 0) EQ 1) 
    109113   IF liste[0] NE -1 THEN BEGIN 
    110 ; on recupere lx et ly qui sont les indices ds un tableau 2d des 
    111 ; points donnes par liste 
     114; We recuperate lx and ly which are indexes in a 2d array of points given by list 
    112115      ly = liste/nx & lx = temporary(liste)-nx*ly 
    113       indice = where(ly NE 0)   ; on ne prend pas les points concernant  
    114 ; la premiere ligne car ds ce cas le pt j-1 n''est pas definit 
     116      indice = where(ly NE 0) ; We do not take points concernining  
     117; the first line because in this case, the point j-1 is not defined  
    115118      if indice[0] NE -1 then begin 
    116119         lx = lx[indice] & ly = ly[temporary(indice)] 
    117120         IF testvar(var = key_performance) EQ 2 THEN $ 
    118121          print, 'temps tracemask: liste traits verticaux', systime(1)-tempdeux 
    119          tempdeux = systime(1)  ; pour key_performance =2 
    120 ; boucle sur les points concernes et trace du segment 
    121 ; rq: on utilise plots au lieu de plot car plots est bcp plus rapide. 
     122         tempdeux = systime(1)  ; For key_performance =2 
     123; loop on concerned points and drawing of the segment. 
     124; comments: we use plots instead of plot because plots is faster. 
    122125         for pt = 0L, n_elements(lx)-1 do BEGIN  
    123126            i = lx[pt] & j = ly[pt] 
     
    132135   ENDIF 
    133136; 
    134 ; on trace les segments horizontaux: 
     137; We trace horizontal segments: 
    135138; 
    136    tempdeux = systime(1)        ; pour key_performance =2 
     139   tempdeux = systime(1)        ; For key_performance =2 
    137140   liste = where(mask+shift(mask, 0, -1) EQ 1) 
    138141   IF liste[0] NE -1 THEN BEGIN 
    139142      ly = liste/nx & lx = temporary(liste)-nx*ly 
    140       indice = where(lx NE 0)   ; on ne prend pas les points de la  premiere colonne 
     143      indice = where(lx NE 0)   ; We do not take point sof the first column. 
    141144      if indice[0] EQ -1 then return 
    142145      lx = lx[indice] & ly = ly[temporary(indice)] 
    143146      IF testvar(var = key_performance) EQ 2 THEN $ 
    144147       print, 'temps tracemask: liste traits horizontaux', systime(1)-tempdeux 
    145       tempdeux = systime(1)     ; pour key_performance =2 
     148      tempdeux = systime(1)     ; For key_performance =2 
    146149      for pt = 0L, n_elements(lx)-1 do BEGIN  
    147150         i = lx[pt] & j = ly[pt] 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/triangule.pro

    r134 r150  
     1;------------------------------------ 
     2;+ 
     3; 
     4; @todo 
     5; seb 
     6; 
     7;- 
     8;------------------------------------ 
    19FUNCTION triangule, maskentree, BASIC = basic, COINMONTE = coinmonte, COINDESCEND = coindescend, _extra = ex 
    210; 
  • trunk/SRC/ToBeReviewed/TRIANGULATION/triangule_c.pro

    r134 r150  
    33;------------------------------------------------------------ 
    44;+ 
    5 ; NAME:triangule_c 
    6 ; 
    7 ; PURPOSE:construit le tableau de triangulation. 
    8 ; 
    9 ; L''idee est de 
    10 ; construire une liste de triangles qui relient les points entre 
    11 ; eux. Ceci est fait automatiquement avec la fonction TRIANGULATE. 
    12 ;  ICI: 
    13 ; on tient compte du fait que les points sont disposes sur une grille 
    14 ; (reguliere ou pas, mais pas destructuree, cad que les points sont 
    15 ; ecrits suivant une matrice rectangulaire). Un moyen tres simple de 
    16 ; faire des triangles entre tous les points est alors: 
    17 ; 
    18 ;     pour chaque point (i,j) de la matrice -sauf ceux de la derniere 
    19 ;     ligne et de la derniere colonne- on on appelle le rectangle 
    20 ;     (i,j) le rectangle forme par les 4 points (i,j), (i+1,j), 
    21 ;     (i,j+1), (i+1,j+1). Pour tracer tous les triangles, il suffit de 
    22 ;     tracer les 2 triangles contenus ds les rectangles (i,j) 
    23 ; 
    24 ; au passage on remarque que chaque rectangle (i,j) possede 2 diagonales (si 
    25 ; si faites un dessin c''est vrai), il y a donc 2 choix possibles pour 
    26 ; chaque rectangles qd on veut le couper en 2 triangles... 
     5; 
     6; @file_comments 
     7; Construct the triangulation array. 
     8; 
     9; The idea is: construct a list of triangle which link points between them.  
     10; This is automatically done by the function TRIANGULATE 
     11;  Here: 
     12; we consider the fact that points are disposed on a grid (regular or not,  
     13; but not unstructured, that is to say that points are written following a  
     14; rectangular matrix). A easy way to do triangles between all points is then:  
     15; 
     16;     for each point (i,j) of the matrix -exept those of the last line and of 
     17;     the last column- we call rectangle (i,j) the rectangle made of the four 
     18;     points (i,j), (i+1,j), (i,j+1), (i+1,j+1). To trace all triangle, we just 
     19;     have to trace the 2 triangles contained in rectangles (i,j) 
     20; 
     21; We notice that each rectangle (i,j) have 2 diagonals (it is true... Make a 
     22; drawing to make sure!!), so there are two possible choice for each rectangle 
     23; we want to cut in 2 triangles... 
    2724;  
    28 ; C''est grace a ce choix que l''on va pouvoir tracer les cotes avec 
    29 ; des angles droits. A chaque angle de cote remarquable par 
    30 ; l''existance d''un unique point terre ou d''un unique point mer sur 
    31 ; les 4 cotes d''un rectangle (i,j), il faut couper le rectangle 
    32 ; suivant la diagonale qui qui passe par le point singulier. 
     25; It is thanks to this choice that we will be able to trace coast with right 
     26; angles. At each angle of coast remarkable by the existence of an unique land 
     27; point or of an unique ocean point on one of the four summit of a rectangle (i,j), 
     28; we have to cut the rectangle following the diagonal passing by this point. 
    3329;  
    34 ; CATEGORY:pour faire de beaux graphiques masques 
    35 ; 
    36 ; CALLING SEQUENCE:res=triangule([mask]) 
    37 ; 
    38 ; INPUTS:optionnel:mask c''est le tableau 2d qui sevira a masquer le 
    39 ; champ que l''on tracera apres avec CONTOUR, 
     30; @categories 
     31; graphic 
     32; 
     33; @param MASKENTREE {in}{optional} 
     34; It is a 2d array which will serve to mask the field we will trace after with CONTOUR,  
    4035; ...TRIANGULATION=triangule(mask) 
    41 ; si cet argument n''est pas specifie, la function utilise tmask. 
    42 ; 
    43 ; KEYWORD PARAMETERS: 
    44 ; 
    45 ;       /BASIC: specifie que le masque est sur une grille basice 
    46 ;       (utiliser pour la triangulation ds les coupes verticales et 
    47 ;       des hovmoellers) 
    48 ; 
    49 ;       /KEEP_CONT: to keep the triangulation even on the continents 
    50 ; 
    51 ;       COINMONTE=tableau, pour obtenir le tableau de "coins de terre 
    52 ;       montant" a traiter avec completecointerre.pro ds la variable 
    53 ;       tableau plutot que de la faire passer par la variable globale 
    54 ;       twin_corners_up. 
    55 ; 
    56 ;       COINDESCEND=tableau cf COINMONTE 
    57 ; 
    58 ; OUTPUTS: 
    59 ;       res: tableau 2d (3,nbre de triangles). 
    60 ;    chaque ligne de res represente les indices des points 
    61 ;    constituants les sommets d''un triangle. 
    62 ;    cf. comment on trace les triangles ds dessinetri.pro 
    63 ; 
    64 ; COMMON BLOCKS: 
    65 ;       common.pro different.pro definetri.pro 
    66 ; 
    67 ; SIDE EFFECTS: 
    68 ; 
    69 ; RESTRICTIONS:les donnees dont un veut ensuite faire le contour 
    70 ; doivent etre disposees dans une matrice. Par contre dans la matrice, 
    71 ; la disposition des points peut ne pas etre irreguliere. Si les 
    72 ; donnees sont disposees completement de facon irreguliere, utiliser 
    73 ; TRIANGULE. 
    74 ; 
    75 ; EXAMPLE: 
    76 ; 
    77 ; MODIFICATION HISTORY: Sebastien Masson (smasson@lodyc.jussieu.fr) 
     36; If this argument is not specified, the function use tmask 
     37; 
     38; @keyword BASIC 
     39; Specify that the mask is on a basic grid (use the triangulation for vertical cuts and hovmoellers) 
     40; 
     41; @keyword KEEP_CONT 
     42; To keep the triangulation even on the continents 
     43; 
     44; @keyword COINMONTE 
     45; It is an array. To obtain the array of "ascending land corner" to be treated with  
     46; completecointerre.pro in the variable array instead of make it pass by the global  
     47; variable twin_corners_up. 
     48; 
     49; @keyword COINDESCEND 
     50; It is an array. See COINMONTE 
     51; 
     52; @returns 
     53; res: tableau 2d (3,nbre de triangles). 
     54; Each line of res represent indexes of points constituing summits of a triangle.  
     55; See how we trace triangles in definetri.pro 
     56; 
     57; @uses 
     58; common.pro 
     59; different.pro 
     60; definetri.pro 
     61; 
     62; @restrictions 
     63; Datas whose we want to do the contour must be disposed in a matrix.  
     64; On the other hand, in the matrix, the points's arrangement can not be  
     65; irregular. If it is, use TRIANGULE. 
     66; 
     67; @history 
     68; Sebastien Masson (smasson@lodyc.jussieu.fr) 
    7869;                       26/4/1999 
     70; 
     71; @version 
     72; $Id$ 
     73; 
     74; @todo 
     75; seb L.267->268 je ne pense pas que ce soit ce que tu voulais dire mais  
     76; c'est la traduction de ce qu'il y avait écrit. Correction si besoin. 
    7977;- 
    8078;------------------------------------------------------------ 
     
    8381FUNCTION triangule_c, maskentree, COINMONTE = coinmonte, COINDESCEND = coindescend, BASIC = basic, KEEP_CONT = keep_cont 
    8482; 
    85   compile_opt idl2, strictarrsubs 
    86 ; 
    87    tempsun = systime(1)         ; pour key_performance 
     83compile_opt idl2, strictarrsubs 
     84; 
     85tempsun = systime(1)            ; For key_performance 
    8886;--------------------------------------------------------- 
    8987@cm_4mesh 
    90   IF NOT keyword_set(key_forgetold) THEN BEGIN 
     88IF NOT keyword_set(key_forgetold) THEN BEGIN 
    9189@updatenew 
    92   ENDIF 
    93 ;------------------------------------------------------------ 
    94 ; le masque est donne ou il faut prendre tmask? 
    95 ;------------------------------------------------------------ 
    96 ; 
    97    msk = maskentree 
    98    taille = size(msk) 
    99    nx = taille[1] 
    100    ny = taille[2] 
    101 ; 
    102    IF n_elements(keep_cont) EQ 0 THEN keep_cont = 1-key_irregular 
    103 ;------------------------------------------------------------ 
    104    if keyword_set(key_periodic)*(nx EQ jpi) $ 
    105     AND NOT keyword_set(basic) then BEGIN  
    106       msk = [msk, msk[0, *]] 
    107       nx = nx+1 
    108    ENDIF 
    109 ;------------------------------------------------------------ 
    110 ; on va trouver la liste des rectangles (i,j) (reperes par leur coin 
    111 ; en bas a gauche) qu''il faut couper suivant une diagonale descendante 
    112 ; on appellera cette liste : pts_downward 
     90ENDIF 
     91;------------------------------------------------------------ 
     92; Is the mask given or do we have to take tmask? 
     93;------------------------------------------------------------ 
     94; 
     95msk = maskentree 
     96taille = size(msk) 
     97nx = taille[1] 
     98ny = taille[2] 
     99; 
     100IF n_elements(keep_cont) EQ 0 THEN keep_cont = 1-key_irregular 
     101;------------------------------------------------------------ 
     102if keyword_set(key_periodic)*(nx EQ jpi) $ 
     103  AND NOT keyword_set(basic) then BEGIN  
     104    msk = [msk, msk[0, *]] 
     105    nx = nx+1 
     106ENDIF 
     107;------------------------------------------------------------ 
     108; We will find the list of rectangles (i,j)(located by their left  
     109; bottom corner) we have to cut folowing a descendant diagonal.  
     110; We will call this list : pts_downward 
    113111;  
    114    pts_downward = 0 
    115  
    116 ; on construit le test qui permet de trouver un tel triangle: 
     112pts_downward = 0 
     113 
     114; We construct the test which allow to find this triangle : 
    117115; 
    118116; 
     
    124122;             msk---------------------shift(msk, -1,  0) 
    125123; 
    126    sum1 = msk+shift(msk, -1, 0)+shift(msk, -1, -1) ;pts qui entourrent le pt en haut a gauche 
    127    sum2 = msk+shift(msk, 0, -1)+shift(msk, -1, -1) ;pts qui entourrent le pt en bas a droite 
    128  
    129  
    130    tempdeux = systime(1)        ; pour key_performance =2 
    131 ; pt terre en haut a gauche entoure de pts mer 
    132    liste = where( (4-sum1)*(1-shift(msk, 0, -1)) EQ 1 ) 
    133    if liste[0] NE -1 THEN pts_downward = [pts_downward,liste ] 
    134 ; pt mer en haut a gauche entoure de pts terre 
    135    liste = where( (1-sum1)*shift(msk, 0, -1) EQ 1) 
    136    if liste[0] NE -1 THEN pts_downward = [pts_downward,liste ] 
    137 ; pt terre en bas a droite entoure de pts mer 
    138    liste = where( (4-sum2)*(1-shift(msk, -1,  0)) EQ 1) 
    139    if liste[0] NE -1 THEN pts_downward = [pts_downward,liste ] 
    140 ; pt mer en bas a droite entoure de pts terre 
    141    liste = where( (1-sum2)*shift(msk, -1,  0) EQ 1) 
    142    if liste[0] NE -1 THEN pts_downward = [pts_downward,liste ] 
    143    undefine, liste 
    144 ; 
    145    IF testvar(var = key_performance) EQ 2 THEN $ 
    146     print, 'temps triangule: trouver pts_downward', systime(1)-tempdeux 
    147 ; 
    148    if (NOT keyword_set(basic)) OR keyword_set(coinmonte) OR keyword_set(coindescend) then begin 
    149       tempdeux = systime(1)     ; pour key_performance =2 
    150 ;2 points terre en diagonale montante avec 2 points mer sur la diagonale descendante 
    151       coinmont = where( (1-msk)*(1-shift(msk, -1, -1)) $ 
    152                         *(shift(msk, 0, -1)*shift(msk, -1,  0) EQ 1) ) 
    153       if coinmont[0] NE -1 THEN pts_downward = [pts_downward, coinmont] 
    154 ; 
    155       IF testvar(var = key_performance) EQ 2 THEN $ 
    156        print, 'temps triangule: trouver coinmont', systime(1)-tempdeux 
    157       tempdeux = systime(1)     ; pour key_performance =2 
    158 ; 
    159 ;2 points terre en diagonale descendante avec 2 points mer sur la diagonale montante 
    160       coindesc = where( ((1-shift(msk,  0, -1))*(1-shift(msk, -1, 0)) $ 
    161                          *msk*shift(msk, -1, -1) EQ 1) ) 
    162 ; 
    163       IF testvar(var = key_performance) EQ 2 THEN $ 
    164        print, 'temps triangule: trouver coindesc', systime(1)-tempdeux 
    165 ; 
    166     ENDIF 
    167 ; 
    168    if n_elements(pts_downward) EQ 1 then BEGIN  
    169       tempdeux = systime(1)     ; pour key_performance =2 
    170 ; 
    171       triang = definetri(nx, ny) 
    172 ; 
    173       IF testvar(var = key_performance) EQ 2 THEN $ 
    174        print, 'temps triangule: definetri', systime(1)-tempdeux 
    175       coinmont = -1 
    176       coindesc = -1 
    177    ENDIF ELSE BEGIN  
    178       tempdeux = systime(1)     ; pour key_performance =2 
    179       pts_downward = pts_downward[1:n_elements(pts_downward)-1] 
    180       pts_downward = pts_downward[uniq(pts_downward, sort(pts_downward))] 
    181 ; aucun rectangle ne peut avoir comme coin en bas a gauche un element 
    182 ; de la derniere colonne ou de la derniere ligne. 
    183 ; il faut donc enlever ces points si ils ont ete selectionnes dans 
    184 ; pts_downward. 
    185       derniere_colonne = (lindgen(ny)+1)*nx-1  
    186       derniere_ligne = lindgen(nx)+(ny-1)*nx  
    187       pts_downward =different(pts_downward,derniere_colonne ) 
    188       pts_downward =different(pts_downward,derniere_ligne ) 
    189       if (NOT keyword_set(basic)) OR keyword_set(coinmonte) OR keyword_set(coindescend) then begin 
    190          if coinmont[0] NE -1 then begin 
     124sum1 = msk+shift(msk, -1, 0)+shift(msk, -1, -1)    ;points which surround the left top point. 
     125sum2 = msk+shift(msk, 0, -1)+shift(msk, -1, -1)    ;points which surround the right bottom point. 
     126 
     127 
     128tempdeux = systime(1)           ; For key_performance =2 
     129; The left top land point surrounded by ocean points 
     130liste = where( (4-sum1)*(1-shift(msk, 0, -1)) EQ 1 ) 
     131if liste[0] NE -1 THEN pts_downward = [pts_downward,liste ] 
     132; The left top ocean point surrounded by land points 
     133liste = where( (1-sum1)*shift(msk, 0, -1) EQ 1) 
     134if liste[0] NE -1 THEN pts_downward = [pts_downward,liste ] 
     135; The right bottom land point surrounded by ocean points 
     136liste = where( (4-sum2)*(1-shift(msk, -1,  0)) EQ 1) 
     137if liste[0] NE -1 THEN pts_downward = [pts_downward,liste ] 
     138; The right bottom ocean point surrounded by land points 
     139liste = where( (1-sum2)*shift(msk, -1,  0) EQ 1) 
     140if liste[0] NE -1 THEN pts_downward = [pts_downward,liste ] 
     141undefine, liste 
     142; 
     143IF testvar(var = key_performance) EQ 2 THEN $ 
     144  print, 'temps triangule: trouver pts_downward', systime(1)-tempdeux 
     145; 
     146if (NOT keyword_set(basic)) OR keyword_set(coinmonte) OR keyword_set(coindescend) then begin 
     147    tempdeux = systime(1)       ; For key_performance =2 
     148;2 land points in ascendant diagonal with 2 ocean points in descendant diagonal. 
     149    coinmont = where( (1-msk)*(1-shift(msk, -1, -1)) $ 
     150                      *(shift(msk, 0, -1)*shift(msk, -1,  0) EQ 1) ) 
     151    if coinmont[0] NE -1 THEN pts_downward = [pts_downward, coinmont] 
     152; 
     153    IF testvar(var = key_performance) EQ 2 THEN $ 
     154      print, 'temps triangule: trouver coinmont', systime(1)-tempdeux 
     155    tempdeux = systime(1)       ; pour key_performance =2 
     156; 
     157    coindesc = where( ((1-shift(msk,  0, -1))*(1-shift(msk, -1, 0)) $ 
     158                       *msk*shift(msk, -1, -1) EQ 1) ) 
     159; 
     160;2 land points in descendant diagonal with 2 ocean points in ascendant diagonal. 
     161    IF testvar(var = key_performance) EQ 2 THEN $ 
     162      print, 'temps triangule: trouver coindesc', systime(1)-tempdeux 
     163; 
     164ENDIF 
     165; 
     166if n_elements(pts_downward) EQ 1 then BEGIN  
     167    tempdeux = systime(1)       ; For key_performance =2 
     168; 
     169    triang = definetri(nx, ny) 
     170; 
     171    IF testvar(var = key_performance) EQ 2 THEN $ 
     172      print, 'temps triangule: definetri', systime(1)-tempdeux 
     173    coinmont = -1 
     174    coindesc = -1 
     175ENDIF ELSE BEGIN  
     176    tempdeux = systime(1)       ; For key_performance =2 
     177    pts_downward = pts_downward[1:n_elements(pts_downward)-1] 
     178    pts_downward = pts_downward[uniq(pts_downward, sort(pts_downward))] 
     179; None rectangle can have an element of the last column or of the  
     180; last line as left bottom corner. 
     181; so we have to remove these points if they has been selected in pts_downward. 
     182    derniere_colonne = (lindgen(ny)+1)*nx-1  
     183    derniere_ligne = lindgen(nx)+(ny-1)*nx  
     184    pts_downward =different(pts_downward,derniere_colonne ) 
     185    pts_downward =different(pts_downward,derniere_ligne ) 
     186    if (NOT keyword_set(basic)) OR keyword_set(coinmonte) OR keyword_set(coindescend) then begin 
     187        if coinmont[0] NE -1 then begin 
    191188            coinmont =different(coinmont,derniere_colonne ) 
    192189            coinmont =different(coinmont,derniere_ligne ) 
    193          endif 
    194          if coindesc[0] NE -1 then begin 
     190        endif 
     191        if coindesc[0] NE -1 then begin 
    195192            coindesc =different(coindesc,derniere_colonne ) 
    196193            coindesc =different(coindesc,derniere_ligne ) 
    197          endif 
    198       ENDIF ELSE BEGIN  
    199          coinmont = -1 
    200          coindesc = -1 
    201       ENDELSE  
    202       IF testvar(var = key_performance) EQ 2 THEN $ 
    203        print, 'temps triangule: menage ds pts_downward coinmont et coindesc', systime(1)-tempdeux 
    204 ; 
    205       tempdeux = systime(1)     ; pour key_performance =2 
    206       if  pts_downward[0] EQ -1 then triang = definetri(nx, ny) $ 
    207       ELSE triang = definetri(nx, ny, pts_downward) 
    208       IF testvar(var = key_performance) EQ 2 THEN $ 
    209        print, 'temps triangule: definetri', systime(1)-tempdeux 
    210    ENDELSE  
    211 ;------------------------------------------------------------ 
    212 ; on vire les triangles qui ne contiennent que des points terre 
    213 ;------------------------------------------------------------ 
    214 ; 
    215 ;  tres bonne idee qui ne marche pas encore a 200% avec IDL 5.2 
    216 ;  ca devrait aller mieux dans les prochaines versions d''IDL... 
    217 ; 
    218    if (NOT keyword_set(basic)) AND (NOT keyword_set(keep_cont)) then begin 
    219       tempdeux = systime(1)     ; pour key_performance =2 
    220 ; on enleve les rectangles qui sont entierement dans la terre 
    221       recdsterre = where((1-msk)*(1-shift(msk, -1, 0))*(1-shift(msk, 0, -1))*(1-shift(msk, -1, -1)) EQ 1) 
    222       IF testvar(var = key_performance) EQ 2 THEN $ 
    223        print, 'temps triangule: tous les recdsterre', systime(1)-tempdeux 
    224  
    225 ; en attendant une version qui marche parfaitement, on est contraint 
    226 ; de faire un nouveau tri: 
    227 ; il ne faut pas enlever les rectangles qui n''ont qu''un sommet en 
    228 ; commun. 
     194        endif 
     195    ENDIF ELSE BEGIN  
     196        coinmont = -1 
     197        coindesc = -1 
     198    ENDELSE  
     199    IF testvar(var = key_performance) EQ 2 THEN $ 
     200      print, 'temps triangule: menage ds pts_downward coinmont et coindesc', systime(1)-tempdeux 
     201; 
     202    tempdeux = systime(1)       ; For key_performance =2 
     203    if  pts_downward[0] EQ -1 then triang = definetri(nx, ny) $ 
     204    ELSE triang = definetri(nx, ny, pts_downward) 
     205    IF testvar(var = key_performance) EQ 2 THEN $ 
     206      print, 'temps triangule: definetri', systime(1)-tempdeux 
     207ENDELSE  
     208;------------------------------------------------------------ 
     209; We delete land points which only contain land points. 
     210;------------------------------------------------------------ 
     211; 
     212 
     213if (NOT keyword_set(basic)) AND (NOT keyword_set(keep_cont)) then begin 
     214    tempdeux = systime(1)       ; For key_performance =2 
     215; We delete rectangles which are entirely in the land. 
     216    recdsterre = where((1-msk)*(1-shift(msk, -1, 0))*(1-shift(msk, 0, -1))*(1-shift(msk, -1, -1)) EQ 1) 
     217    IF testvar(var = key_performance) EQ 2 THEN $ 
     218      print, 'temps triangule: tous les recdsterre', systime(1)-tempdeux 
     219 
     220; We do an other sort : 
     221; We have to do not remove rectangles which only have one common summit. 
    229222; t1 = systime(1) 
    230       indice = intarr(nx, ny) 
    231       trimask = intarr(nx, ny) 
    232       trimask[0:nx-2, 0:ny-2] = 1 
    233       IF recdsterre[0] NE -1 then BEGIN  
    234          tempdeux = systime(1)  ; pour key_performance =2 
    235          indice[recdsterre] = 1 
    236 ;      if NOT keyword_set(basic) then begin 
    237          vire1 = 0 
    238          vire2 = 0 
    239          while (vire1[0] NE -1 OR vire2[0] NE -1) ne 0 do begin 
    240 ; vire sont les rectangles qu''il faut retirer de recsterre (en fait 
    241 ; qu''il faut garder bien qu''ils soient entirement dans la terre)   
    242             vire1 = where( (indice*shift(indice, -1, -1) $ 
    243                             *(1-shift(indice, 0, -1))*(1-shift(indice, -1, 0))*trimask) EQ 1) 
    244             if vire1[0] NE -1 THEN BEGIN  
    245                indice[vire1] = 0 
     223    indice = intarr(nx, ny) 
     224    trimask = intarr(nx, ny) 
     225    trimask[0:nx-2, 0:ny-2] = 1 
     226    IF recdsterre[0] NE -1 then BEGIN  
     227        tempdeux = systime(1)   ; For key_performance =2 
     228        indice[recdsterre] = 1 
     229        if NOT keyword_set(basic) then begin 
     230            vire1 = 0 
     231            vire2 = 0 
     232            while (vire1[0] NE -1 OR vire2[0] NE -1) ne 0 do begin 
     233; Delete rectangles we have to remove from recsterre (in fact those we have  
     234; to keep although they ar eentirely in the land). 
     235                vire1 = where( (indice*shift(indice, -1, -1) $ 
     236                                *(1-shift(indice, 0, -1))*(1-shift(indice, -1, 0))*trimask) EQ 1) 
     237                if vire1[0] NE -1 THEN BEGIN  
     238                    indice[vire1] = 0 
    246239;               indice[vire1+nx+1] = 0 
    247             endif 
    248              
    249             vire2 = where( ((1-indice)*(1-shift(indice, -1, -1)) $ 
    250                             *shift(indice, 0, -1)*shift(indice, -1, 0)*trimask) EQ 1) 
    251             if vire2[0] NE -1 THEN BEGIN  
    252                indice[vire2+1] = 0 
     240                endif 
     241                 
     242                vire2 = where( ((1-indice)*(1-shift(indice, -1, -1)) $ 
     243                                *shift(indice, 0, -1)*shift(indice, -1, 0)*trimask) EQ 1) 
     244                if vire2[0] NE -1 THEN BEGIN  
     245                    indice[vire2+1] = 0 
    253246;               indice[vire2+nx] = 0 
    254             endif 
    255          endwhile 
    256          IF testvar(var = key_performance) EQ 2 THEN $ 
    257           print, 'temps triangule: trier les recdsterre', systime(1)-tempdeux 
    258 ;      endif 
    259          indice[*, ny-1] = 1    ; la deriere colonne te la derniere ligne 
    260          indice[nx-1, *] = 1    ; ne peuvent definir de rectangle. 
    261 ; 
    262          tempdeux = systime(1)  ; pour key_performance =2 
    263          recgarde = where(indice EQ 0) 
    264 ; on recupere les numeros des triangles que l'' on va garder 
    265          trigarde = 2*[recgarde-recgarde/nx] 
    266          trigarde = transpose(temporary(trigarde)) 
    267          trigarde = [trigarde, trigarde+1] 
     247                endif 
     248            endwhile 
     249            IF testvar(var = key_performance) EQ 2 THEN $ 
     250              print, 'temps triangule: trier les recdsterre', systime(1)-tempdeux 
     251        endif 
     252        indice[*, ny-1] = 1     ; The last column and the last line 
     253        indice[nx-1, *] = 1     ; can not define any rectangle. 
     254; 
     255        tempdeux = systime(1)   ; For key_performance =2 
     256        recgarde = where(indice EQ 0) 
     257; We recuperate numbers of triangles we will keep. 
     258        trigarde = 2*[recgarde-recgarde/nx] 
     259        trigarde = transpose(temporary(trigarde)) 
     260        trigarde = [trigarde, trigarde+1] 
    268261;  
    269          triang = triang[*, temporary(trigarde[*])] 
    270          IF testvar(var = key_performance) EQ 2 THEN $ 
     262        triang = triang[*, temporary(trigarde[*])] 
     263        IF testvar(var = key_performance) EQ 2 THEN $ 
    271264          print, 'temps triangule: virer les triangle de la liste', systime(1)-tempdeux 
    272       endif 
    273    endif 
     265    endif 
     266endif 
    274267; print, 'temps tri triangles', systime(1)-t1  
    275268;------------------------------------------------------------ 
    276 ; quand key_periodic eq 1, triang est une liste d''indice d'un 
    277 ; tableau qui a une colonne de trop. 
    278 ; il faut ramener ca a la matrice initiale en mettant les indivces de 
    279 ; la derniere colonne egaux a ceux de la derniere colonne... 
    280 ;------------------------------------------------------------ 
    281    tempdeux = systime(1)        ; pour key_performance =2 
    282    if keyword_set(key_periodic)*(nx-1 EQ jpi) $ 
    283     AND NOT keyword_set(basic) then BEGIN  
    284       indicey = triang/nx 
    285       indicex = triang-indicey*nx 
    286       nx = nx-1 
    287       liste = where(indicex EQ nx) 
    288       if liste[0] NE -1 then indicex[liste] = 0 
    289       triang = indicex+nx*indicey 
    290       nx = nx+1 
    291       if coinmont[0] NE -1 then begin 
    292          indicey = coinmont/nx 
    293          indicex = coinmont-indicey*nx 
    294          nx = nx-1 
    295          liste = where(indicex EQ nx) 
    296          if liste[0] NE -1 THEN indicex[liste] = 0 
    297          coinmont = indicex+nx*indicey 
    298          nx = nx+1 
    299       endif 
    300       if coindesc[0] NE -1 then begin 
    301          indicey = coindesc/nx 
    302          indicex = coindesc-indicey*nx 
    303          nx = nx-1 
    304          liste = where(indicex EQ nx) 
    305          if liste[0] NE -1 THEN indicex[liste] = 0 
    306          coindesc = indicex+nx*indicey 
    307          nx = nx+1 
    308       endif 
    309    endif 
    310    IF testvar(var = key_performance) EQ 2 THEN $ 
    311     print, 'temps triangule: finitions', systime(1)-tempdeux 
    312  
    313 ;------------------------------------------------------------ 
    314    if keyword_set(coinmonte) THEN coinmonte = coinmont ELSE twin_corners_up = coinmont 
    315    if keyword_set(coindescend) THEN coindescend = coindesc ELSE twin_corners_dn = coindesc 
    316 ;------------------------------------------------------------ 
    317   IF NOT keyword_set(key_forgetold) THEN BEGIN 
     269; When key_periodic equal 1, triang is a list of indexes's array which  
     270; have a surplus column. 
     271; We have to put it back to the initial matrix by putting indexes of  
     272; the last column equal to these of the last column... 
     273;------------------------------------------------------------ 
     274tempdeux = systime(1)           ; For key_performance =2 
     275if keyword_set(key_periodic)*(nx-1 EQ jpi) $ 
     276  AND NOT keyword_set(basic) then BEGIN  
     277    indicey = triang/nx 
     278    indicex = triang-indicey*nx 
     279    nx = nx-1 
     280    liste = where(indicex EQ nx) 
     281    if liste[0] NE -1 then indicex[liste] = 0 
     282    triang = indicex+nx*indicey 
     283    nx = nx+1 
     284    if coinmont[0] NE -1 then begin 
     285        indicey = coinmont/nx 
     286        indicex = coinmont-indicey*nx 
     287        nx = nx-1 
     288        liste = where(indicex EQ nx) 
     289        if liste[0] NE -1 THEN indicex[liste] = 0 
     290        coinmont = indicex+nx*indicey 
     291        nx = nx+1 
     292    endif 
     293    if coindesc[0] NE -1 then begin 
     294        indicey = coindesc/nx 
     295        indicex = coindesc-indicey*nx 
     296        nx = nx-1 
     297        liste = where(indicex EQ nx) 
     298        if liste[0] NE -1 THEN indicex[liste] = 0 
     299        coindesc = indicex+nx*indicey 
     300        nx = nx+1 
     301    endif 
     302endif 
     303IF testvar(var = key_performance) EQ 2 THEN $ 
     304  print, 'temps triangule: finitions', systime(1)-tempdeux 
     305 
     306;------------------------------------------------------------ 
     307if keyword_set(coinmonte) THEN coinmonte = coinmont ELSE twin_corners_up = coinmont 
     308if keyword_set(coindescend) THEN coindescend = coindesc ELSE twin_corners_dn = coindesc 
     309;------------------------------------------------------------ 
     310IF NOT keyword_set(key_forgetold) THEN BEGIN 
    318311   @updateold 
    319  ENDIF  
    320  
    321    IF keyword_set(key_performance) THEN print, 'temps triangule', systime(1)-tempsun  
    322  
    323    return, triang 
     312ENDIF  
     313 
     314IF keyword_set(key_performance) THEN print, 'temps triangule', systime(1)-tempsun  
     315 
     316return, triang 
    324317 
    325318END  
  • trunk/SRC/ToBeReviewed/TRIANGULATION/triangule_e.pro

    r134 r150  
    33;------------------------------------------------------------ 
    44;+ 
    5 ; NAME:triangule_e 
    6 ; 
    7 ; PURPOSE:buid the triangulation for a E-grid type 
    8 ; 
    9 ; CATEGORY: 
    10 ; 
    11 ; CALLING SEQUENCE: 
    12 ;  
    13 ; INPUTS: 
    14 ; 
    15 ; KEYWORD PARAMETERS: 
    16 ; 
    17 ; OUTPUTS: 
    18 ; 
    19 ; COMMON BLOCKS:common.pro 
    20 ; 
    21 ; SIDE EFFECTS: 
    22 ; 
    23 ; RESTRICTIONS: 
    24 ; 
    25 ; EXAMPLE: 
    26 ; 
    27 ; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr) 
     5; @file_comments 
     6; Buid the triangulation for a E-grid type 
     7; 
     8; @categories 
     9; graphic 
     10;  
     11; @param MASKENTREE {in}{optional} 
     12; It is a 2d array which will serve to mask the field we will trace after with CONTOUR,  
     13; ...TRIANGULATION=triangule(mask) 
     14; If this argument is not specified, the function use tmask 
     15;  
     16; @keyword BASIC 
     17; Specify that the mask is on a basic grid (use the triangulation for vertical cuts and hovmoellers) 
     18;  
     19; @keyword COINMONTE 
     20; It is an array. To obtain the array of "ascending land corner" to be treated with  
     21; completecointerre.pro in the variable array instead of make it pass by the global  
     22; variable twin_corners_up. 
     23; 
     24; @keyword COINDESCEND 
     25; It is an array. See COINMONTE 
     26; 
     27; @keyword SHIFTED 
     28;  
     29; @uses 
     30; common.pro 
     31;  
     32; @history 
     33; Sebastien Masson (smasson@lodyc.jussieu.fr) 
    2834;                      june 2001 
     35;  
     36; @version  
     37; $Id$ 
     38;  
     39; @todo 
     40; seb L.152->153 je ne pense pas que ce soit ce que tu voulais dire mais  
     41; c'est la traduction de ce qu'il y avait écrit. Correction si besoin. 
    2942;- 
    3043;------------------------------------------------------------ 
     
    4255  ENDIF 
    4356;--------------------------------------------------------- 
    44    tempsun = systime(1)         ; pour key_performance 
    45 ;------------------------------------------------------------ 
    46 ; le masque est donne ou il faut prendre tmask? 
     57   tempsun = systime(1)         ; For key_performance 
     58;------------------------------------------------------------ 
     59; Is the mask given or do we have to take tmask? 
    4760;------------------------------------------------------------ 
    4861; 
     
    138151; ; 
    139152;------------------------------------------------------------ 
    140 ; quand key_periodic eq 1, triang est une liste d''indice d'un 
    141 ; tableau qui a une colonne de trop. 
    142 ; il faut ramener ca a la matrice initiale en mettant les indivces de 
    143 ; la derniere colonne egaux a ceux de la derniere colonne... 
    144 ;------------------------------------------------------------ 
    145    tempdeux = systime(1)        ; pour key_performance =2 
     153; When key_periodic equal 1, triang is a list of indexes's array which  
     154; have a surplus column. 
     155; We have to put it back to the initial matrix by putting indexes of  
     156; the last column equal to these of the last column... 
     157;------------------------------------------------------------ 
     158   tempdeux = systime(1)        ; For key_performance =2 
    146159   if keyword_set(key_periodic)*(nx-1 EQ jpi) $ 
    147160    AND NOT keyword_set(basic) then BEGIN  
     
    176189 
    177190;------------------------------------------------------------ 
    178 ;    if arg_present(coinmonte) THEN coinmonte = coinmont ELSE twin_corners_up = coinmont 
    179 ;    if arg_present(coindescend) THEN coindescend = coindesc ELSE twin_corners_dn = coindesc 
    180 ; 
    181 ;   IF NOT keyword_set(key_forgetold) THEN BEGIN 
    182 ;    @updateold 
    183 ;   ENDIF 
     191    if arg_present(coinmonte) THEN coinmonte = coinmont ELSE twin_corners_up = coinmont 
     192    if arg_present(coindescend) THEN coindescend = coindesc ELSE twin_corners_dn = coindesc 
     193 
     194   IF NOT keyword_set(key_forgetold) THEN BEGIN 
     195    @updateold 
     196   ENDIF 
    184197;------------------------------------------------------------ 
    185198 
Note: See TracChangeset for help on using the changeset viewer.