!**********************************************

  BASIC GEOMETRY MODULE sectarea(
   REF   l_id  >"ID of left hull surface";
   REF   wp_id >"ID of waterline b_plane";
   FLOAT ywl   >"Y-coordinate of waterline plane";
   FLOAT x     >"X-coordinate of requested section";
   VAR FLOAT beam;
   VAR FLOAT area);

!*  Calculates the area of the submerged part
!*  of a hull section.
!*
!*  This program is free software; you can redistribute it and/or
!*  modify it under the terms of the GNU General Public License
!*  as published by the Free Software Foundation; either version 2
!*  of the License, or (at your option) any later version.
!*
!*  This program is distributed in the hope that it will be useful,
!*  but WITHOUT ANY WARRANTY; without even the implied warranty of
!*  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!*  GNU General Public License for more details.
!*
!*  You should have received a copy of the GNU General Public License
!*  along with this program; if not, write to the Free Software
!*  Foundation, Inc., 59 Temple Place - Suite 330, Boston,
!*  MA  02111-1307, USA.
!*
!*  (C)Microform AB 1998-09-23, J.Kjellander
!*
!**********************************************

 VECTOR pa,pb,tmp,pwl;
 FLOAT  t,t1,t2,dl;
 INT    font,ns_geo,ns_gra;

BEGINMODULE

!*
!***Intersect between hull and XY-plane of #1.
!*
   csys_1p(#1,"Cutplane",vec(x,0,0),0,90,0:BLANK=1,PEN=5);
   cur_int(#2,l_id,#1,1:BLANK=1,PEN=5);
!*
!***Intersect between curve and waterline plane.
!*
   pwl:=intersect(#2,wp_id);
!*
!***Cut #2 and trim it to pwl. Compare Y-coordinates of
!***both ends of #2 to select the right part to keep.
!*
   t:=inv_on(#2,pwl,1);
   pa:=startp(#2);
   pb:=endp(#2);

   if pa.y < pb.y then
     t1:=0;
     t2:=t;
   else
     t1:=t;
     getcurh(#2,font,dl,ns_geo,ns_gra);
!     t2:=inv_on(#2,pb);  Detta ger felstopp !!!!!
      t2:=ns_geo;
   endif;

   cur_trim(#3,#2,t1,t2:BLANK=0,PEN=2);
!*
!***Connect ends of intersect curve with lines
!***and calculate section area..
!*
   pa:=startp(#3);
   pb:=endp(#3);
   if pa.y > pb.y then
     tmp:=pa;
     pa:=pb;
     pb:=tmp;
   endif;

   beam:=2.0*pb.z;
   
   if abs(pa.z) > 0.1 then
     lin_free(#7,pa,vec(pa.x,pa.y,0.0):BLANK=0,PEN=9);
     pa:=endp(#7);
     lin_free(#4,pa,vec(pa.x,ywl,0):BLANK=0,PEN=3);
     lin_free(#5,endp(#4),pb:BLANK=0,PEN=4);
     cur_comp(#6,#7,#4,#5,#3:BLANK=1,PEN=5);
     area:=2.0*section_area(#6,1.0);
   else
     lin_free(#8,pa,vec(pa.x,ywl,0):BLANK=0,PEN=3);
     lin_free(#9,endp(#8),pb:BLANK=0,PEN=4);
     cur_comp(#10,#8,#9,#3:BLANK=1,PEN=5);
     area:=2.0*section_area(#10,1.0);
   endif;

ENDMODULE

!*****************************************************