samer@0: /*
samer@0:  *	Copyright (c) 2002, Samer Abdallah, King's College London.
samer@0:  *	All rights reserved.
samer@0:  *
samer@0:  *	This software is provided AS iS and WITHOUT ANY WARRANTY;
samer@0:  *	without even the implied warranty of MERCHANTABILITY or
samer@0:  *	FITNESS FOR A PARTICULAR PURPOSE.
samer@0:  */
samer@0: 
samer@0: package samer.models;
samer@0: 
samer@0: import samer.core.*;
samer@0: import samer.core.types.*;
samer@0: import samer.maths.*;
samer@0: import samer.maths.opt.*;
samer@0: import samer.tools.*;
samer@0: 
samer@0: /**
samer@0: 	Differential scaler: scales and offsets each element of a vector
samer@0: 	independently aiming to match a given prior model. This is
samer@0: 	like ICA using a diagonal weight matrix, with a 'zero-mean'
samer@0: 	normalisation built in, though it doesn't actually the mean to
samer@0: 	do this, but some statistic based on the prior model.
samer@0: */
samer@0: 
samer@0: public class AlignedGaussian extends AnonymousTask implements Model
samer@0: {
samer@0: 	private int			n;
samer@0: 	private Vec			x;
samer@0: 	private VVector	s;
samer@0: 	private VVector	w;		// vector of multipliers
samer@0: 	private VDouble	E;
samer@0: 	private VDouble	logA;
samer@0: 
samer@0: 	double []		_x, _s, _g, _w, _m, phi;
samer@0: 
samer@0: 	public AlignedGaussian( Vec input) { this(input.size()); setInput(input); }
samer@0: 	public AlignedGaussian( int N)
samer@0: 	{
samer@0: 		n = N;
samer@0: 
samer@0: 		x = null;
samer@0: 		w = new VVector("w",n);		w.addSaver();
samer@0: 		s = new VVector("s",n);
samer@0: 		E = new VDouble("E");
samer@0: 		logA = new VDouble("log|A|",0);
samer@0: 		
samer@0: 
samer@0: 		_s = s.array();
samer@0: 		_w = w.array();
samer@0: 		_g = new double[n];
samer@0: 		phi = null;
samer@0: 		reset();
samer@0: 	}
samer@0: 
samer@0: 	public int getSize() { return n; }
samer@0: 	public VVector output() { return s; }
samer@0: 	public VVector weights() { return w; }
samer@0: 	public void setInput(Vec in) { x=in; _x=x.array(); }
samer@0: 	public void reset() { reset(1.0); }
samer@0: 	public void reset(double k) {
samer@0: 		Mathx.setAll(_w,k); w.changed();
samer@0: 		logA.set(-sumlog(_w));
samer@0: 	}
samer@0: 
samer@0: 	public String toString() { return "AlignedGaussian("+x+")"; }
samer@0: 	private static double sumlog(double [] x) {
samer@0: 		double S=0;
samer@0: 		for (int i=0; i<x.length; i++) S += Math.log(x[i]);
samer@0: 		return S;
samer@0: 	}
samer@0: 
samer@0: 	public void dispose()
samer@0: 	{
samer@0: 		logA.dispose();
samer@0: 		w.dispose();
samer@0: 		s.dispose();
samer@0: 		super.dispose();
samer@0: 	}
samer@0: 
samer@0: 	public void infer() { Mathx.mul(_s,_x,_w); s.changed(); }
samer@0: 	public void compute() { 
samer@0: 		Mathx.mul(_g,_s,_w); 
samer@0: 		E.set(0.5*Mathx.norm2(_s)+logA.value);
samer@0: 	}
samer@0: 
samer@0: 	public double	getEnergy() { return E.value; }
samer@0: 	public double [] getGradient() { return _g; }
samer@0: 
samer@0: 	public Functionx functionx() {
samer@0: 		return new Functionx() {
samer@0: 			double [] s=new double[n];
samer@0: 
samer@0: 			public void dispose() {}
samer@0: 			public void evaluate(Datum P) { P.f=evaluate(P.x,P.g); }
samer@0: 			public double evaluate(double [] x, double [] g) {
samer@0: 				Mathx.mul(s,x,_w);
samer@0: 				Mathx.mul(g,s,_w);
samer@0: 				return 0.5*Mathx.norm2(s); // +logA.value;
samer@0: 			}
samer@0: 		};
samer@0: 	}
samer@0: 
samer@0: 	public void starting() { logA.set(-sumlog(_w)); }
samer@0: 	public void stopping() {}
samer@0: 	public void run() { infer(); compute(); }
samer@0: }
samer@0: