source: trunk/SRC/ToBeReviewed/TRIANGULATION/tracemask.pro @ 76

;------------------------------------------------------------
;------------------------------------------------------------
;------------------------------------------------------------
;+
; NAME:tracemask
;
; PURPOSE:dessiner des contour d''un mask
;
; CATEGORY:plus simple que tracecote, car ne s''occuppe pas du type de
; projection et de la periodicite de la grille
;
; CALLING SEQUENCE: tracemask, maskentree, xentree, yentree
; 
; INPUTS:maskentree, xentree, yentree tableaux 2d specifiant le mask
; et ses coordonees en longitude te latitude.
;
; KEYWORD PARAMETERS:
;
;        COAST_COLOR: the color of the coastline.
;                     defaut value is 0 => black
;
;        COAST_THICK: l''epaisseur du trait pour tracer les
;        continents. par defaut c''est 1.
;
; OUTPUTS: none
;
; COMMON BLOCKS:common.pro
;
; SIDE EFFECTS:
;
; RESTRICTIONS:
;
; EXAMPLE:
;
; MODIFICATION HISTORY:Sebastien Masson (smasson@lodyc.jussieu.fr)
;
;-
;------------------------------------------------------------
;------------------------------------------------------------
;------------------------------------------------------------
PRO tracemask, maskentree, xin, yin, COAST_COLOR = coast_color, COAST_THICK = coast_thick, OVERPLOT = overplot, _extra = ex
;
   if keyword_set(overplot) then return
;---------------------------------------------------------
@cm_general
  IF NOT keyword_set(key_forgetold) THEN BEGIN
@updatekwd
  ENDIF
;--------------------------------------------------------- 
   tempsun = systime(1)         ; pour key_performance
; on s''afranchit des problemes de bord:
   tempdeux = systime(1)        ; pour key_performance =2
   tailleentree = size(maskentree)
   nx = tailleentree[1]+1
   ny = tailleentree[2]+1
; we check the input axis
  IF n_elements(xin) EQ 0 THEN xentree = findgen(nx-1) ELSE xentree = xin
  IF (size(xentree))[0] EQ 1 THEN xentree = xentree#replicate(1,ny-1)
  IF n_elements(yin) EQ 0 THEN yentree = findgen(ny-1) ELSE yentree = yin
  IF (size(yentree))[0] EQ 1 THEN yentree = replicate(1,nx-1)#yentree
; on agrandi le mask de une colonne a gauche et de une colonne en bas
   mask = intarr(tailleentree[1]+1, tailleentree[2]+1)
   mask[1:tailleentree[1], 1:tailleentree[2]] = maskentree
; les 2 premieres colonnes sont identiques
   mask[0, 1:tailleentree[2]] = maskentree[0, *]
; les 2 premieres lignes sont identiques
   mask[1:tailleentree[1], 0] = maskentree[*, 0]
; on calcul la position suivant x des points qui seviront a tracer le
; masque. ils sont situes entre chaque points du masque, sauf pour la
; derniere colonne que l''on ne peut pas calculer et que l''on met
; donc a max(!x.range)
   xrange = !x.range[sort(!x.range)] ; si reverse_x est utilise!
   xentree = .5*(xentree+shift(xentree, -1, 0))
   IF not keyword_set(overplot) THEN xentree[nx-2, *] = xrange[1] $
   ELSE xentree[nx-2, *] = xentree[nx-3, *]
; on seuil
   xentree = xrange[0] > xentree < xrange[1]
; on agrandit le tableau
   xf = fltarr(nx, ny)
   xf[1:nx-1, 1:ny-1] = xentree
   IF not keyword_set(overplot) THEN xf[0, *] = xrange[0] $
   ELSE xf[0, *] = xf[1, *]
   xf[1:nx-1, 0] = xentree[*, 0]
;
   yinverse = yentree[0, 0] GT yentree[0, ny-2]
   yrange = !y.range[sort(!y.range)]
   yentree = .5*(yentree+shift(yentree, 0, -1))
   IF not keyword_set(overplot) THEN BEGIN 
      if yinverse then yentree[*, ny-2] = yrange[0] ELSE yentree[*, ny-2] = yrange[1]
   ENDIF ELSE yentree[*, ny-2] = yentree[*, ny-3]
   yentree = yrange[0] > yentree < yrange[1]
   yf = fltarr(nx, ny)
   yf[1:nx-1, 1:ny-1] = yentree
   yf[0, 1:ny-1] = yentree[0, *]
   IF not keyword_set(overplot) THEN BEGIN 
      if yinverse then yf[*, 0] = yrange[1] ELSE yf[*, 0] = yrange[0]
   ENDIF ELSE yentree[*, 0] = yentree[*, 1]
;
   IF testvar(var = key_performance) EQ 2 THEN $
    print, 'temps tracemask: determination du mask et des ses coordonnes', systime(1)-tempdeux
;
; on trace les segments verticaux:
;
   tempdeux = systime(1)        ; pour key_performance =2
   liste = where(mask+shift(mask, -1, 0) EQ 1)
   IF liste[0] NE -1 THEN BEGIN
; on recupere lx et ly qui sont les indices ds un tableau 2d des
; points donnes par liste
      ly = liste/nx & lx = temporary(liste)-nx*ly
      indice = where(ly NE 0)   ; on ne prend pas les points concernant 
; la premiere ligne car ds ce cas le pt j-1 n''est pas definit
      if indice[0] NE -1 then begin
         lx = lx[indice] & ly = ly[temporary(indice)]
         IF testvar(var = key_performance) EQ 2 THEN $
          print, 'temps tracemask: liste traits verticaux', systime(1)-tempdeux
         tempdeux = systime(1)  ; pour key_performance =2
; boucle sur les points concernes et trace du segment
; rq: on utilise plots au lieu de plot car plots est bcp plus rapide.
         for pt = 0L, n_elements(lx)-1 do BEGIN 
            i = lx[pt] & j = ly[pt]
            plots, [xf[i, j-1], xf[i, j]], [yf[i, j-1], yf[i, j]] $
              , color = coast_color, thick = coast_thick, _extra = ex
            if pt LT 5 then begin
            endif
         endfor
         IF testvar(var = key_performance) EQ 2 THEN $
          print, 'temps tracemask: trace traits verticaux', systime(1)-tempdeux
      endif
   ENDIF
;
; on trace les segments horizontaux:
;
   tempdeux = systime(1)        ; pour key_performance =2
   liste = where(mask+shift(mask, 0, -1) EQ 1)
   IF liste[0] NE -1 THEN BEGIN
      ly = liste/nx & lx = temporary(liste)-nx*ly
      indice = where(lx NE 0)   ; on ne prend pas les points de la  premiere colonne
      if indice[0] EQ -1 then return
      lx = lx[indice] & ly = ly[temporary(indice)]
      IF testvar(var = key_performance) EQ 2 THEN $
       print, 'temps tracemask: liste traits horizontaux', systime(1)-tempdeux
      tempdeux = systime(1)     ; pour key_performance =2
      for pt = 0L, n_elements(lx)-1 do BEGIN 
         i = lx[pt] & j = ly[pt]
         plots, [xf[i-1, j], xf[i, j]], [yf[i-1, j], yf[i, j]] $
           , color = coast_color, thick = coast_thick, _extra = ex
      endfor
      IF testvar(var = key_performance) EQ 2 THEN $
       print, 'temps tracemask: trace traits horizontaux', systime(1)-tempdeux
   endif
;
;
   if keyword_set(key_performance) THEN print, 'temps tracemask', systime(1)-tempsun 
   
   return
end