1 files changed, 86 insertions, 0 deletions
diff --git a/pkg/proto/vol/src/imjoin.gx b/pkg/proto/vol/src/imjoin.gx
new file mode 100644
index 00000000..04d2cc93
--- /dev/null
+++ b/pkg/proto/vol/src/imjoin.gx
@@ -0,0 +1,86 @@
+include <imhdr.h>
+
+define	VPTR		Memi[$1+$2-1]	# Array of axis vector pointers
+
+
+# IMJOIN -- Join the set of input images into an output image along the
+# specified axis, any dimension up to 7 (NOT necessarily IM_MAXDIM!).
+
+$for (silrdx)
+procedure imjoin$t (inptr, nimages, out, joindim, outtype)
+pointer	inptr[nimages]		# Input IMIO pointers
+int	nimages			# Number of input images
+pointer	out			# Output IMIO pointer
+int	joindim			# Dimension along which to join images
+int	outtype			# Output datatype
+
+pointer	in, inbuf, outbuf, sp, vin
+int	stat, image, cum_len, line, band, dim4, dim5, dim6, dim7
+long	vout[IM_MAXDIM]
+
+pointer	imgnl$t()
+pointer	impnl$t()
+
+begin
+	call smark (sp)
+	call salloc (vin, nimages, TY_INT)
+
+	call amovkl (long(1), vout, IM_MAXDIM)
+	do image = 1, nimages {
+	    call salloc (VPTR(vin,image), IM_MAXDIM, TY_LONG)
+	    call amovkl (long(1), Meml[VPTR(vin,image)], IM_MAXDIM)
+	}
+
+	# Join input images along specified dimension.  Joins along columns
+	# and lines require processing in special order, all others in the
+	# same order.  In the first two cases we process all input images
+	# in inner loops, so we have to keep all those imdescriptors open.
+
+	switch (joindim) {
+	case 1:						# join columns
+	    do band = 1, IM_LEN(out,3)
+		do line = 1, IM_LEN(out,2) {
+		    stat = impnl$t (out, outbuf, vout)
+		    cum_len = 0
+		    do image = 1, nimages {
+			in = inptr[image]
+			stat = imgnl$t (in, inbuf, Meml[VPTR(vin,image)])
+			call amov$t (Mem$t[inbuf], Mem$t[outbuf+cum_len],
+			    IM_LEN(in,1))
+			cum_len = cum_len + IM_LEN(in,1)
+		    }
+		}
+	case 2:						# join lines
+	    do band = 1, IM_LEN(out,3)
+		do image = 1, nimages {
+		    in = inptr[image]
+		    do line = 1, IM_LEN(in,2) {
+			stat = impnl$t (out, outbuf, vout)
+			stat = imgnl$t (in, inbuf, Meml[VPTR(vin,image)])
+			call amov$t (Mem$t[inbuf], Mem$t[outbuf], IM_LEN(in,1))
+		    }
+		}
+	case 3,4,5,6,7:					# join bands or higher
+	    do image = 1, nimages {
+		in = inptr[image]
+		do dim7 = 1, IM_LEN(in,7)
+		    do dim6 = 1, IM_LEN(in,6)
+			do dim5 = 1, IM_LEN(in,5)
+			    do dim4 = 1, IM_LEN(in,4)
+				do band = 1, IM_LEN(in,3)
+				    do line = 1, IM_LEN(in,2) {
+					stat = impnl$t (out, outbuf, vout)
+					stat = imgnl$t (in, inbuf,
+					    Meml[VPTR(vin,image)])
+					call amov$t (Mem$t[inbuf],
+					    Mem$t[outbuf], IM_LEN(in,1))
+				    }
+		# Unmap last image to free resources.
+		call imunmap (in)
+	    }
+	}
+
+	call sfree (sp)
+end
+
+$endfor