;------------------------------------------------------------
;------------------------------------------------------------
;------------------------------------------------------------
;+
; NAME:
;
; PURPOSE:vire les tableaux qui ne doivent pas etre dessines grace a 2
; tests:
;     1) les coins du tableau doivent etre ds la fenetre
;     2) les clongeurs des cotes des triangfles exprimes en
;     coordonnees normalisesne doivent pas depasser une certaine
;     longueur seuil. 
;
; CATEGORY:
;
; CALLING SEQUENCE:
;
; INPUTS:
;
; KEYWORD PARAMETERS:
;
; OUTPUTS:
;
; COMMON BLOCKS:
;	common.pro
;
; SIDE EFFECTS:
;
; RESTRICTIONS:
;
; EXAMPLE:
;
; MODIFICATION HISTORY: Sebastien Masson (smasson@lodyc.jussieu.fr)
;                       20/2/99
;-
;------------------------------------------------------------
;------------------------------------------------------------
;------------------------------------------------------------
FUNCTION ciseauxtri, triang, glam, gphi, TOUT = tout, _EXTRA = ex
@common
   tempsun = systime(1)         ; pour key_performance
   taille = size(glam)
   nx = taille[1]
   ny = taille[2]

   tempdeux = systime(1)        ; pour key_performance =2
   z = convert_coord(glam(*),gphi(*),/data,/to_normal) 
   x = z(0, *)
   y = z(1, *)
   tempvar = SIZE(TEMPORARY(z)) ; delete z
   IF testvar(var = key_performance) EQ 2 THEN $
    print, 'temps ciseauxtri: convert_coord data to_normal', systime(1)-tempdeux
;
; attention, suivant la projection certains points x ou y peuvent
; devenir NaN (cf. points deriere la terre ds une projection
; orthographique) il faut dans ce cas enlever tous les triangles qui
; contiennent un de ces points
;
   if (!map.projection LE 7 AND !map.projection NE 0) $
    OR !map.projection EQ 14 OR !map.projection EQ 15 OR !map.projection EQ 18 then begin
      tempdeux = systime(1)     ; pour key_performance =2

      ind = where(( finite((x*y)[triang[0, *]])$ ; points NaN ds la premiere colonne
                    *finite((x*y)[triang[1, *]]) $ ; points NaN ds la 2eme colonne
                    *finite((x*y)[triang[2, *]]) ) EQ 1) ; points NaN ds la 3eme colonne

      if ind[0] NE -1 then triang = triang[*, ind] ELSE return,  -1
      IF testvar(var = key_performance) EQ 2 THEN $
       print, 'temps ciseauxtri: recherche points a NAN', systime(1)-tempdeux
   endif
;

   seuil = 5.
;
; maintenant on vire les triangles dont un des cotes a une taille
; superieure a 1/seuil du domaine reserve au dessin.
;

   if keyword_set(key_periodique)  then begin
      tempdeux = systime(1)     ; pour key_performance =2
      ind=where( (abs(x[triang[1,*]]-x[triang[0,*]]) LE (!p.position[2]-!p.position[0])/seuil) $
                 AND (abs(x[triang[2,*]]-x[triang[1,*]]) LE (!p.position[2]-!p.position[0])/seuil) $
                 AND (abs(x[triang[0,*]]-x[triang[2,*]]) LE (!p.position[2]-!p.position[0])/seuil) $
;
      AND (abs(y[triang[0,*]]-y[triang[1,*]]) LE (!p.position[3]-!p.position[1])/seuil) $
       AND (abs(y[triang[1,*]]-y[triang[2,*]]) LE (!p.position[3]-!p.position[1])/seuil) $ 
       AND (abs(y[triang[2,*]]-y[triang[0,*]]) LE (!p.position[3]-!p.position[1])/seuil) ) 

      IF testvar(var = key_performance) EQ 2 THEN $
       print, 'temps ciseauxtri: trouver les triangles trop grands', systime(1)-tempdeux
;
      tempdeux = systime(1)     ; pour key_performance =2
      if ind[0] NE -1 then triang = triang[*, ind] ELSE return,  -1
      IF testvar(var = key_performance) EQ 2 THEN $
       print, 'temps ciseauxtri: virer les triangles trop grands', systime(1)-tempdeux
   endif
;
; on supprime les triangles qui sont un peu trop hors du cadre
; valable qd /TOUT est active
;
;    if keyword_set(key_irregular) then begin

;       tempdeux = systime(1)     ; pour key_performance =2
;       seuil = 0
;       seuil = .1*max([!p.position[2]-!p.position[0],!p.position[3]-!p.position[1]])
;       ind=where(x[triang[0,*]] GE !p.position[0]-seuil AND x[triang[0,*]] LE !p.position[2]+seuil $
;                 AND x[triang[1,*]] GE !p.position[0]-seuil AND x[triang[1,*]] LE !p.position[2]+seuil $
;                 AND x[triang[2,*]] GE !p.position[0]-seuil AND x[triang[2,*]] LE !p.position[2]+seuil $
; ;
;       AND y[triang[0,*]] GE !p.position[1]-seuil AND y[triang[0,*]] LE !p.position[3]+seuil $
;        AND y[triang[1,*]] GE !p.position[1]-seuil AND y[triang[1,*]] LE !p.position[3]+seuil $
;        AND y[triang[2,*]] GE !p.position[1]-seuil AND y[triang[2,*]] LE !p.position[3]+seuil )
;       IF testvar(var = key_performance) EQ 2 THEN $
;        print, 'temps ciseauxtri: trouver les triangles hors du cadre', systime(1)-tempdeux
; ;
;       tempdeux = systime(1)     ; pour key_performance =2
;       if ind[0] NE -1 then triang = triang[*, ind] ELSE return,  -1
;       IF testvar(var = key_performance) EQ 2 THEN $
;        print, 'temps ciseauxtri: virer les triangles hors du cadre', systime(1)-tempdeux
;    ENDIF
;
;        dernier tri
;
   if n_elements(ind) GT 1 then BEGIN
;
; il faut verifier que les triangles que l''on garde ne
; forment pas une triangulation dans laquelle 2 triangles ont un
; sommet en commun sans avoir de cote de commun.
;
; on constitue la liste des rectangles que l''on souhaite garder (on
; garde un rectangle des qu''il y a un triangle dedans)
;
; dans definetri, on a construit les triangles tels que le premier et 
; le dernier sommets soient ceux de la diagonale du rectangle definit
; par le maillage.
; pour retouver de quel rectangle provient un triangle, on cherche le
; min de l''indice suivant x et suivant y de chaque triangle. Apres on
; repasse cette indissage suivant x et y en indicage suivant nx*ny
;
      tempdeux = systime(1)     ; pour key_performance =2
; indices en y des rectangles
      indytriang = (triang[0, *]/nx) < (triang[2, *]/nx)
; indices en x des rectangles
      indxtriang0 = triang[0, *]-nx*(triang[0, *]/nx)
      indxtriang2 = triang[2, *]-nx*(triang[2, *]/nx)
      indxmin = indxtriang0 < indxtriang2
; attention dans le cas ou la grille est periodique et ou on a tous
; les points suivant x, les triangles qui assurent la periodicite en x
; ont des indices suivant x qui sont 0 et nx-1. Ils appatienent au
; rectangles de la colonne nx-1 et non de la colonne 0
      if keyword_set(key_periodique) AND nx EQ jpi then BEGIN
      indxmax = indxtriang0 > indxtriang2
      indxtriang = indxmin + (nx-1)*(indxmin EQ 0 AND indxmax EQ (nx-1))
      ENDIF ELSE indxtriang = indxmin 
      listrect = nx*indytriang+indxtriang
      IF testvar(var = key_performance) EQ 2 THEN $
       print, 'temps ciseauxtri: liste des rectangles', systime(1)-tempdeux
;
; maintenant qu''on a cette liste on va s''assuter que l''on a pas de
; triangles qui n''ont qu''on sommet en commun.
;
      test = bytarr(nx, ny)
      test[listrect] = 1
      dejavire = 1b-test
;
      tempdeux = systime(1)     ; pour key_performance =2
      vire1 = 0
      vire2 = 0
      while (vire1[0] NE -1 OR vire2[0] NE -1) ne 0 do begin
         vire1 = where( (test*shift(test, -1, -1) $
                         *(1-shift(test, 0, -1))*(1-shift(test, -1, 0))) EQ 1)
         if vire1[0] NE -1 THEN test[vire1] = 0 ; on vire le rectangle 
         vire2 = where( ((1-test)*(1-shift(test, -1, -1)) $
                         *shift(test, 0, -1)*shift(test, -1, 0)) EQ 1)
; on vire le rectangle du dessus (meme indice x mais egale a 1)
         if vire2[0] NE -1 THEN test[vire2+nx] = 0 
      ENDWHILE
;stop
      test = test+temporary(dejavire)
;
      avirer = where(temporary(test) EQ 0)
      IF testvar(var = key_performance) EQ 2 THEN $
       print, 'temps ciseauxtri: determinationdes rectangles a virer', systime(1)-tempdeux
;
      if avirer[0] NE -1 then begin
      tempdeux = systime(1)     ; pour key_performance =2
      indnx = n_elements(listrect)
      indny = n_elements(avirer)
      ind = listrect[*]#replicate(1l, indny) 
      ind = ind EQ replicate(1, indnx)#avirer[*]
      if indny GT 1 then ind = total(ind, 2)
      ind = where(ind EQ 0)
      if ind[0] NE -1 then triang = triang[*, ind] ELSE return,  -1
      endif
      IF testvar(var = key_performance) EQ 2 THEN $
       print, 'temps ciseauxtri: derniere retouche de la triangulation', systime(1)-tempdeux
   endif
;    
;
   if keyword_set(key_performance) THEN print, 'temps ciseauxtri', systime(1)-tempsun 
   return,  triang
end