Index: /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/compute_sections.f90
===================================================================
--- /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/compute_sections.f90	(revision 2876)
+++ /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/compute_sections.f90	(revision 2877)
@@ -301,6 +301,8 @@
            DO ji=1,jpi  
               coord_t=COORD_SECTION(glamf(ji,jj),gphif(ji,jj))
-              zdistFirst = distance(coord_t,coordFirst)
-              zdistLast = distance(coord_t,coordLast)
+!cbr              zdistFirst = distance(coord_t,coordFirst)
+              zdistFirst = distance2(coord_t,coordFirst)
+!cbr              zdistLast = distance(coord_t,coordLast)
+              zdistLast = distance2(coord_t,coordLast)
               IF( zdistFirst .LT. zdistmesh .AND. zdistFirst .LT. zdistante )THEN
                  sec%listPoint(1) = POINT_SECTION(ji,jj)
@@ -391,5 +393,6 @@
 
         !initialize distnew value (with distance between section's extremities) 
-        zdistnew  = distance(coordFirst,coordLast,sec%ll_date_line) 
+!cbr        zdistnew  = distance(coordFirst,coordLast,sec%ll_date_line) 
+        zdistnew  = distance2(coordFirst,coordLast,sec%ll_date_line) 
         prevPoint  = POINT_SECTION(0,0)
         jseg       = 1
@@ -707,4 +710,5 @@
 
      !now compute new slopeSection with ij-coordinates of first and last point 
+     sec%slopeSection = 0  ! default value
      IF( sec%nb_point .ne. 0 )THEN
         IF (  sec%listPoint(sec%nb_point)%I .NE.  sec%listPoint(1)%I ) THEN
Index: /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/readcoordmesh.f90
===================================================================
--- /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/readcoordmesh.f90	(revision 2876)
+++ /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/readcoordmesh.f90	(revision 2877)
@@ -22,4 +22,5 @@
   !!
   !! ** Purpose :   Read a coordinate and a meshmask file in NetCDF format 
+  !!
   !!---------------------------------------------------------------------      
   PRINT*,'              '
@@ -54,7 +55,8 @@
   SUBROUTINE getdim(cdfile)
   !!----------------------------------------------------------------------
-  !!              ***  ROUTINE coord_mesh_read  ***
+  !!              ***  ROUTINE getdim  ***
   !!
-  !! ** Purpose :   Read a coordinate and a meshmask file in NetCDF format
+  !! ** Purpose :   get dimsensions of a netcdf file
+  !!
   !!----------------------------------------------------------------------
   !! * Arguments
@@ -110,4 +112,5 @@
   !!
   !! ** Purpose :   Read a coordinate and a meshmask file in NetCDF format
+  !!
   !!----------------------------------------------------------------------
   !! * Arguments
Index: /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/readsections.f90
===================================================================
--- /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/readsections.f90	(revision 2876)
+++ /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/readsections.f90	(revision 2877)
@@ -22,13 +22,7 @@
      !!         ***  ROUTINE read_list_sections ***
      !!
-     !! ** Purpose
+     !! ** Purpose: read ascii file 'list_sections.ascii' that contains
+     !!             sections description
      !!
-     !! ** Method
-     !!
-     !! ** Input
-     !!
-     !! ** Action
-     !!
-     !! History
      !!---------------------------------------------------------------------
      !! * arguments
@@ -45,7 +39,7 @@
      CHARACTER(len=9)    :: cdstrpond
      CHARACTER(len=110)  :: clname,cdsecname,cltmp
-     REAL,DIMENSION(nb_type_class):: zclass_value
-     TYPE(COORD_SECTION) :: coord_point1,coord_point2,coordTemp
-     TYPE(COORD_SECTION), DIMENSION(2)::coord_sec
+     REAL,DIMENSION(nb_type_class)     :: zclass_value
+     TYPE(COORD_SECTION)               :: coord_point1,coord_point2,coordTemp
+     TYPE(COORD_SECTION), DIMENSION(2) :: coord_sec
      !!---------------------------------------------------------------------
      PRINT*,'              '
Index: /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/sections_tools.f90
===================================================================
--- /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/sections_tools.f90	(revision 2876)
+++ /branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/TOOLS/SECTIONS_DIADCT/src/sections_tools.f90	(revision 2877)
@@ -13,13 +13,13 @@
 
    !! * Routine accessibility
-   PUBLIC pointToCoordF
-   PUBLIC distance
-   PUBLIC distance2
-   PUBLIC distance3
-   PUBLIC intersec
-   PUBLIC slope_coeff
-   PUBLIC qcksrt
-   PUBLIC write_debug   
-   PUBLIC file_open
+   PUBLIC pointToCoordF         ! define a point with geographical coordinates
+   PUBLIC distance              ! compute distance between 2 points   
+   PUBLIC distance2             ! compute distance between 2 points
+   PUBLIC distance3             ! compute distance between a point and a line
+   PUBLIC intersec              ! intersection between 2 lines
+   PUBLIC slope_coeff           ! slope coefficient of a line
+   PUBLIC qcksrt                ! organize a list of point
+   PUBLIC write_debug           ! write debug messages
+   PUBLIC file_open             ! open a file ( ascii, bin)
  
 CONTAINS
@@ -28,12 +28,17 @@
      !!---------------------------------------------------------------------
      !!         ***  FUNCTION pointToCoordF  ***
-     !! ** Purpose
+     !!           
+     !! ** Purpose: define a point with geographical coordinates
      !!
      !!---------------------------------------------------------------------
      !! * arguments
      TYPE(POINT_SECTION), INTENT(IN) ::  p
-     !!---------------------------------------------------------------------
+
+     !!---------------------------------------------------------------------
+
      pointToCoordF = COORD_SECTION(glamf(p%I,p%J),gphif(p%I,p%J))
+
      RETURN
+
   END FUNCTION pointToCoordF
 
@@ -41,20 +46,29 @@
      !!---------------------------------------------------------------------
      !!         ***  FUNCTION distance  ***
-     !! ** Purpose:compute distance between coordA and coordB
-     !!            We add 360] to coordB%long if the line (coordA,coordB)
-     !!            crosses the date-line.
+     !!
+     !! ** Purpose: Compute distance between coordA and coordB
+     !!             We add 360° to coordB%long if the line (coordA,coordB)
+     !!             crosses the date-line.
+     !!
      !!---------------------------------------------------------------------
      !! * arguments
      TYPE(COORD_SECTION), INTENT(IN) :: coordA,coordB
      LOGICAL,INTENT(IN),OPTIONAL     :: ld_date_line
+
      !! * Local declarations 
-     INTEGER ::idateline    
+     INTEGER ::idateline             ! local integer
+
      !-----------------------------------------------
+     !
      idateline=0
      IF( PRESENT( ld_date_line))THEN
        IF( ld_date_line )idateline=1
      ENDIF
+     !
+     !function output
      distance = (sqrt((coordA%lon-coordB%lon-360*idateline )**2 + (coordA%lat-coordB%lat )**2))
+     !
      RETURN
+     !
   END FUNCTION distance
 
@@ -62,24 +76,29 @@
      !!---------------------------------------------------------------------
      !!         ***  FUNCTION distance  ***
+     !!
      !! ** Purpose:compute distance between coordA and coordB
      !!            We add 360 to coordB%long if the line (coordA,coordB)
      !!            crosses the date-line.
+     !!
      !!---------------------------------------------------------------------
      !! * arguments
      TYPE(COORD_SECTION), INTENT(IN) :: coordA,coordB
      LOGICAL,INTENT(IN),OPTIONAL     :: ld_date_line
+
      !! * Local declarations 
-     REAL(wp) :: rad, zrayon
-     REAL(wp) :: rpi = 3.141592653589793
-     REAL(wp) :: lam1,lam2,phi1,phi2
+     REAL(wp) :: zrad, zrayon
+     REAL(wp) :: zpi = 3.141592653589793
+     REAL(wp) :: zlam1, zlam2, zphi1, zphi2
+
      !-----------------------------------------------
-     zrayon=6367000.0
-     rad=rpi/180.
-     lam1=coordA%lon*rad
-     lam2=coordB%lon*rad
-     phi1=coordA%lat*rad
-     phi2=coordB%lat*rad
-     !
-     distance2=0.001*2.*zrayon *ASIN(SQRT(SIN((phi2-phi1)/2.0)**2.0+COS(phi1)*COS(phi2)*SIN((lam2-lam1)/2.0)**2.0) )
+     zrayon = 6367000.0
+     zrad   = zpi/180.
+     zlam1  = coordA%lon*zrad
+     zlam2  = coordB%lon*zrad
+     zphi1  = coordA%lat*zrad
+     zphi2  = coordB%lat*zrad
+     !
+     !function output
+     distance2=0.001*2.*zrayon *ASIN(SQRT(SIN((zphi2-zphi1)/2.0)**2.0+COS(zphi1)*COS(zphi2)*SIN((zlam2-zlam1)/2.0)**2.0) )
      !
      RETURN
@@ -89,5 +108,8 @@
      !!---------------------------------------------------------------------
      !!         ***  FUNCTION cosinus  ***
-     !!         théorème d'Al-Kashi
+     !!
+     !! ** Purpose: compute distance between a point B and a line (AC)
+     !!
+     !! ** Methode: use Al-Kashi theorem
      !!
      !!              B                   A first point of the section
@@ -105,7 +127,9 @@
      TYPE(COORD_SECTION), INTENT(IN) :: coordA,coordB,coordC
      LOGICAL                         :: ll_debug
+
      !! * Local declarations
      REAL(wp)                        :: ztmp
-     REAL(wp)                        :: za,zb,zc,angle
+     REAL(wp)                        :: za, zb, zc, zangle
+
      !-----------------------------------------------
      za = SQRT( (coordB%lon-coordC%lon)**2 + (coordB%lat-coordC%lat)**2 )
@@ -114,18 +138,16 @@
      !
      IF( za /= 0. .AND. zb /= 0. )THEN
-        ztmp=( za**2 + zb**2  - zc**2 ) / ( 2.0*za*zb )
-        ztmp=MIN(ztmp,1.00_wp)
-        !angle=ABS(acos( ( za**2 + zb**2  - zc**2 ) / ( 2.0*za*zb ) ) )
-        IF( ztmp==1.00_wp )THEN 
-             angle=0.0 
-        ELSE 
-             angle=ABS(acos(ztmp)) 
+        ztmp = ( za**2 + zb**2  - zc**2 ) / ( 2.0*za*zb )
+        ztmp = MIN(ztmp,1.00_wp)
+        IF( ztmp==1.00_wp )THEN ; zangle = 0.0 
+        ELSE                    ; zangle = ABS(acos(ztmp)) 
         ENDIF
      ELSE
-        PRINT*,'za zb =',za,zb ; STOP
-     ENDIF
-     !
-     distance3=za*ABS(sin(angle)) 
-     !
+        PRINT*,'You should not have this situation: za zb =',za,zb ; STOP
+     ENDIF
+     !
+     distance3 = za*ABS(sin(zangle)) 
+     !
+     !function output
      RETURN
   END FUNCTION distance3
@@ -146,7 +168,9 @@
      TYPE(COORD_SECTION), INTENT(IN) :: coordA, coordB
      LOGICAL,INTENT(IN),OPTIONAL     :: ld_dateline
+
      !! * Local declarations
      REAL(wp)                        :: zcoeff
      INTEGER                         :: idateline
+
      !!---------------------------------------------------------
      !initialization
@@ -194,4 +218,5 @@
      TYPE(SECTION), INTENT(IN)        :: sec
      TYPE(COORD_SECTION)  , INTENT(IN):: coord_a, coord_b
+
      !! * Local declarations
      TYPE(COORD_SECTION) :: coordInter
@@ -200,4 +225,5 @@
                             zX,zY     ! coordinates of intersection
      LOGICAL             :: ll_date_line=.FALSE. !segment [a,b] crosses the date-line ?
+
      !----------------------------------------------------------------------------
 
@@ -391,114 +417,104 @@
      !!         ***  SUBROUTINE qcksrt  ***
      !!
-     !! ** Purpose
-     !!
-     !! ** Method
-     !!
-     !! ** Input
-     !!
-     !! ** Action
-     !!
-     !! History
-     !!---------------------------------------------------------------------
-     !! * arguments
-     INTEGER,INTENT(IN)                 :: n
-     REAL(wp),DIMENSION(n),INTENT(INOUT):: arr1,arr2
+     !! ** Purpose : organize a list point by latitude , and after by longitude
+     !!
+     !!     History
+     !!---------------------------------------------------------------------
+     !! * arguments
+     INTEGER,INTENT(IN)                 :: n           ! number of points
+     REAL(wp),DIMENSION(n),INTENT(INOUT):: arr1,arr2   ! longitude and latitude
 
      !! * Local declarations
-     REAL(wp)                           :: fx,a,b
-     INTEGER,PARAMETER                  :: m=7,nstack=500
-     REAL(wp),PARAMETER                 :: fm=7875.,fa=211.,fc=1663.,fmi=1./fm
-     INTEGER                            :: l,ir,j,i,iq,jstack
-     INTEGER,DIMENSION(nstack)          :: istack(nstack)
+     INTEGER , PARAMETER                :: m = 7, nstack = 500
+     REAL(wp), PARAMETER                :: fm=7875.,fa=211.,fc=1663.,fmi=1./fm
+     INTEGER                            :: il, ir, i, iq, ist    !local integers
+     INTEGER                            :: jj                    !loop indice
+     REAL(wp)                           :: zfx, za, zb           !local real
+     INTEGER,DIMENSION(nstack)          :: istack(nstack)        !temp array
      !---------------------------------------------------
-      jstack=0
-      l=1
-      ir=n
-      fx=0.
- 10   if(ir-l .lt.m)then
-          do j=l+1,ir
-            a=arr1(j)
-            b=arr2(j)
-            do i=j-1,1,-1
-              if (arr1(i).le.a) goto 12
-              arr1(i+1)=arr1(i)
-              arr2(i+1)=arr2(i)
-            enddo
-            i=0
- 12         arr1(i+1)=a
-            arr2(i+1)=b
-          enddo
-          if(jstack.eq.0)return
-          ir=istack(jstack)
-          l=istack(jstack-1)
-          jstack=jstack-2
-      else
-          i=l
-          j=ir
-          fx=MOD(fx*fa+fc,fm)
-          iq=l+(ir-l+1)*(fx*fmi)
-          a=arr1(iq)
-          arr1(iq)=arr1(l)
-          b=arr2(iq)
-          arr2(iq)=arr2(l)
- 20       continue
- 21       if (j.gt.0) then
-              if(a.lt.arr1(j))then
-                  j=j-1
-                  goto 21
-              endif
-          endif
-          if (j.le.i)then
-              arr1(i)=a
-              arr2(i)=b
-              goto 30
-          endif
-          arr1(i)=arr1(j)
-          arr2(i)=arr2(j)
-          i=i+1
- 22       if (i.le.n) then
-              if (a.gt.arr1(i)) then
-                  i=i+1
-                  goto 22
-              endif
-          endif
-          if (j.le.i)then
-              arr1(j)=a
-              arr2(j)=b
-              i=j
-              goto 30
-          endif
-          arr1(j)=arr1(i)
-          arr2(j)=arr2(i)
-          j=j-1
-          goto 20
- 30       jstack=jstack+2
-          if(jstack.gt.nstack)pause 'nstack trop petit'
-          if (ir-i.ge.i-1) then
-              istack(jstack)=ir
-              istack(jstack-1)=i+1
-              ir=i-1
-          else
-              istack(jstack)=i-1
-              istack(jstack-1)=l
-              l=i+1
-          endif
-      endif
-      goto 10
+     ist=0
+     il=1
+     ir=n
+     zfx=0.
+     !
+ 10  IF( ir-il .LT. m )THEN
+
+        DO jj = il+1,ir
+           za = arr1(jj)
+           zb = arr2(jj)
+           DO i = jj-1,1,-1
+              IF( arr1(i) .LE. za )GOTO 12
+              arr1(i+1) = arr1(i)
+              arr2(i+1) = arr2(i)
+           ENDDO
+           i=0
+ 12        arr1(i+1) = za
+           arr2(i+1) = zb
+        ENDDO
+        IF( ist .EQ. 0 )RETURN
+        ir  = istack(ist)
+        il  = istack(ist-1)
+        ist = ist-2
+     ELSE
+        i  = il
+        jj = ir
+        zfx = MOD(zfx*fa+fc,fm)
+        iq = il+(ir-il+1)*(zfx*fmi)
+        za  = arr1(iq)
+        arr1(iq) = arr1(il)
+        zb = arr2(iq)
+        arr2(iq) = arr2(il)
+ 20     CONTINUE
+ 21     IF( jj .GT. 0 )THEN
+           IF( za .LT. arr1(jj) )THEN
+              jj = jj-1
+              GOTO 21
+           ENDIF
+        ENDIF
+        IF( jj .LE. i )THEN
+           arr1(i) = za
+           arr2(i) = zb
+           GOTO 30
+        ENDIF
+        arr1(i) = arr1(jj)
+        arr2(i) = arr2(jj)
+        i=i+1
+ 22     IF( i .LE. n )THEN
+           IF( za .GT. arr1(i) )THEN
+              i = i+1
+              GOTO 22
+           ENDIF
+        ENDIF 
+        IF( jj .LE.i )THEN
+           arr1(jj) = za
+           arr2(jj) = zb
+           i        = jj
+           GOTO 30
+        ENDIF
+        arr1(jj) = arr1(i)
+        arr2(jj) = arr2(i)
+        jj = jj - 1
+        GOTO 20
+ 30     ist = ist + 2
+        IF( ist .GT. nstack )PAUSE 'nstack too small'
+        IF( ir-i .GE. i-1 )THEN
+           istack(ist  ) = ir
+           istack(ist-1) = i + 1
+           ir = i-1
+        ELSE
+           istack(ist  ) = i-1
+           istack(ist-1) = il
+           il = i+1
+        ENDIF
+     ENDIF
+     GOTO 10
   END SUBROUTINE qcksrt
 
   SUBROUTINE write_debug(ksec,cd_write)
      !!---------------------------------------------------------------------
-     !!         ***  SUBROUTINE qcksrt  ***
-     !!
-     !! ** Purpose
-     !!
-     !! ** Method
-     !!
-     !! ** Input
-     !!
-     !! ** Action
-     !!
-     !! History
+     !!         ***  SUBROUTINE write_debug  ***
+     !!
+     !! ** Purpose: write debug messages
+     !!
      !!---------------------------------------------------------------------
      !! * arguments
@@ -545,13 +561,6 @@
      !!         ***  ROUTINE file_open ***
      !!
-     !! ** Purpose
-     !!
-     !! ** Method
-     !!
-     !! ** Input
-     !!
-     !! ** Action
-     !!
-     !! Hisotry
+     !! ** Purpose: open a file
+     !!
      !!---------------------------------------------------------------------
      !! * arguments
