using Reactive,Interact
using GaussianProcess
using GPT_SGLD
using PyPlot

close()
length_scale=0.5; # length_scale for RBF kernel of GP
xylim=3; # range of x,y is (-xylim,xylilm)
resolution=50; # number of samples per dimension

x,y,X=createmesh(-xylim,xylim,resolution); # X=[x[1] y[1]; x[1] y[2]; ...]
xgrid=float([i for i in x, j in y])
ygrid=float([j for i in x, j in y])

figure(figsize=(10,2))
cov=SECov(length_scale,1);
mygp=GP(0,cov,2); #2D example
for i=1:5
    z=GPrand(mygp,X); 
    z=reshape(z,resolution,resolution)' #z[i,j]=value of GP at (x[i],y[j])
    subplot(1,5,i)
    pcolormesh(xgrid,ygrid,z,cmap=ColorMap("coolwarm"),vmin=-3,vmax=3)
    if i>1
        ax=gca()
        setp(ax[:get_xticklabels](),visible=false)
        setp(ax[:get_yticklabels](),visible=false)
    end
    #colorbar()
end
suptitle("GP draws")
tight_layout()

f=figure(figsize=(10,2))
@manipulate for r=1:2, n=1:2, Q=1:2; withfig(f) do
        for i=1:5
            fhat=fhatdraw(X,n,length_scale,r,Q,i)
            fhat=reshape(fhat,resolution,resolution)'
            subplot(1,5,i)
            pcolormesh(xgrid,ygrid,fhat,cmap=ColorMap("coolwarm"),vmin=-3,vmax=3)
            if i>1
                ax=gca()
                setp(ax[:get_xticklabels](),visible=false)
                setp(ax[:get_yticklabels](),visible=false)
            end
            #colorbar()
        end
        suptitle("Tensor model draws")
        tight_layout()
    end
end

#### Need r < n, Q < r^2 for plot to update