Trailing-Edge
-
PDP-10 Archives
-
decuslib10-06
-
43,50416/linpak.sai
There are 2 other files named linpak.sai in the archive. Click here to see a list.
Entry;
COMMENT
.SEC(LINPAK.SAI - line geometry package)
.index(LINPAK.SAI - line geometry package)
.;
Begin "LINPAK.SAI"
COMMENT
Bruce Shapiro
Image Processing Unit
Division of Cancer Biology and Diagnosis
National Cancer Institute
National Institutes of Health
Bethesda, Maryland 20014
301-496-2394
Revised August 6, 1976 - Lemkin, Shapiro removed debugs
Revised August 2, 1976 - Shapiro fixed infinite slope in
routine LINE
Reved May 17, 1976 - Lemkin, made procedures SIMPLE
.<< January 9, 1976
. Feb 26, 1976 - Lemkin, Shapiro added PUB statements
. March 3, 1976 - Lemkin, edited long lines for PUB
. March 10, 1976 - Lemkin, PUB edit>>
Abstract
--------
This set of routines will manipulate lines and vectors.
It currently deals with 2-D space with the following
capabilities:
. datetitle_(DATE&" "&TIME)
.NEXT PAGE
.SS(List of procedure calls)
.INDEX(List of procedure calls)
The following is a list of procedures contained in
LINPAK. The actual argument data types is contained in a list
in the next section.
1. LINE(X1,Y1,X2,Y2,M,B)- requires 2 points and
returns the slope and intercept of the
line passing through these points.
2. MIDPOINT(X1,Y1,X2,Y2,Xm,Ym)-requires 2 points and
returns the mid-point of the line joining
these two points.
3. PERPENDICULAR(M,B,X,Y,Mp,Bp)- requires the slope
intercept and a point on the line. It returns
the slope and intercept of the line PERPENDICULAR
to the given line a the specified point.
4. INTERSECTION(M1,M2,B1,B2,Xi,Yi)- requires the slopes
and intercepts of two lines. It returns the
point of INTERSECTION of the two lines.
5. DCOSX(X1,Y1,X2,Y2)- given a directed line segment
this function returns the direction cosine
of the line measured from the positive X axis.
6. DCOSY(X1,Y1,X2,Y2)- given a directed line segment, this
function returns the direction cosine of this
segment measured from the positive Y axis.
7. EDISTANCE(X1,Y1,X2,Y2)- given a pair of points this function
returns the euclidian distance between these
points.
8. SIGN(Num)- this function returns either -1. or +1. depending
on the SIGN of the argument (Num).
9. VECTORANGLE(X1h,Y1h,X1t,Y1t,X2h,Y2h,X2t,Y2t)- this function
returns the angle (in radians) between two directed
line segments. The head and tail of the segments
must be specified in the correct order. The first
segment indicated in the arguments (distinguish by 1)
is being considered as the base line. An angle
between 0 and 2 PI radians is returned measured in
a counter-clockwise direction from the directed
directed base line.
10. DEFLECTVECTOR(Xh,Yh,Xc,Yc,Dangle,NewXh,NewYh)- this
procedure when presented with a base vector
head and tail coordinates (Xh,Yh) and (Xc,Yc)
respectively and a deflection angle will return
the head coordinates of a new vector at Dangle
radians from the base vector. The new vector's tail
will remain at (Xc,Yc) and its head will be
unit distance from (Xc,Yc).
11. SCALEVECTOR(Xh,Yh,Xt,Yt,Scale,NewXh,NewYh)-this
procedure when presented with a vector
defined by its head and tail coordinates,
(Xh,Yh) and (Xt,Yt) respectively and a
scale factor (Scale) will return the head
coordinates of a vector of magnitude Scale
and direction the same as the given vector.
Note that if a negative scale value is given
the direction of the new vector will be opposite
the given vector. The tail coordinates will
not be changed.
12. FIND!REC(Boundary,Size,Xmin,Xmax,Ymin,Ymax)- this
procedure when presented with a boundary file
and its size will compute the extrema rectangle
for this file. The results are returned in
Xmin, Ymin, Xmax, Ymax.
.NEXT PAGE
.SS(Table of Procedure Arguments)
.INDEX(Table of Procedure Arguments)
x1 - Real x1
x2 - Real x2
y1 - Real y1
y2 - Real y2
m1 - Real m1
m2 - Real m2
b1 - Real b1
b2 - Real b2
num - Real num
scale - Real scale
x1h - Real x1h
x2h - Real x2h
y1h - Real y1h
y2h - Real y2h
x1t - Real x1t
x2t - Real x2t
y1t - Real y1t
y2t - Real y2t
xh - Real xh
yh - Real yh
xc - Real xc
yc - Real yc
xt - Real xt
yt - Real yt
dangle - Real dangle
b - Reference Real b
m - Reference Real m
xm - Reference Real xm
ym - Reference Real ym
bp - Reference Real bp
mp - Reference Real mp
xi - Reference Real xi
yi - Reference Real yi
newxh - Reference Real newxh
newyh - Reference Real newyh
xmin - Reference Real xmin
xmax - Reference Real xmax
ymin - Reference Real ymin
ymay - Reference Real ymay
;
COMMENT
.NEXT PAGE
.ss(REQUIRE files)
.;
Require "DEFINE.REQ" Source!file;
Require "SOLVER.REQ" Source!file;
COMMENT
.NEXT PAGE
.SS(LINE)
.INDEX(LINE)
.;
Internal Simple Procedure LINE(Real X1,Y1,X2,Y2;Reference Real M,B);
Begin "Line"
Real
Denom;
#--------------------------------------------------------------
# This procedure when given two points (X1,Y1) and (X2,Y2)
# will return the slope and intercept of a line passing through
# the two points.
#---------------------------------------------------------------;
If (X2-X1)=0 then
Begin "Infinite Slope"
Outstr("Infinite slope in Line"&crlf);
Denom_.00000001;
End "Infinite Slope"
Else
Denom_(X2-X1);
M_(Y2-Y1)/Denom;
B_Y1-M*X1;
Return;
End "Line";
COMMENT
.NEXT PAGE
.SS(MIDPOINT)
.INDEX(MIDPOINT)
.;
Internal Simple Procedure MIDPOINT(Real X1,Y1,X2,Y2;
Reference Real Xm,Ym);
Begin "Midpoint"
#-----------------------------------------------------------
# This procedure when passed the end points of a line
# segment will return its mid-point.
#-----------------------------------------------------------;
Xm_(X2+X1)/2.;
Ym_(Y2+Y1)/2.;
Return;
End "Midpoint";
COMMENT
.NEXT PAGE
.SS(PERPENDICULAR)
.INDEX(PERPENDICULAR)
.;
Internal Simple Procedure PERPENDICULAR(Real M,B,X,Y;
Reference Real Mp,Bp);
Begin "PERPENDICULAR"
#------------------------------------------------------------
# This procedure when passed the slope (M) and intercept (B)
# if a line and a point (X,Y) on the line will return
# the slope (Mp)
# and intercept (Bp) of a line PERPENDICULAR to the given line
# at the given point.
#--------------------------------------------------------------------;
# If M=0 then Outstr("0 slope in PERPENDICULAR"&crlf);
if Y neq M*X+B
then Outstr("Point not on given line in PERPENDICULAR"&crlf);
# Mp_-1./M;
Mp_If M=0 then 9999999. else -1./m;
Bp_Y-Mp*X;
Return;
End "PERPENDICULAR";
COMMENT
.NEXT PAGE
.SS(INTERSECTION)
.NEXT PAGE
.SS(INTERSECTION)
.;
Internal Procedure INTERSECTION(Real M1,M2,B1,B2;
Reference Real Xi,Yi);
Begin "INTERSECTION"
#------------------------------------------------------------------
# This procedure when passed the slope of two lines (M1 and M2)
# and the Y intercepts (B1 and B2) will return the point of
# INTERSECTION of the two line (Xi,Yi).
#-----------------------------------------------------------------;
Real Array
A[1:2,1:2],
B[1:2];
Integer Array
Ip[1:2],
Jp[1:2];
A[2,1]_A[1,1]_1.;
A[1,2]_-M1;
A[2,2]_-M2;
DECOMP(A,2,Ip,Jp);
If Ip[2]=0 then Outstr("Singular Matrix in INTERSECTION"&crlf);
B[1]_B1;
B[2]_B2;
SOLVE(A,2,Ip,B);
Yi_B[1];
Xi_B[2];
Return;
End "INTERSECTION";
COMMENT
.NEXT PAGE
.SS(DCOSX)
.INDEX(DCOSX)
.;
Internal Real Simple Procedure DCOSX(Real X1,Y1,X2,Y2);
Begin "DCOSX"
#--------------------------------------------------------------
# This function when passed a directed line segment will
# return the direction cosine of the segment relative to
# positive X axis. The head of the segment is the point
# (X2,Y2) and the tail is the point (X1,Y1).
#
#-----------------------------------------------------------;
Return((X2-X1)/Sqrt((X2-X1)^2+(Y2-Y1)^2));
End "DCOSX";
COMMENT
.NEXT PAGE
.SS(DCOSY)
.INDEX(DCOSY)
.;
Internal Real Simple Procedure DCOSY(Real X1,Y1,X2,Y2);
Begin "DCOSY"
#--------------------------------------------------------------
# This function when passed a directed line segment will
# return the direction cosine of the segment relative to
# positive Y axis. The head of the segment is the point
# (X2,Y2) and the tail is the point (X1,Y1).
#
#-----------------------------------------------------------;
Return((Y2-Y1)/Sqrt((X2-X1)^2+(Y2-Y1)^2));
End "DCOSY";
COMMENT
.NEXT PAGE
.SS(EDISTANCE)
.INDEX(EDISTANCE)
.;
Internal Real Simple Procedure EDISTANCE(Real X1,Y1,X2,Y2);
#---------------------------------------------------------
# This function returns the euclidian distance between
# a pair of points in the cartisian coordinate system.
#----------------------------------------------------------;
Return(Sqrt((X1-X2)^2+(Y1-Y2)^2));
COMMENT
.NEXT PAGE
.SS(SIGN)
.INDEX(SIGN)
.;
Internal Real Simple Procedure SIGN(Real Num);
Begin "SIGN"
#-------------------------------------------------------------
# This function returns the SIGN +1. or -1. of its
# argument.
#--------------------------------------------------------------;
Real
Value;
Value_Num/Abs(Num);
Value_If Value=0 then 1. else Value;
Return(Value);
End "SIGN";
COMMENT
.NEXT PAGE
.SS(VECTORANGLE)
.INDEX(VECTORANGLE)
.;
Internal Real Simple Procedure VECTORANGLE(Real X1h, Y1h, X1t, Y1t,
X2h, Y2h, X2t, Y2t);
Begin "VECTORANGLE"
#----------------------------------------------------------------
# This function returns the angle between two directed
# line segments. The head and tail of the segments are indicated
# in the argument list by the letters h and t repectively. The
# first segment indicated is taken as the base segment. The
# returned angle is computed by moving in a counter-clockwise
# direction from this base segment. The angle will lie
# within the closed interval 0 to 2 PI.
#-----------------------------------------------------------------;
Real
Fangle,
Bangle,
SIGNf,
SIGNb,
Anglediff,
Cosfx,
Cosfy,
Cosbx,
Cosby;
#----------------------------------------------------------------
# Determine the direction cosines for the base segment.
#----------------------------------------------------------------;
Cosbx_DCOSX(X1t,Y1t,X1h,Y1h);
Cosby_DCOSY(X1t,Y1t,X1h,Y1h);
#----------------------------------------------------------------
# Determine the direction cosines for the non-base segment.
#---------------------------------------------------------------;
Cosfx_DCOSX(X2t,Y2t,X2h,Y2h);
Cosfy_DCOSY(X2t,Y2t,X2h,Y2h);
#-------------------------------------------------------------------
# Determine the absolute coordinate system quadrant
# in which the segments lie based upon the direction cosine
# measured relative to the absolute Y axis. and the
# arc cosine of the direction cosine measured relative to
# the absolute X axis.
#---------------------------------------------------------------;
SIGNf_SIGN(Cosfy);
SIGNb_SIGN(Cosby);
Fangle_Acos(Cosfx)*SIGNf;
Bangle_Acos(Cosbx)*SIGNb;
#--------------------------------------------------------------
# Now determine the angular difference.
# Make sure the resultant angles are
# between 0 and + 360 degrees.
#---------------------------------------------------------------;
Anglediff_Fangle-Bangle;
Anglediff_If Anglediff<0. then 6.2832-Abs(Anglediff) else Anglediff;
Return(Anglediff);
End "VECTORANGLE";
COMMENT
.NEXT PAGE
.SS(DEFLECTVECTOR)
.INDEX(DEFLECTVECTOR)
.;
Internal Simple Procedure DEFLECTVECTOR(Real Xh,Yh,Xc,Yc,Dangle;
Reference Real NewXh,NewYh);
Begin "DEFLECTVECTOR"
#--------------------------------------------------------------
# This procedure when presented with a base vector with
# head and tail coordinates (Xh,Yh) and (Xc,Yc) repectively
# and a deflection angle (Dangle) will return the
# head coordinates of a new vector. The new vector's tail
# will remain at (Xc,Yc). The angle, Dangle may be positive
# or negative. It is measured in a counter-clockwise direction
# from the base vector. The head of the new vector will
# be unit distance from the coordinates (Xc,Yc).
#------------------------------------------------------------------;
Real
Dx,
Dy,
Angle,
Newangle,
DeltaX,
DeltaY;
#-----------------------------------------------------------------
# Compute the direction cosines of the base vector.
#-----------------------------------------------------------------;
Dx_DCOSX(Xc,Yc,Xh,Yh);
Dy_DCOSY(Xc,Yc,Xh,Yh);
#------------------------------------------------------------------
# Determine the absolute angular position of the base vector
# and make sure it lies within 0 to 2 PI radians.
#------------------------------------------------------------------;
Angle_Acos(Dx)*SIGN(Dy);
Angle_If Angle <0 then twopi+Angle else Angle;
#----------------------------------------------------------
# Get the absolute angular position of the deflected
# vector.
#---------------------------------------------------------;
Newangle_Dangle+Angle;
#----------------------------------------------------------
# Determine the X and Y displacements for the new vector
# and the position of the head of the new vector.
# Note that head position is unit distance from
# (Xc,Yc).
#------------------------------------------------------------;
DeltaX_Cos(Newangle);
DeltaY_Sin(Newangle);
NewXh_Xc+DeltaX;
NewYh_Yc+DeltaY;
Return;
End "DEFLECTVECTOR";
COMMENT
.NEXT PAGE
.SS(SCALEVECTOR)
.INDEX(SCALEVECTOR)
.;
Internal Simple Procedure SCALEVECTOR(Real Xh,Yh,Xt,Yt,Scale;
Reference Real NewXh,NewYh);
Begin "SCALEVECTOR"
#------------------------------------------------------------------
# This procedure when presented with the head and tail
# coordinates of a vector (Xh,Yh) and (Xt,Yt) respectively,
# and a scale factor (Scale) will return the head coordinates
# of a new vector with magnitude scale and direction the
# same as the given vector. Note that a negative scale
# value will point the vector in the opposite direction.
# The tail coordinates remain the same.
#---------------------------------------------------------------;
Real
Dx,
Dy,
Deltax,
Deltay;
#-------------------------------------------------------------------
Find direction cosines of given vector.
#-------------------------------------------------------------------;
Dx_DCOSX(Xt,Yt,Xh,Yh);
Dy_DCOSY(Xt,Yt,Xh,Yh);
#------------------------------------------------------------
# Determine delta values to determine the new head
# vector coordinates.
#-------------------------------------------------------------;
Deltax_Scale*Dx;
Deltay_Scale*Dy;
NewXh_Xt+Deltax;
NewYh_Yt+Deltay;
End "SCALEVECTOR";
COMMENT
.NEXT PAGE
.SS(FIND!REC)
.INDEX(FIND!REC)
.;
Internal Simple Procedure FIND!REC(Integer Array Boundary;
Integer Size;
Reference Integer Xmin,Xmax,Ymin,Ymax);
Begin "FIND!REC"
#------------------------------------------------------
# This procedure when presented with a boundary file
# will compute the rectangular extrema of the file.
# The values are returned in Xmin,Xmax,Ymin,Ymax.
#-------------------------------------------------------;
Integer
X,
Y,
I;
Xmin_maxinteger;
Ymin_maxinteger;
Xmax_mininteger;
Ymax_mininteger;
For I_0 step 1 until size-1 do
Begin
X_X!BND!FETCH(Boundary,I);
Y_Y!BND!FETCH(Boundary,I);
Xmax_X max Xmax;
Ymax_Y max Ymax;
Xmin_X min Xmin;
Ymin_ Y min Ymin;
End;
End "FIND!REC";
End "LINPAK.SAI";