Xiaolin Wu's line algorithm: Difference between revisions
→{{header|FreeBASIC}}
Line 1,843:
im2.savePPM6("xiaolin_lines2.ppm");
}</syntaxhighlight>
=={{header|Fortran}}==
{{trans|Common Lisp}}
The program outputs a Portable Gray Map representing a transparency mask. The mask is full of straight lines, drawn by a solution of this Rosetta Code task. Some of the lines draw a piecewise approximation of an ellipse, and others draw the normals of the ellipse. The normals form an envelope that is the evolute of the ellipse.
<syntaxhighlight lang="fortran">
program xiaolin_wu_line_algorithm
use, intrinsic :: ieee_arithmetic
implicit none
type :: drawing_surface
integer :: u0, v0, u1, v1
real, allocatable :: pixels(:)
end type drawing_surface
interface
subroutine point_plotter (s, x, y, opacity)
import drawing_surface
type(drawing_surface), intent(inout) :: s
integer, intent(in) :: x, y
real, intent(in) :: opacity
end subroutine point_plotter
end interface
real, parameter :: pi = 4.0 * atan (1.0)
integer, parameter :: u0 = -499
integer, parameter :: v0 = -374
integer, parameter :: u1 = 500
integer, parameter :: v1 = 375
real, parameter :: a = 200.0 ! Ellipse radius on x axis.
real, parameter :: b = 350.0 ! Ellipse radius on y axis.
type(drawing_surface) :: s
integer :: i, step_size
real :: t, c, d
real :: x, y
real :: xnext, ynext
real :: u, v
real :: rhs, ad, bc
real :: x0, y0, x1, y1
s = make_drawing_surface (u0, v0, u1, v1)
! Draw both an ellipse and the normals of that ellipse.
step_size = 2
do i = 0, 359, step_size
! Parametric representation of an ellipse.
t = i * (pi / 180.0)
c = cos (t)
d = sin (t)
x = a * c
y = b * d
! Draw a piecewise linear approximation of the ellipse. The
! endpoints of the line segments will lie on the curve.
xnext = a * cos ((i + step_size) * (pi / 180.0))
ynext = b * sin ((i + step_size) * (pi / 180.0))
call draw_line (s, x, y, xnext, ynext)
! The parametric equation of the normal:
!
! (a * sin (t) * xnormal) - (b * cos (t) * ynormal)
! = (a**2 - b**2) * cos (t) * sin (t)
!
! That is:
!
! (a * d * xnormal) - (b * c * ynormal) = (a**2 - b**2) * c * d
!
rhs = (a**2 - b**2) * c * d
ad = a * d
bc = b * c
if (abs (ad) < abs (bc)) then
x0 = -1000.0
y0 = ((ad * x0) - rhs) / bc
x1 = 1000.0
y1 = ((ad * x1) - rhs) / bc
else
y0 = -1000.0
x0 = (rhs - (bc * y0)) / ad
y1 = 1000.0
x1 = (rhs - (bc * y1)) / ad
end if
call draw_line (s, x0, y0, x1, y1)
end do
call write_transparency_mask (s)
contains
function make_drawing_surface (u0, v0, u1, v1) result (s)
integer, intent(in) :: u0, v0, u1, v1
type(drawing_surface) :: s
integer :: w, h
if (u1 < u0 .or. v1 < v0) error stop
s%u0 = u0; s%v0 = v0
s%u1 = u1; s%v1 = v1
w = u1 - u0 + 1
h = v1 - v0 + 1
allocate (s%pixels(0:(w * h) - 1), source = 0.0)
end function make_drawing_surface
function drawing_surface_ref (s, x, y) result (c)
type(drawing_surface), intent(in) :: s
integer, intent(in) :: x, y
real :: c
! In calls to drawing_surface_ref and drawing_surface_set, indices
! outside the drawing_surface are allowed. Such indices are
! treated as if you were trying to draw on the air.
if (s%u0 <= x .and. x <= s%u1 .and. s%v0 <= y .and. y <= s%v1) then
c = s%pixels((x - s%u0) + ((s%v1 - y) * (s%u1 - s%u0 + 1)))
else
c = ieee_value (s%pixels(0), ieee_quiet_nan)
end if
end function drawing_surface_ref
subroutine drawing_surface_set (s, x, y, c)
type(drawing_surface), intent(inout) :: s
integer, intent(in) :: x, y
real, intent(in) :: c
! In calls to drawing_surface_ref and drawing_surface_set, indices
! outside the drawing_surface are allowed. Such indices are
! treated as if you were trying to draw on the air.
if (s%u0 <= x .and. x <= s%u1 .and. s%v0 <= y .and. y <= s%v1) then
s%pixels((x - s%u0) + ((s%v1 - y) * (s%u1 - s%u0 + 1))) = c
end if
end subroutine drawing_surface_set
subroutine write_transparency_mask (s)
type(drawing_surface), intent(in) :: s
! Write a transparency mask, in plain (ASCII or EBCDIC) Portable
! Gray Map format, but representing opacities rather than
! whitenesses. (Thus there will be no need for gamma corrections.)
! See the pgm(5) manpage for a discussion of this use of PGM
! format.
integer :: w, h
integer :: i
w = s%u1 - s%u0 + 1
h = s%v1 - s%v0 + 1
write (*, '("P2")')
write (*, '("# transparency mask")')
write (*, '(I0, 1X, I0)') w, h
write (*, '("255")')
write (*, '(15I4)') (nint (255 * s%pixels(i)), i = 0, (w * h) - 1)
end subroutine write_transparency_mask
subroutine draw_line (s, x0, y0, x1, y1)
type(drawing_surface), intent(inout) :: s
real, intent(in) :: x0, y0, x1, y1
real :: xdiff, ydiff
xdiff = abs (x1 - x0)
ydiff = abs (y1 - y0)
if (ydiff <= xdiff) then
if (x0 <= x1) then
call drawln (s, x0, y0, x1, y1, plot_shallow)
else
call drawln (s, x1, y1, x0, y0, plot_shallow)
end if
else
if (y0 <= y1) then
call drawln (s, y0, x0, y1, x1, plot_steep)
else
call drawln (s, y1, x1, y0, x0, plot_steep)
end if
end if
end subroutine draw_line
subroutine drawln (s, x0, y0, x1, y1, plot)
type(drawing_surface), intent(inout) :: s
real, intent(in) :: x0, y0, x1, y1
procedure(point_plotter) :: plot
real :: dx, dy, gradient
real :: yend, xgap
real :: first_y_intersection, intery
integer :: xend
integer :: xpxl1, ypxl1
integer :: xpxl2, ypxl2
integer :: x
dx = x1 - x0; dy = y1 - y0
if (dx == 0.0) then
gradient = 1.0
else
gradient = dy / dx
end if
! Handle the first endpoint.
xend = iround (x0)
yend = y0 + (gradient * (xend - x0))
xgap = rfpart (x0 + 0.5)
xpxl1 = xend
ypxl1 = ipart (yend)
call plot (s, xpxl1, ypxl1, rfpart (yend) * xgap)
call plot (s, xpxl1, ypxl1 + 1, fpart (yend) * xgap)
first_y_intersection = yend + gradient
! Handle the second endpoint.
xend = iround (x1)
yend = y1 + (gradient * (xend - x1))
xgap = fpart (x1 + 0.5)
xpxl2 = xend
ypxl2 = ipart (yend)
call plot (s, xpxl2, ypxl2, (rfpart (yend) * xgap))
call plot (s, xpxl2, ypxl2 + 1, fpart (yend) * xgap)
! Loop over the rest of the points.
intery = first_y_intersection
do x = xpxl1 + 1, xpxl2 - 1
call plot (s, x, ipart (intery), rfpart (intery))
call plot (s, x, ipart (intery) + 1, fpart (intery))
intery = intery + gradient
end do
end subroutine drawln
subroutine plot_shallow (s, x, y, opacity)
type(drawing_surface), intent(inout) :: s
integer, intent(in) :: x, y
real, intent(in) :: opacity
real :: combined_opacity
! Let us simply add opacities, up to the maximum of 1.0. You might
! wish to do something different, of course.
combined_opacity = opacity + drawing_surface_ref (s, x, y)
call drawing_surface_set (s, x, y, min (combined_opacity, 1.0))
end subroutine plot_shallow
subroutine plot_steep (s, x, y, opacity)
type(drawing_surface), intent(inout) :: s
integer, intent(in) :: x, y
real, intent(in) :: opacity
call plot_shallow (s, y, x, opacity)
end subroutine plot_steep
elemental function ipart (x) result (i)
real, intent(in) :: x
integer :: i
i = floor (x)
end function ipart
elemental function iround (x) result (i)
real, intent(in) :: x
integer :: i
i = ipart (x + 0.5)
end function iround
elemental function fpart (x) result (y)
real, intent(in) :: x
real :: y
y = modulo (x, 1.0)
end function fpart
elemental function rfpart (x) result (y)
real, intent(in) :: x
real :: y
y = 1.0 - fpart (x)
end function rfpart
end program xiaolin_wu_line_algorithm
</syntaxhighlight>
Here is a shell script that runs the program and creates a PNG. The background of the image is one pattern, and the foreground is another. The foreground, however, is masked by the transparency mask, and so only all those straight lines we drew show up in the PNG.
<syntaxhighlight lang="shell">
#!/bin/sh
# Using the optimizer, even at low settings, avoids trampolines and
# executable stacks.
gfortran -std=f2018 -g -O1 xiaolin_wu_line_algorithm.f90
./a.out > alpha.pgm
ppmpat -anticamo -randomseed=36 1000 750 | pambrighten -value=-60 -saturation=50 > fg.pam
ppmpat -poles -randomseed=57 1000 750 | pambrighten -value=+200 -saturation=-80 > bg.pam
pamcomp -alpha=alpha.pgm fg.pam bg.pam | pamtopng > image.png
</syntaxhighlight>
{{out}}
[[File:Xiaolin wu line algorithm Fortran.png|thumb|none|alt=An ellipse and its normals (forming an evolute as their envelope), all of it in goofy colors.]]
=={{header|FreeBASIC}}==
| |||