// This function returns the value of a general quadratic classifier at
// a point defined by the vector x.  The populations are defined by
// w1 and w2, and the probabilities of the various classes are defined
// in the vector P.
 
function [r]=quad_disc (x, w1, w2, P)

   C1=cov (w1);
   C2=cov (w2);
   U1 = mean (w1, "r");
   U1=U1';
   U2 = mean (w2, "r");
   U2=U2';
   I1=inv(C1);
   I2=inv(C2);
   
   g1= -0.5 * x' * I1 * x;
   g2 = 0.5 * x' * I1 * U1;
   g3 = - 0.5 * U1' * I1 * U1;
   g4 = 0.5 * x' * I1 * U1;
   c1 = log (P(1)) - log (6.28) - 0.5 * log (det (C1));
   g = g1 + g2 + g3 + g4 + c1;

   h1= -0.5 * x' * I2 * x;
   h2 = 0.5 * x' * I2 * U2;
   h3 = -0.5 * U2' * I2 * U2;
   h4 = 0.5 * x' * I2 * U2;
   c2 = log (P(2)) - log (6.28) - 0.5 * log (det (C2));
   h = h1 + h2 + h3 + h4 + c2;
   
   r = g - h;