# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.

include	<imhdr.h>

# IMA_DIV -- Image arithmetic division.

$for (silrd)
procedure ima_div$t (im_a, im_b, im_c, a, b, c)

pointer	im_a, im_b, im_c
PIXEL	a, b, c

int	len
pointer	im[3], buf[3]
long	v[IM_MAXDIM, 3]

int	ima_nl$t()
PIXEL	ima_efnc$t()
extern	ima_efnc$t

PIXEL	divzero
common	/imadcom$t/ divzero

begin
	# Loop through all of the image lines.
	divzero = c
	im[1] = im_c
	len = IM_LEN (im[1], 1)
	call amovkl (long(1), v, 3 * IM_MAXDIM)

	# If imagea is constant then read imageb and do a vector
	# reciprical to imagec.
	if (im_a == NULL) {
	    im[2] = im_b
	    while (ima_nl$t (im, buf, v, 2) != EOF)
		call arcz$t (a, Mem$t[buf[2]], Mem$t[buf[1]], len,
		    ima_efnc$t)

	# If imageb is constant then read imagea.  If the constant
	# is 1 do a vector move to imagec otherwise do a vector/scalar
	# divide to imagec.
	} else if (im_b == NULL) {
	    im[2] = im_a
	    while (ima_nl$t (im, buf, v, 2) != EOF) {
		if (b == 0$f)
		    call amovk$t (divzero, Mem$t[buf[1]], len)
		else if (b == 1$f)
		    call amov$t (Mem$t[buf[2]], Mem$t[buf[1]], len)
		else
		    call adivk$t (Mem$t[buf[2]], b, Mem$t[buf[1]], len)
	    }

	# Read imagea and imageb and do the vector divide to imagec.
	} else {
	    im[2] = im_a
	    im[3] = im_b
	    while (ima_nl$t (im, buf, v, 3) != EOF)
		call advz$t (Mem$t[buf[2]], Mem$t[buf[3]], Mem$t[buf[1]],
		    len, ima_efnc$t)
	}
end


# IMA_EFNC -- Error function for division by zero.

PIXEL procedure ima_efnc$t (a)

PIXEL	a
PIXEL	divzero
common	/imadcom$t/ divzero

begin
	return (divzero)
end
$endfor