/*Copyright (C) 1993, 1996 by Thomas Glen Smith.  All Rights Reserved.*/
/* reduces.h APL2 V1.0.0 ***********************************************
* The is the heart of both subroutines reducesb.c and ireduces.c.      *
* After using #define for the following values, this code is included  *
* in those subroutines to do both reduce and scan processing for the   *
* two data types integer and real, respectively:                       *
*   SUBRTNE -> Subroutine name.                                        *
*   DATATYPE-> Data type, i.e. int or double.                          *
*   APLTYPE -> APL data type, i.e. APLINT or APLNUMB.                  *
*   CONVERT -> Conversion subroutine, i.e. integer or real.            *
*   COPYSUB -> Copy subroutine, i.e. intcopy or dblcopy.               *
*   DATAPTR -> union name to be used with aplptr, i.e. aplint, apldata.*
***********************************************************************/
struct aplcb *SUBRTNE(id,oper,identity,rite,axis)
int id; /* 1=reduce, 0=scan */
DATATYPE (*oper)(); /* operator */
DATATYPE *identity; /* identity value */
Aplcb rite; /* operand */
int axis; /* axis of reduction */
{
	Errinit; Errstop; Reducecm;
	extern int aplerr;
	Aplcb CONVERT(), out;
	int axicnt,botcnt,topcnt;
	int i,j,k,m,n,p;
	DATATYPE *dataout,*COPYSUB(),*ip,*kp,wrk;

	if (errinit())
		return(errstop(0,NULL,rite,NULL));
	if (!(rite->aplflags & APLTYPE)) {
		rite=CONVERT(rite); /* convert to desired input data type */
		if (aplerr)
			return(NULL);
	}
	out=reducecm(id,identity,rite,&axis,&axicnt,&botcnt,&topcnt,APLTYPE);
	if (aplerr)
		return(NULL);
	n = (id) ? 1 : axicnt; /* n == 1 if reduce, axicnt if scan */
	if (out->aplcount) { /* 1 or more elements of output */
		dataout = out->aplptr.DATAPTR;
		if (0 == rite->aplcount) /* is input empty? */
			dataout=COPYSUB(dataout,identity,out->aplcount,0); 
		else {
			for (i = 0; i < topcnt; i++) {
				ip=rite->aplptr.DATAPTR+(p=i*botcnt*axicnt);
				for (j=0; j<botcnt; j++)
					for (m = n; m > 0; m--) {
						wrk=*(kp=ip+j+(axicnt-m)*botcnt);
						for (k = 1; k < axicnt-m+1; k++)
							wrk=(*oper)(*(kp-=botcnt),wrk);
						if (id) *dataout++=wrk; /* reduce */
						else *(dataout+p+j+(axicnt-m)*botcnt)=wrk;
					}
			}
		}
	}
	return(errstop(0,NULL,rite,out));
}