/* FFNet_Eigen.c
 *
 * Copyright (C) 1994-2002 David Weenink
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or (at
 * your option) any later version.
 *
 * This program is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
 */

/*
 djmw 20020712 GPL header
*/

#include "FFNet_Eigen.h"
#include "Graphics.h"

void FFNet_Eigen_drawIntersection (FFNet me, Eigen eigen, Graphics g, long pcx, long pcy,
    double xmin, double xmax, double ymin, double ymax)
{
    long i, ix = labs (pcx), iy = labs (pcy); double x1, x2, y1, y2;
    long numberOfEigenvalues = eigen -> numberOfEigenvalues;
    long dimension = eigen -> dimension;
    
    if (ix > numberOfEigenvalues ||
    	iy > numberOfEigenvalues || my nInputs != dimension) return;
    Melder_assert (ix > 0 && iy > 0);
    if (xmax <= xmin || ymax <= ymin) Graphics_inqWindow (g, & x1, & x2, & y1, & y2);
    if (xmax <= xmin) { xmin = x1; xmax = x2; }
    if (ymax <= ymin) { ymin = y1; ymax = y2; }
    Graphics_setInner (g);
    Graphics_setWindow (g, xmin, xmax, ymin, ymax);
    for (i=1; i <= my nUnitsInLayer[1]; i++)
    {
		long j, unitOffset = my nInputs + 1;
		double c1 = 0.0, c2 = 0.0, bias = my w[ my wLast[unitOffset+i] ];
		double x[6], y[6], xs[3], ys[3]; int ns = 0;
		for (j=1; j <= my nInputs; j++)
		{
	    	c1 += my w[ my wFirst[unitOffset+i] + j - 1 ] * eigen->eigenvectors[ix][j];
	    	c2 += my w[ my wFirst[unitOffset+i] + j - 1 ] * eigen->eigenvectors[iy][j];
		}
		x[1] = x[2] = x[5] = xmin; x[3] = x[4] = xmax;
		y[1] = y[4] = y[5] = ymin; y[2] = y[3] = ymax;
		for (j=1; j <= 4; j++)
		{
	    	double p1 = c1 * x[j  ] + c2 * y[j  ] + bias;
	    	double p2 = c1 * x[j+1] + c2 * y[j+1] + bias;
	    	double r = fabs (p1) / ( fabs (p1) + fabs (p2));
	    	if (p1*p2 > 0 || r == 0.0) continue;
	    	if (++ns > 2) break;
	    	xs[ns] = x[j] + ( x[j+1] - x[j]) * r;
	    	ys[ns] = y[j] + ( y[j+1] - y[j]) * r;
		}
		if (ns < 2) Melder_casual ("intersection for unit %ld outside range", i);
		else Graphics_line (g, xs[1], ys[1], xs[2], ys[2]);
    }
    Graphics_unsetInner (g);
}

/* End of file FFNet_Eigen.c */