
% Generate data from a 2D Gaussian process
clear all
close all

Nsamples=100;
u=2;  % threshold

% Plotting options
plot_cov=0;
plot_cross=1;
plot_level=1;
plot_data=1;

gp1.cfun='sqexp';
gp1.N=64;
gp1.s=0.01; % Smoothness
gp1.sigma=1; % Variance
gp1=init_gp_1d(gp1);

if plot_cov
    % Plot covariance function
    figure
    plot(gp1.v,gp1.C(1,:),'Linewidth',3);
    set(gca,'Fontsize',18);
    xlabel('d');
    ylabel('r(d)');
    title('Covariance function');
end

if plot_data
    figure
    hold on
    set(gca,'Fontsize',18);
    xlabel('v');
    ylabel('x');
end
u=2;
if plot_level
    
end

for n=1:Nsamples,
    clf;
    plot([0 1],[u u],'k-','Linewidth',3);
    set(gca,'Fontsize',18);
    xlabel('v');
    ylabel('x');
    hold on
    
    z=randn(gp1.N,1);
    x=gp1.A*z;
    xmax=max(x);
    
    s=sign(x-u);
    s_prev=s(1:end-1);
    s_next=s(2:end);
    cross=s_prev+s_next;
    crossing=find(cross==0);
    Nc=length(crossing);
    if plot_cross
        % plot crossings
        for i=1:Nc,
            v=gp1.v(crossing(i));
            plot([v v],[0 xmax],'r-','Linewidth',3); 
        end
    end
    Nu=length(find(s(crossing)==-1));
    
    if plot_data
        plot(gp1.v,x,'Linewidth',3);
        axis([0 1 -2*u 2*u]);
        title(sprintf('Number of upcrossings =%d',Nu));
    end
    
    pause
end






