function [y,x] = threenorm (x, opts, output_mgn)
%THREENORM  The threenorm target in Breiman'96 [TR-460] and
%           Breiman'98 [Arcing Classifiers].
%
%   The default dimension is 20, so
%   [Y,X] = THREENORM(500) gives 500 20-D examples.
%
%   Y = THREENORM(X) gives output (1/-1) with dimension from X.
%   Y = THREENORM(X, .6) also specifies the probability of Class 1.

if nargin < 3, output_mgn = false; end
if nargin < 2, opts = []; end
opts = options([0.5, 20], opts);
p1 = opts(1); dim = opts(2);
if length(x) > 1, dim = size(x, 2); end

if output_mgn, error('Margin part is not implemented.'); end

a = 2 / sqrt(dim);
if length(x) == 1, %% Generate random inputs
    assert(x == round(x));
    r = rand(x, 1);
    i1 = (r < p1); i2 = ~i1;
    i11 = (r < p1 / 2); i12 = ~i11 & i1;
    y = ones(x, 1); y(i2) = -1;
    x = randn(x, dim);
    x(i11,:) = x(i11,:) + a;
    x(i12,:) = x(i12,:) - a;
    x(i2,1:2:end) = x(i2,1:2:end) + a;
    x(i2,2:2:end) = x(i2,2:2:end) - a;
else
    d11 = x - a; d11 = sum(d11.^2, 2);
    d12 = x + a; d12 = sum(d12.^2, 2);
    d2 = x - repmat(a * (mod(1:dim,2)*2-1), size(x,1), 1);
    d2 = sum(d2.^2, 2);
    post =  p1 ./ (p1 + 2*(1-p1) ./ (exp((d2-d11)/2) + exp((d2-d12)/2)));
    y = (rand(size(x,1),1) < post) * 2 - 1;
end