text:
% general script for non-systematic LDPC simulation

% Things to Remember - Steps to follow ALWAYS !!!!!!!!
% 0. This file is good for non-systematic H only (errors calculated for all n)
% 1. Load H
% 2. Define SNR Range
% 3. Set Maximum number of iterations
% 4. Set Maximum number of codeword-errors for which to run simulation
% 5. Select decoder - Matlab or C-based

clear all
close all
clc
tic

%get H from somewhere , change R accordingly !!
load 128x256regular.mat H
% load 128x256regular_v6.mat  %for users of Matlab 6.5

ind=find(H==1);
[r,c]=ind2sub(size(H),ind);
[rows,cols]=size(H);
h=sparse(H);                         % for use with Matlab-based LDPC Decoder
n=cols;
k=n-rows;

% Find
% 1: maximum check degree
% 2: column indeces in each row which contain '1'
% 3: maximum variable degree
% 4: find column indeces in each row which contain '1'
[max_check_degree,check_node_ones,BIGVALUE_COLS,max_variable_degree,variable_node_ones,BIGVALUE_ROWS]=one_finder(H);

rand('seed',584);
randn('seed',843);

dB=[0:3];                           % range of SNR values in dB
SNRpbit=10.^(dB/10);                % Eb/No conversion from dB to decimal
No_uncoded=1./SNRpbit;              % since Eb=1
R=k/n;                              % code rate
No=No_uncoded./R;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
max_iter=10;                   %maximum number of decoder iterations
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
maximum_blockerror=30;          % maximum blockerrors per SNR point
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

blockerrors=0;
biterrors=0;
block=0;
FER=zeros(1,length(dB));        % array for Frame Error Rate
BER=zeros(1,length(dB));        % array for Channel Error Rate
block_array=zeros(1,length(dB));

for z=1:length(SNRpbit), % loop for testing over range of SNR values

    biterrors=0;
    blockerrors=0;
    block=0;

    while(blockerrors<maximum_blockerror)  %while loop

        %%%%%%%%%%%%%%% u is the codeword to be transmitted %%%%%%%%%%%%%%%%%%%
        u=zeros(1,cols);  %%all zero codeword
        tx_waveform=bpsk(u);

        sigma=sqrt(No(z)/2);
        rx_waveform=tx_waveform + sigma*randn(1,length(tx_waveform));

        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        %%%%%%%%%%%%%%%%% Use C - based LDPC Decoder %%%%%%%%%%%%%%
        %%% Fastest C-based LDPC Decoder  %%%%

        gamma_n=(4/No(z))*rx_waveform;
        vhat=decode_ldpc_new(max_iter,gamma_n,check_node_ones,max_check_degree,BIGVALUE_COLS-1,variable_node_ones,max_variable_degree,BIGVALUE_ROWS-1,rows,cols);
        uhat=vhat;

        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        %%% Slightly slower C-based LDPC Decoder  %%%%

%         gamma_n=(4/No(z))*rx_waveform;
%         sg=zeros(rows,cols);
%         sg(ind)=gamma_n(c);
%
%         [vhat,iteration]=decode_ldpc(No(z),H,max_iter,sg,gamma_n);
%         uhat=vhat;
%
        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        %%%%%%%%%%%%%%%%% Use Matlab-based LDPC Decoder %%%%%%%%%%%%%%
%
%         [vhat,iteration]=decode_ldpc_matlab(rx_waveform,No(z),h,rows,cols,ind,r,c,max_iter);
%         uhat=vhat;

        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

        Errors=zeros(1,length(vhat));
        Errors(find(u~=vhat))=1;

        if sum(Errors)~=0
            blockerrors=blockerrors+1;

            if rem(blockerrors,5)==0; %save statistics after every 5 blockerrors

                ber=biterrors/(block*n);
                fer=blockerrors/block;

                s=sprintf('%dx%d_regular_results.txt',n-k,n);
                fid = fopen(s,'a');
                fprintf(fid,'\n');
                fprintf(fid,'%s %2.1EdB\n','SNR =',dB(z));
                fprintf(fid,'%s %6.3E\n','BER =',ber);
                fprintf(fid,'%s %6.3E\n','FER =',fer);
                fprintf(fid,'%s %d\n','blocks =',block);
                fprintf(fid,'%s %d\n','blockerrors =',blockerrors);
                fprintf(fid,'%s %d\n','biterrors =',biterrors);
                fclose(fid);
            end

        end

        biterrors=biterrors+sum(Errors);
        block=block+1;

    end     %while loop

    block_array(z)=block; %counting blocks per SNR point
    BER(z)=biterrors/(block*n);
    FER(z)=blockerrors/block;
    fprintf(1,'\n\n Simulation finished for SNR: %d \n',dB(z))

%     save results.mat dB BER FER block_array

end      % loop for testing over range of SNR values

semilogy(dB,BER,'b-o')
title('Bit Error Rate')
ylabel('BER')
xlabel('Eb/No (dB)')
grid
figure
semilogy(dB,FER,'b-o')
title('Frame Error Rate')
ylabel('FER')
xlabel('Eb/No (dB)')
grid

toc

%   Copyright (c) 2005 by Shaikh Faisal Zaheer, faisal.zaheer@gmail.com
%   $Revision: 3.0 $  $Date: 26/12/2005 $


text:
function y=flip(x)
%   function y=flip(x)
%   reverses the order of the elements of a vector, first->last, etc.
%   ddj 02 08 26
n=length(x);
for i=1:n
  y(i)=x(n+1-i);
end

author: CelsoReyes


text:
function [theta] = normalEqn(X, y)
%NORMALEQN Computes the closed-form solution to linear regression
%   NORMALEQN(X,y) computes the closed-form solution to linear
%   regression using the normal equations.

theta = zeros(size(X, 2), 1);

% ====================== YOUR CODE HERE ======================
% Instructions: Complete the code to compute the closed form solution
%               to linear regression and put the result in theta.
%

% ---------------------- Sample Solution ----------------------

theta = pinv(X' * X) * X' * y;

% -------------------------------------------------------------

% ============================================================

end

author: vugsus

repo: coursera-machine-learning


text:
function [opts]=test_rnn(net,opts)

    opts.training=0;

    opts.MiniBatchError=[];
    opts.MiniBatchLoss=[];

    for mini_b=1:opts.n_test_batch
        idx=1+(mini_b-1)*opts.parameters.batch_size:mini_b*opts.parameters.batch_size;
        opts.input_data=opts.test(:,idx,:);
        opts.input_labels=opts.test_labels(idx,:);

        %forward
        if strcmp(opts.network_name,'lstm')
            [ net,res,opts ] = lstm_ff( net,opts );
        end
        if strcmp(opts.network_name,'gru')
            [ net,res,opts ] = gru_ff( net,opts );
        end
        if strcmp(opts.network_name,'rnn')
            [ net,res,opts ] = rnn_ff( net,opts );
        end
        if strcmp(opts.network_name,'qrnn')
            [ net,res,opts ] = qrnn_ff( net,opts );
        end

        if isfield(opts,'err')
            opts.MiniBatchError=[opts.MiniBatchError;gather( opts.err(1))];
        end
        opts.MiniBatchLoss=[opts.MiniBatchLoss;gather( opts.loss)];

    end

    opts.results.TestEpochError=[opts.results.TestEpochError;mean(opts.MiniBatchError(:))];
    opts.results.TestEpochLoss=[opts.results.TestEpochLoss;mean(opts.MiniBatchLoss(:))];

end

author: yechengxi

repo: LightNet

save_path: github-repos/MATLAB/yechengxi-LightNet


text:
function h = plus(f, g)
%+   Plus for SEPARABLEAPPROX objects.
%
% F + G adds F and G. F and G can be scalars or SEPARABLEAPPROX objects.

% Copyright 2017 by The University of Oxford and The Chebfun Developers.
% See http://www.chebfun.org/ for Chebfun information.

if ( ~isa(f, 'separableApprox') ) % ??? + SEPARABLEAPPROX

    h = plus(g, f);

elseif ( isempty(g) ) % SEPARABLEAPPROX + []

    h = f;

elseif ( isempty(f) ) % [] + SEPARABLEAPPROX

    h = g;

elseif ( isa( g, 'double' ) )           % SEPARABLEAPPROX + DOUBLE

    g = compose( 0*f,@plus, g);   % promote double to object class.
    h = plus(f, g);

elseif ( ~isa(g, 'separableApprox') )          % SEPARABLEAPPROX + ???

    error( 'CHEBFUN:SEPARABLEAPPROX:plus:unknown', ...
        ['Undefined function ''plus'' for input arguments of type %s ' ...
        'and %s.'], class(f), class(g));

else                                     % SEPARABLEAPPROX + SEPARABLEAPPROX

    % Domain Check:
    if ( ~domainCheck(f, g) )
        error('CHEBFUN:SEPARABLEAPPROX:plus:domain', 'Inconsistent domains.');
    end

    % Type check:
    if ( ~strcmp( class(f),class(g) ) )
        error( 'CHEBFUN:SEPARABLEAPPROX:plus:unknown', ...
            ['Undefined function ''plus'' for input arguments of type %s ' ...
            'and %s. Try converting to the same type.'], class(f), class(g));
    end

    % Check for zero SEPARABLEAPPROX objects:
    if ( iszero(f) )
        h = g;
    elseif ( iszero(g) )
        h = f;
    else
        % Add together two nonzero SEPARABLEAPPROX objects:
        h = compression_plus(f, g);

    end

end

end

function h = compression_plus(f, g)
% Add SEPARABLEAPPROX objects together by a compression algorithm.

% The algorithm is as follows:
% If A = XY^T and B = WZ^T, then A + B = [X W]*[Y Z]^T,
% [Qleft, Rleft] = qr([X W])
% [Qright, Rright] = qr([Y Z])
% A + B = Qleft * (Rleft * Rright') * Qright'
% [U, S, V] = svd( Rleft * Rright' )
% A + B = (Qleft * U) * S * (V' * Qright')     -> new low rank representation

% Hack: Ensure g has the smaller pivot values.
if ( norm(f.pivotValues, -inf) < norm(g.pivotValues, -inf) )
    % [TODO]: Understand why this works!
    h = compression_plus(g, f);
    return
end

fScl = diag(1./f.pivotValues);
gScl = diag(1./g.pivotValues);
cols = [f.cols, g.cols];
rows = [f.rows, g.rows];

[Qcols, Rcols] = qr(cols);
[Qrows, Rrows] = qr(rows);

Z = zeros(length(fScl), length(gScl));
D = [ fScl, Z ; Z.', gScl ];
[U, S, V] = svd(Rcols * D * Rrows.');
% Take diagonal from SIGMA:
s = diag(S);

% Compress the format if possible.
% [TODO]: What should EPS be in the tolerance check below?
vf = vscale(f);
vg = vscale(g);

vscl = 2*max(vf, vg);
% Remove singular values that fall below eps*vscale:
idx = find( s > 10*eps * vscl, 1, 'last');

if ( isempty(idx) )
    % Return 0 separableApprox
    h = 0*f;
else
    U = U(:,1:idx);
    V = V(:,1:idx);
    s = s(1:idx);
    h = f;
    h.cols = Qcols * U;
    h.rows = Qrows * conj(V);
    % [TODO]: PivotValues have very little meaning after this compression step.
    % For now we assign the singular values as the pivot values.
    h.pivotValues = 1./s;
end

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
% This functions, when given the detectionsQ and detectionsT
% (after propagation step) gives the 3D3D and 3D2D intesection
% area score
function [scoreMatrix_3D3D, scoreMatrix_3D2D] = generateScoreMatrices(detectionsQ, detectionsT)

    % score_3D2D   = get3D2DMatchScore(detectionsT(1).bbox, detectionsQ(1).bvolume_proj);
    % score_3D3D = get3D3DMatchScore(detectionsT(2).bvolume, detectionsQ(1).bvolume);

    scoreMatrix_3D3D = zeros(length(detectionsQ), length(detectionsT));
    scoreMatrix_3D2D = zeros(length(detectionsQ), length(detectionsT));

    % 3D 3D score
    for i = 1:length(detectionsQ)
        for j = 1:length(detectionsT)

            % remember all the shapes are in F2, and we are matching shapes
            % i.e. 3D bvolumes and their bboxes in F2 with their
            % corresponding search bounding volumes and bounding volume
            % projections which are the propagated volumes of cars in F1.
            % Hence here Q is T and vice-versa
            scoreMatrix_3D3D(i,j) = get3D3DMatchScore(detectionsT{j}.bvolume, detectionsQ{i}.bvolume);
            scoreMatrix_3D2D(i,j) = get3D2DMatchScore(detectionsT{j}.bbox, detectionsQ{i}.bvolume_proj);
        end
    end

end

repo: MOTBeyondPixels


text:
function [cellU, S, cellV] = separableFormat(A, xorder, yorder, domain)
%SEPARATBLEFORMAT  Compute separable expression for a linear PDO.
%
% Calculate a separable representation of a partial differential
% operator. These representations can then be using to derive a 2D spectral
% method from 1D ideas.  The linear PDO can have variable coefficients.
%
% This uses the tensor-train decomposition and Proposition
% 4.2 from [1].
%
% [1] A. Townsend and S. Olver, The automatic solution of partial differential
% equations using a global spectral method, submitted, 2014.
%
% Author: Alex Townsend September 2014.

if ( nargin == 1 )
    N = A;
    A = N.coeffs;
    xorder = N.xorder;
    yorder = N.yorder;
    domain = N.domain;
end

% Loop over coefficients of A. Find coefficient of highest degree. We will
% need this to recover the variable coefficients:
n = 10;
for jj = 1:size(A, 1)
    for kk = 1:size(A, 2);
        if ( isa(A{jj,kk}, 'chebfun2') )
            [xdeg, ydeg] = length(A{jj,kk});    % get degrees
            n = max([xdeg, ydeg, n]) + 1;       % take maximum degree we find
        end
    end
end

% Set up Chebyshev points that we need:
x = chebpts(xorder+1);
s = chebpts(n, domain(1:2));
y = chebpts(yorder+1);
t = chebpts(n, domain(3:4));
[xx, ss, yy, tt] = ndgrid(x, s, y, t);
[newx, newy] = meshgrid(s, t);

% We need to apply Proposition 4.2 from [1]. Use the linear operator
% T motivated by umbral calculus. That is, convert
%
%       d^j/dx^j-> x^j
%       d^j/dy^j-> y^j
%         x     -> s
%         y     -> t
% We obtain a function of 4 variables, H(x,s,y,t). See [1].
H = @(x,s,y,t) 0*x;
for jj = 1:size(A, 1)
    for kk = 1:size(A, 2)
        if ( isa(A{jj,kk}, 'double') && ~(A{jj,kk} == 0) )
            H = @(x,s,y,t) H(x,s,y,t) + A{jj,kk}*x.^(kk - 1).*y.^(jj - 1);
        elseif ( isa(A{jj,kk}, 'chebfun2') )
            v = zeros(1, n, 1, n);
            v(1,:,1,:) = feval(A{jj,kk}, newx, newy).';
            out = repmat(v, [xorder + 1, 1, yorder + 1, 1]);
            H = @(x,s,y,t) H(x,s,y,t) + out.*x.^(kk-1).*y.^(jj-1);
        end
    end
end
H = H(xx, ss, yy, tt);

% Using tensor-train ideas. Calculate the splitting rank of the function
% H(x,s,y,t):
A = reshape(H, n*(xorder+1), n*(yorder+1));
[U, S, V] = svd(A);
rk = find(abs(diag(S)/S(1,1)) > 1000*eps, 1, 'last' );  % splitting rank

% Restrict to singular vectors of interest.
S = S(1:rk, 1:rk);
U = U(:, 1:rk);
V = V(:, 1:rk);

% We have the splitting rank of the PDO, now we want the corresponding
% separable representation. The following is tricky to get right...
cellU = cell(yorder+1, rk);
cellV = cell(xorder+1, rk);
c1 = cell(xorder, 1);
c2 = cell(yorder, 1);
% Matrices to convert ChebT -> monomials:
converty = fliplr(poly(chebpoly(0:yorder)));
convertx = fliplr(poly(chebpoly(0:xorder)));

% Figure out the separable representation:
for jj = 1:rk

    % This is giving us the 1D ODEs that go on the right in the generalized
    % Sylvester matrix equation:
    f1 = chebfun2.vals2coeffs( reshape(U(:,jj), xorder+1, n) );  % @(x, s)
    f1 = f1.' * convertx;
    for kk = 1:xorder+1
        c1{kk} = chebfun(f1(:,kk), domain(1:2), 'coeffs');
        cellV(kk,jj) = c1(kk);
    end

    % This is giving us the 1D ODEs that go on the left in the generalized
    % Sylvester matrix equation:
    f2 = chebfun2.vals2coeffs( reshape(conj(V(:,jj)), yorder+1, n) );  % @(y, t)
    f2 = f2.' * converty;
    for kk = 1:yorder+1
        c2{kk} = chebfun(f2(:,kk), domain(3:4), 'coeffs');
        cellU(kk,jj) = c2(kk);
    end
end

end


text:
function precision = Precision(SEG, GT)
    % SEG, GT are the binary segmentation and ground truth areas, respectively.
    % precision
    precision = double(sum(uint8(SEG(:) & GT(:)))) / double(sum(uint8(SEG(:))));
end

repo: cluster


text:
function info = muteinf(A, Y, method)
	n = size(A,1);
	Z = [A Y];
	if(n/10 > 20)
		nbins = 20;
	else
		nbins = max(floor(n/10),10);
	end;
	pA = hist(A, nbins);
	pA = pA ./ n;

	i = find(pA == 0); % that's a hack!
	pA(i) = 0.00001;
	if(strcmp(method,'regression'))
		pY = hist(Y, nbins);
		pY = pY ./ n;
		j = find(pY == 0);
		pY(j) = 0.00001;

		rx = abs(max(A) - min(A)) / nbins;
		ry = abs(max(Y) - min(Y)) / nbins;
		yl = min(Y);
		p = zeros(nbins, nbins);
		for i = 1:nbins
			xl = min(A);
			for j = 1:nbins
				%disp(['intervals [' num2str(xl) ',' num2str(xl+rx) '|' num2str(yl) ',' num2str(yl+ry) ']'])
				interval = (xl <= Z(:,1)) & (yl <= Z(:,2));
				if(j < nbins)
					interval = interval & (Z(:,1) < xl + rx);
				end;
				if(i < nbins)
					interval = interval & (Z(:,2) < yl + ry);
				end;
				%find(interval)
				p(i,j) = length(find(interval));

				if p(i,j) == 0 % hack!
					p(i,j) = 0.00001;
				end

				xl = xl + rx;
			end;
			yl = yl + ry;
		end;
		HA = -sum(pA .* log(pA));
		HY = -sum(pY .* log(pY));
		pA = repmat(pA,nbins,1);
		pY = repmat(pY',1,nbins);
	else
		od = size(Y,2);
		cl = od;
		if(od == 1)
			pY = [length(find(Y==+1)) length(find(Y==-1))] / n;
			cl = 2;
		else
			pY = zeros(1,od);
			for i=1:od
				pY(i) = length(find(Y==+1));
			end;
			pY = pY / n;
		end;
		p = zeros(cl,nbins);
		rx = abs(max(A) - min(A)) / nbins;
		for i = 1:cl
			xl = min(A);
			for j = 1:nbins
				if(i == 2) & (od == 1)
					interval = (xl <= Z(:,1)) & (Z(:,2) == -1);
				else
					interval = (xl <= Z(:,1)) & (Z(:,i+1) == +1);
				end;
				if(j < nbins)
					interval = interval & (Z(:,1) < xl + rx);
				end;
				%find(interval)
				p(i,j) = length(find(interval));

				if p(i,j) == 0 % hack!
					p(i,j) = 0.00001;
				end

				xl = xl + rx;
			end;
		end;
		HA = -sum(pA .* log(pA));
		HY = -sum(pY .* log(pY));
		pA = repmat(pA,cl,1);
		pY = repmat(pY',1,nbins);
	end;
	p = p ./ n;

	info = sum(sum(p .* log(p ./ (pA .* pY))));
	info = 2 * info ./ (HA + HY);

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
function [Fimage,objimage,x,y] = imagemodel(F,obj,options)
% IMAGEMODEL Explicitely removes equality constraints from model.
%
% [Fi,hi,x,y] = imagemodel(F,h)
%
% Input
%  F   : Constraint in form F(x)>0, Ax=b
%  h   : objective function h(x)
%
% Output
%  Fi  : Constraints in form F(y)>0
%  hi  : objective function h(y)
%  x   : original variables
%  z   : old variables expressed in basis y
%
% To obtain a solution in the original variables x, use assign(x,double(y))
%
% Note that this reduction is automatically done when you call solvesdp and use
% a solver that cannot handle equality constraints. Hence, there is typically no
% reason to use this command, unless some further manipulations are going
% to be done.
%
% See also DUALIZE, PRIMALIZE

% Check for unsupported problems
err = 0;
p1 = ~(isreal(F) & isreal(obj));
p2 = ~(islinear(F) & islinear(obj));
p3 = any(is(F,'integer')) | any(is(F,'binary'));
if p1 | p2 | p3
    if nargout == 5
        Fdual = ([]);objdual = [];y = []; X = []; t = []; err = 1;
    else
        problems = {'Cannot imagalize complex-valued problems','Cannot imagalize nonlinear problems','Cannot imagalize discrete problems'};
        error(problems{min(find([p1 p2 p3]))});
    end
end

if nargin < 3
    options = sdpsettings('solver','sedumi','remove',1);
end

if any(is(F,'equality'))
    [model,recoverdata,solver,diagnostic,F] = compileinterfacedata(F,[],[],obj,options,0);
    if isfield(diagnostic,'problem')
        if diagnostic.problem == 1
            warning('Problem is infeasible');
            Fimage = [];
            objimage = [];
            x = [];
            y = [];
            return
        end
    end
    if isempty(model)
            Fimage = [];
            objimage = [];
            x = [];
            y = [];
            warning('Reduced problem does not have free variables. Optimal solution computed');
            return
    end
else
    Fimage = F;
    objimage = obj;
    x  = recover(unique([depends(obj) depends(F)]));
    y = x;
    return
end

y = sdpvar(length(model.c),1);

vecF = model.F_struc*[1;y];
K = model.K;
Fimage = ([]);
start = 1;
if any(model.K.l)
    Fimage = Fimage + (vecF(start:start+K.l-1) >= 0);
    start = start + K.l;
end

if any(model.K.q)
    for i = 1:length(model.K.q)
        z = vecF(start:start+K.q(i)-1)
        Fimage = Fimage + (cone(z(2:end),z(1)));
        start = start + K.q(i);
    end
end

if any(model.K.s)
    for i = 1:length(model.K.s)
        z = vecF(start:start+K.s(i)^2-1);
        Fimage = Fimage + (reshape(z,K.s(i),K.s(i)));
        start = start + K.s(i)^2;
    end
end

objimage = model.c'*y+y'*model.Q*y+model.f;

x = recover(recoverdata.used_variables);
y = recoverdata.H*y+recoverdata.x_equ;

author: yalmip

repo: YALMIP

save_path: github-repos/MATLAB/yalmip-YALMIP


text:
% DeJong_f3.m
% De Jong's f3 function, ND, also called STEP
% from: http://www.cs.rpi.edu/~hornda/pres/node4.html
%
% f(x) = sum( floor(x) )
%
% x = N element row vector containing [ x0, x1,..., xN ]
%   each row is processed independently,
%   you can feed in matrices of timeXN no prob
%
% example: cost = DeJong_f3([1,2;3,4;5,6])
% note minimum @ x= -Inf or whatever lower bound you've chosen

% Brian Birge
% Rev 1.0
% 9/12/04

function [out]=DeJong_f3(in)
 out = sum( floor(in) , 2);

author: HuangCongQing


text:
classdef FeatureVocabularyForSelectFeatures < handle

  % instance variables
  properties
%     jld  % instance of JLabelData, a ref
%     scoreFeatures  % the scores-as-input structure, detailing the files that
%                    % hold classifiers to be used as inputs, and their
%                    % per-frame feature names
%     windowFeatureParams
%    featureLexicon  % the dictionary of possible per-frame features
    wfParamsFromAmount  % scalar structure holding the window feature parameters
                        % that go with each of the possible "amount" levels.
                        % e.g. 'normal', 'more', 'less'.
    pfCategoriesFromName  % scalar structure.  Each field is the name of a
                          % per-frame feature, and the value is a cell
                          % array of strings, telling what categories that
                          % per-frame feature is in.  Most PFs are in only
                          % one category, but some are in more than one.
    pfCategoryNames  % a cell array of strings, holding the names of all the
                     % possible per-frame feature categories
    subdialectPFNames
      % A cell array of strings, containing the names of all the
      % to-be-calculated per-frame features.  This includes the
      % scores-as-input PFs.  I like to call this the "subdialect".
    defaultPFTransTypesFromName
      % A structure, each field is a per-frame feature, and each value is a cell
      % array of strings holding the default translation types
      % for that per-frame feature  (formerly transType)
    vocabulary
      % a structure where each field is a per-frame feauture, and the value
      % holds the window feature parameters for that PF.
      % Should have as many fields as there are per-frame features in the
      % to-be-calculated vocabulary (the subdialect).  vocabulary.(perframeFeatureName).enabled is a boolean
      % indicating whether that per-frame feature is in the vocabulary.
      % (formerly data)
  end  % properties

  % things that look like instance vars from the outside, but
  % which are calculated on-the-fly
  properties (Dependent=true)
    wfAmounts
  end

  % class constants
  properties (Access=public,Constant=true)
    % Window features come in different types.  Below is a list of the
    % different types.  (formerly windowComp)
    wfTypes = {'default','mean','min','max','hist','prctile',...
               'change','std','harmonic','diff_neighbor_mean',...
               'diff_neighbor_min','diff_neighbor_max','zscore_neighbors'};
    % Window features have parameters.  Below are the names of those parameters.
    % These parameters exist for each WF type.  (They should arguably be
    % called the "standard" parameters, or the "common" parameters,
    % because they exist for all WF types.  Whereas the extra parameters,
    % below, are each particular to one WF type.
    % (Formerly winParams.)
    wfParamNames = {'max_window_radius','min_window_radius','nwindow_radii', ...
                    'trans_types','window_offsets'};
%     % These are the default values for the window feature parameters.  (Formerly defaultWinParams.)
%     % They are only used if the lexicon fails to specify pre-set
%     % window-feature amounts (e.g. normal, more, less)
%     defaultWFParams = {10,1,3,{'none'},0};
    % A Window feature type can sometimes have an extra parameter.  Below are the names
    % of those parameters, one for each WF type.  WF types without an extra param have the
    % empty string for their entry.  (Formerly winextraParams.)
    wfExtraParamNames = {'','','','','hist_edges','prctiles','change_window_radii', ...
                         '','num_harmonic','','','',''};
%     % Below are the default values for the extra window feature parameters.  (Formerly winextraDefaultParams.)
%     % They are only used if the lexicon fails to specify pre-set
%     % window-feature amounts (e.g. normal, more, less)
%     defaultWFExtraParams = {[],[],[],[],[-400000 0 40000],[5 10 30 50 70 90 95],[1 3],...
%                             [],2,[],[],[],[]};
  end  % class constants

  % class methods
  methods (Access=public,Static=true)
    % ---------------------------------------------------------------------
    function wfParamsFromAmount=computeWFParamsFromAmount(featureLexicon)
      % Specify some constants we need
      % Below are the default values for the window feature parameters.  (Formerly defaultWinParams.)
      % They are only used if the lexicon fails to specify pre-set
      % window-feature amounts (e.g. normal, more, less)
      defaultWFParams = {10,1,3,{'none'},0};
      % Below are the default values for the extra window feature parameters.  (Formerly winextraDefaultParams.)
      % They are only used if the lexicon fails to specify pre-set
      % window-feature amounts (e.g. normal, more, less)
      defaultWFExtraParams = {[],[],[],[],[-400000 0 40000],[5 10 30 50 70 90 95],[1 3],...
                              [],2,[],[],[],[]};

      % Get some class constants we'll need
      wfTypes = ...
        FeatureVocabularyForSelectFeatures.wfTypes;
      wfParamNames = ...
        FeatureVocabularyForSelectFeatures.wfParamNames;
      wfExtraParamNames = ...
        FeatureVocabularyForSelectFeatures.wfExtraParamNames;
      % Read the default parameters for different categories.
      wfAmounts = fieldnames(featureLexicon.defaults);
      if isempty(wfAmounts)
        % If no preset amounts have been specified in the lexicon file,
        % fall back on the built-in defaults.
        wfAmount='default';
        wfParams = struct;
        % Set the parameters for the WFs of type 'default'
        wfParams.default.enabled = true;
        for ndx = 1:numel(wfParamNames)
          wfParamName = wfParamNames{ndx};
          wfParams.default.values.(wfParamName) = defaultWFParams{ndx};
        end
        wfParams.default.values.sanitycheck = false;
        % Set the parameters for the other WF types
        for ndx = 2:numel(wfTypes)
          wfType = wfTypes{ndx};
          wfParams.(wfType).enabled = false;
          for paramIndex = 1:numel(wfParamNames)
            wfParamName = wfParamNames{paramIndex};
            wfParams.(wfType).values.(wfParamName) = defaultWFParams{paramIndex};
          end
          if ~isempty(wfExtraParamNames{ndx})
            extraParamName = wfExtraParamNames{ndx};
            wfParams.(wfType).values.(extraParamName) = defaultWFExtraParams{ndx};
          end
        end
        wfParamsFromAmount.(wfAmount) = wfParams;
      else
        % The usual case --- preset amounts are defined in the lexicon.
        for wfAmountIndex = 1:numel(wfAmounts)
          wfAmount=wfAmounts{wfAmountIndex};
          wfParams = struct;
          wfParamsThisAmountRaw = featureLexicon.defaults.(wfAmount);
          % Set the parameters for the WFs of type 'default'
          wfParams.default.enabled = true;
          for ndx = 1:numel(wfParamNames)
            wfParamName = wfParamNames{ndx};
            wfParams.default.values.(wfParamName) = wfParamsThisAmountRaw.(wfParamName);
          end
          wfParams.default.values.sanitycheck = false;
          % Set the parameters for the other WF types
          for ndx = 2:numel(wfTypes)
            % Copy the standard WF params for this WF type
            wfType = wfTypes{ndx};
            for paramIndex = 1:numel(wfParamNames)
              wfParamName = wfParamNames{paramIndex};
              if isfield(wfParamsThisAmountRaw,wfType) && isfield(wfParamsThisAmountRaw.(wfType),wfParamName)
                wfParams.(wfType).values.(wfParamName) = wfParamsThisAmountRaw.(wfType).(wfParamName);
              else
                % Use the common value
                wfParams.(wfType).values.(wfParamName) = wfParamsThisAmountRaw.(wfParamName);
              end
            end
            % Copy the extra WF parameter, if present for this WF type
            if ~isempty(wfExtraParamNames{ndx})
              extraParamName = wfExtraParamNames{ndx};
              if isfield(wfParamsThisAmountRaw,wfType) && isfield(wfParamsThisAmountRaw.(wfType),extraParamName)
                % Use a value from the preset, if present
                wfParams.(wfType).values.(extraParamName) = wfParamsThisAmountRaw.(wfType).(extraParamName);
              else
                % If not present, fall back to the default default
                wfParams.(wfType).values.(extraParamName) = defaultWFExtraParams{ndx};
              end
            end
            % If this WF type is included in the lexicon preset, enable the WF type
            wfParams.(wfType).enabled=isfield(wfParamsThisAmountRaw,wfType);
            % If this WF type is included in the lexicon preset, check for any WF
            % parameter names besides the usual standard and extra ones,
            % and include them in the processed preset parameters.
            % This seems weird.  Maybe we can leave it out?
            % if isfield(wfParamsThisAmount,wfType)
            %   wfParamNamesThisAmount = fieldnames(wfParamsThisAmount.(wfType));
            %   for dndx = 1:numel(wfParamNamesThisAmount)
            %     if ismember(wfParamNamesThisAmount{dndx},wfParamNames) ,
            %       wfParamName = wfParamNamesThisAmount{dndx};
            %       wfParams.(wfType).values.(wfParamName) = wfParamsThisAmount.(wfType).(wfParamName);
            %     end
            %   end
            % end
          end
          wfParamsFromAmount.(wfAmount) = wfParams;
        end
      end
    end  % method

    % ---------------------------------------------------------------------
    function defaultPFTransTypesFromName= ...
        defaultPFTransTypesFromNameFromFeatureLexicon(featureLexicon,scoreFeatures)
      % Calculate the structure that maps PF names to transformation types.

      % Get the list of all the per-frame features in the lexicon,
      % including the scores-as-inputs ones
      lexiconPFNames=fieldnames(featureLexicon.perframe);
      scorePFNames = {scoreFeatures(:).scorefilename};
      subdialectPFNames=[lexiconPFNames;scorePFNames'];

      % Make one entry in the result structure per PF
      defaultPFTransTypesFromName = struct;
      for i = 1:length(subdialectPFNames)
        pfName = subdialectPFNames{i};
        if ismember(pfName,scorePFNames) ,
          transTypesThis = {'none'};
        else
          transTypesThisRaw = featureLexicon.perframe.(pfName).trans_types;
          if ischar(transTypesThisRaw)
            transTypesThis = {transTypesThisRaw};
          else
            transTypesThis = transTypesThisRaw;
          end
        end
        defaultPFTransTypesFromName.(pfName)=transTypesThis;
      end  % for
    end  % method

    % ---------------------------------------------------------------------
    function pfCategoriesFromName= ...
        pfCategoriesFromNameFromFeatureLexicon(featureLexicon,scoreFeatures)
      % Calculate the structure that maps a PF name to the PF categories
      % that it belongs to.

      % Get the list of all the per-frame features in the lexicon,
      % including the scores-as-inputs ones
      lexiconPFNames=fieldnames(featureLexicon.perframe);
      scorePFNames = {scoreFeatures(:).scorefilename};
      subdialectPFNames=[lexiconPFNames;scorePFNames'];

      % Make one entry in the result structure per PF
      pfCategoriesFromName = struct;
      for i = 1:length(subdialectPFNames)
        pfName = subdialectPFNames{i};
        if ismember(pfName,scorePFNames) ,
          pfCategoriesThis = {'scores'};
        else
          pfCategoriesThisRaw = featureLexicon.perframe.(pfName).type;
          if ischar(pfCategoriesThisRaw)
            pfCategoriesThis = {pfCategoriesThisRaw};
          else
            pfCategoriesThis = pfCategoriesThisRaw;
          end
        end
        pfCategoriesFromName.(pfName)=pfCategoriesThis;
      end  % for
    end  % method

    % ---------------------------------------------------------------------
    function pfCategoryNames=pfCategoryNamesFromFeatureLexicon(featureLexicon)
      fallpf = fieldnames(featureLexicon.perframe);
      pfCategoryNames = {};
      for pfndx = 1:numel(fallpf)
        curpf = fallpf{pfndx};
        curtypes  = featureLexicon.perframe.(curpf).type;
        if ischar(curtypes)
          curT = curtypes;
          if ~any(strcmp(pfCategoryNames,curT))
            pfCategoryNames{end+1} = curT;  %#ok
          end
        else
          for tndx = 1:numel(curtypes)
            curT = curtypes{tndx};
            if ~any(strcmp(pfCategoryNames,curT))
              pfCategoryNames{end+1} = curT;  %#ok
            end
          end
        end
      end
      if ~any(strcmp(pfCategoryNames,'scores')),
        pfCategoryNames{end+1} = 'scores';
      end
    end  % method

  end  % class methods

  % -----------------------------------------------------------------------
  % public instance methods
  methods
    % ---------------------------------------------------------------------
    function self=FeatureVocabularyForSelectFeatures(featureLexicon, ...
                                                     scoreFeatures, ...
                                                     subdialectPFNames, ...
                                                     windowFeatureParams, ...
                                                     maxWindowRadiusCommon)
      % Constructor.  jld an instance of JLabelData, maxWindowRadiusCommon
      % a value to be set for max_window_radius for all PFs, all WF types,
      % and all preset WF amounts.

      % Set a bunch of things that are derived from the constructor input,
      % but stay constant for the life of the object
      self.wfParamsFromAmount= ...
        FeatureVocabularyForSelectFeatures.computeWFParamsFromAmount(featureLexicon);
      self.defaultPFTransTypesFromName= ...
        FeatureVocabularyForSelectFeatures.defaultPFTransTypesFromNameFromFeatureLexicon(featureLexicon, ...
                                                                                scoreFeatures);
      self.pfCategoriesFromName= ...
        FeatureVocabularyForSelectFeatures.pfCategoriesFromNameFromFeatureLexicon(featureLexicon, ...
                                                                                scoreFeatures);
      self.pfCategoryNames= ...
        FeatureVocabularyForSelectFeatures.pfCategoryNamesFromFeatureLexicon(featureLexicon);

      % Populate the list of per-frame feature names
      self.subdialectPFNames = subdialectPFNames;

      % Populate the feature vocabulary proper
      enabledPFs = fieldnames(windowFeatureParams);
      wfParamNames = self.wfParamNames;
      wfAmounts=self.wfAmounts;
      wfAmount=wfAmounts{1};  % by convention, the first one is the default setting (e.g. 'normal')
      wfParamsAmount=self.wfParamsFromAmount.(wfAmount);
      nSubdialectPFNames=length(subdialectPFNames);
      vocabulary = cell(1,nSubdialectPFNames);
      for pfIndex = 1:nSubdialectPFNames
        pfName = subdialectPFNames{pfIndex};
        vocabulary{pfIndex}.name = pfName;
        if ismember(pfName,enabledPFs) ,
          % if this PF is enabled, copy its params out of
          % windowFeatureParams
          vocabulary{pfIndex}.enabled = true;
          vocabulary{pfIndex}.sanitycheck = windowFeatureParams.(pfName).sanitycheck;
          % Fill for the window feature type called 'default'
          for wfParamNdx = 1:numel(wfParamNames)
            wfParamName = wfParamNames{wfParamNdx};
            if isfield(windowFeatureParams.(pfName),wfParamName)
              vocabulary{pfIndex}.default.values.(wfParamName) = windowFeatureParams.(pfName).(wfParamName);
            else % Fill in the value from the default preset amount
              vocabulary{pfNdx}.default.values.(wfParamName) = ...
                wfParamsAmount.default.values.(wfParamName);
                %self.defaultWFParams{wfParamNdx};
            end
            vocabulary{pfIndex}.default.enabled = true;
          end
          % Fill for the rest of the window feature types.
          enabledWFTypes = fieldnames(windowFeatureParams.(pfName));
          for wfTypeNdx = 2:numel(self.wfTypes)
            wfType = self.wfTypes{wfTypeNdx};
            if ismember(wfType,enabledWFTypes) ,
              vocabulary{pfIndex}.(wfType).enabled = true;
              wfParamsThisType = windowFeatureParams.(pfName).(wfType);
              for wfParamNdx = 1:numel(wfParamNames)
                wfParamName = wfParamNames{wfParamNdx};
                if isfield(wfParamsThisType,wfParamName)
                  vocabulary{pfIndex}.(wfType).values.(wfParamName) = wfParamsThisType.(wfParamName);
                else  % fill in the default values
                  vocabulary{pfIndex}.(wfType).values.(wfParamName) = windowFeatureParams.(pfName).(wfParamName);
                end
              end
              if ~isempty(self.wfExtraParamNames{wfTypeNdx})
                extraParamName = self.wfExtraParamNames{wfTypeNdx};
                vocabulary{pfIndex}.(wfType).values.(extraParamName) = wfParamsThisType.(extraParamName);
              end
            else
              % If the window type is disabled, use default values for its
              % parameters
              vocabulary{pfIndex}.(wfType).enabled = false;
              for wfParamNdx = 1:numel(wfParamNames)
                wfParamName = wfParamNames{wfParamNdx};
                vocabulary{pfIndex}.(wfType).values.(wfParamName) = windowFeatureParams.(pfName).(wfParamName);
              end
              if ~isempty(self.wfExtraParamNames{wfTypeNdx})
                extraParamName = self.wfExtraParamNames{wfTypeNdx};
                %vocabulary{pfIndex}.(wfType).values.(extraParamName) = self.defaultWFExtraParams{wfTypeNdx};
                vocabulary{pfIndex}.(wfType).values.(extraParamName) = wfParamsAmount.(wfType).values.(extraParamName);
              end
            end
          end
        else
          % If PF is disabled, use values from the default preset amount
          vocabulary{pfIndex}=wfParamsAmount;
          vocabulary{pfIndex}.enabled = false;
          vocabulary{pfIndex}.sanitycheck = false;
          % Tack on the default translation types, which depend on which
          % per-frame feature this is
          for i=1:length(self.wfTypes)
            wfType=self.wfTypes{i};
            vocabulary{pfIndex}.(wfType).values.trans_types=self.defaultPFTransTypesFromName.(pfName);
          end
%           vocabulary{pfIndex}.default.enabled = true;
%           for wfParamNamesNdx = 1:numel(self.wfParamNames)
%             wfParamName = self.wfParamNames{wfParamNamesNdx};
%             vocabulary{pfIndex}.default.values.(wfParamName) = ...
%               self.defaultWFParams{wfParamNamesNdx};
%           end
%           % Copy the default values into the other window params.
%           for winfnNdx = 2:numel(self.wfTypes)
%             wfType = self.wfTypes{winfnNdx};
%             vocabulary{pfIndex}.(wfType).enabled = false;
%             for wfParamNamesNdx = 1:numel(self.wfParamNames)
%               wfParamName = self.wfParamNames{wfParamNamesNdx};
%               vocabulary{pfIndex}.(wfType).values.(wfParamName) = ...
%                 vocabulary{pfIndex}.default.values.(wfParamName);
%             end
%             if ~isempty(self.wfExtraParamNames{winfnNdx})
%               extraParamName = self.wfExtraParamNames{winfnNdx};
%               vocabulary{pfIndex}.(wfType).values.(extraParamName) = self.defaultWFExtraParams{winfnNdx};
%             end
%           end
        end
      end
      self.vocabulary = vocabulary;

      % If a common max window radius is set, use it
      if exist('maxWindowRadiusCommon','var') && ~isempty(maxWindowRadiusCommon)
        % set the global maxWindowRadius to the given value
        self.setMaxWindowRadiusForAllWFs(maxWindowRadiusCommon);
        self.setMaxWindowRadiusForAllWFAmounts(maxWindowRadiusCommon);
      end
    end  % constructor method

    % ---------------------------------------------------------------------
    function enableAllPFsInCategory(self,pfCategoryIndex)
      %pfCategory=self.pfCategoryNames{pfCategoryIndex};
      %subdialectPFNames=self.subdialectPFNames;
      %pfCategoriesFromName=self.pfCategoriesFromName;
      pfNamesThisCategory=self.getPFNamesInCategory(pfCategoryIndex);
      % Iterate over the per-frame features, looking for ones that are
      % within the selected category.
      for i = 1:length(pfNamesThisCategory)
        pfName=pfNamesThisCategory{i};
        self.enablePerframeFeature(pfName)
      end
    end  % method

    % ---------------------------------------------------------------------
    function setAllPFsInCategoryToWFAmount(self,pfCategoryIndex,wfAmount)
      % This does what it says, but note that it doesn't add or subtract
      % any PFs from the vocabulary.
      %pfCategory=self.pfCategoryNames{pfCategoryIndex};
      %subdialectPFNames=self.subdialectPFNames;
      %pfCategoriesFromName=self.pfCategoriesFromName;
      pfNamesThisCategory=self.getPFNamesInCategory(pfCategoryIndex);
      % Iterate over the per-frame features, looking for ones that are
      % within the selected category.
      for i = 1:length(pfNamesThisCategory)
        pfName=pfNamesThisCategory{i};
        self.setPFToWFAmount(pfName,wfAmount)
      end
    end  % method

    % ---------------------------------------------------------------------
    function setPFEnablement(self,pfSpecifier,enabled)
      % Turn on/off the given per-frame feature, but leave the window feature
      % parameters alone.
      % pfIndicator can be a name or an index.
      if ischar(pfSpecifier)
        % this means pfIndex is really a name
        pfName=pfSpecifier;
        pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      else
        pfIndex=pfSpecifier;
      end
      self.vocabulary{pfIndex}.enabled = enabled;
    end  % method

    % ---------------------------------------------------------------------
    function enablePerframeFeature(self,pfSpecifier)
      % Turn on the given per-frame feature, but leave the window feature
      % parameters alone.
      % This method is deprecated.
      setPFEnablement(self,pfSpecifier,true);
    end  % method

    % ---------------------------------------------------------------------
    function disablePerframeFeature(self,pfSpecifier)
      % Turn off the given per-frame feature.  (I.e. remove all it's window
      % features from the vocabulary)
      % This method is deprecated.
      setPFEnablement(self,pfSpecifier,false);
    end  % method

    % ---------------------------------------------------------------------
    function setPFToWFAmount(self,pfSpecifier,wfAmount)
      % Set the window feature parameters for the named PF to those
      % specified by the named wfAmount (normal, more, less).  This is
      % orthogonal to whether that PF is enabled or not.  Note that this
      % also sets the transformation types for the PF to the default
      % set.
      if isequal(wfAmount,'custom')
        return
      end
      if ischar(pfSpecifier)
        pfName=pfSpecifier;
        pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      else
        pfIndex=pfSpecifier;
        pfName=self.subdialectPFNames{pfIndex};
      end
      defaultPFTransTypes=self.defaultPFTransTypesFromName.(pfName);
      %self.vocabulary{pfIndex}.enabled = true;
      %pfCategories=self.pfCategoriesFromName.(pfName);
      %pfCategory = pfCategories{1};  % why the first one?
      %categoryNdx = find(strcmp(pfCategory,self.pfCategoryNames));
      wfParamTemplate=self.wfParamsFromAmount.(wfAmount);
      %self.vocabulary{pfIndex}.enabled = true;
      for wfTypeIndex = 1:numel(self.wfTypes)
        wfType = self.wfTypes{wfTypeIndex};
        %wfParamsThisType=wfParamTemplate.(wfType);
        %if ~wfParamsThisType.enabled
        %  self.vocabulary{pfIndex}.(wfType).enabled = false;
        %  continue;
        %end
        self.setWindowFeaturesOfSingleType(pfIndex,wfTypeIndex,wfParamTemplate);
        self.vocabulary{pfIndex}.(wfType).values.trans_types = defaultPFTransTypes;
      end
    end  % method

    % ---------------------------------------------------------------------
    function setWFTypeEnablement(self,pfName,wfType,enable)
      % For the given per-frame feature name, enable/disable the given
      % window-feature type
      pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      self.vocabulary{pfIndex}.(wfType).enabled = enable;  %#ok
    end

    % ---------------------------------------------------------------------
    function setWFParam(self,pfName,wfType,wfParamName,newValue)
      % For the given per-frame feature name, window-feature type, and
      % window-feature parameter name, set the parameter to newValue
      if ischar(pfName)
        % the usual case
        pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      else
        % this means pfName is really an index
        pfIndex=pfName;
        %pfName=self.subdialectPFNames{pfIndex};
      end
      if ~ischar(wfType)
        % this means wfType is really an index
        wfIndex=wfType;
        wfType=FeatureVocabulary.wfTypes{wfIndex};
      end
      self.vocabulary{pfIndex}.(wfType).values.(wfParamName) = newValue;
    end

    % ---------------------------------------------------------------------
    function addWFTransformation(self,pfName,wfType,newTransformation)
      % For per-frame feature pfName, window feature type wfType, add a new
      % transformation, newTransformation.  (A transformation is one of 'none',
      % 'flip', 'abs', 'relative'.)
      pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      transformations = self.vocabulary{pfIndex}.(wfType).values.trans_types;
      if ~any(strcmp(newTransformation,transformations))
        % if it's not in the list, add it to the end
        transformations{end+1}=newTransformation;
        self.vocabulary{pfIndex}.(wfType).values.trans_types=transformations;
      end
    end  % method

    % ---------------------------------------------------------------------
    function removeWFTransformation(self,pfName,wfType,transformation)
      % For per-frame feature pfName, window feature type wfType, remove a
      % transformation, transformation.  (A transformation is one of 'none',
      % 'flip', 'abs', 'relative'.)  If transformation is not in the
      % transformation list for that window feature, nothing changes.  If
      % removing the given transformation would leave the transformation
      % list empty, nothing changes.
      pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      transformations = self.vocabulary{pfIndex}.(wfType).values.trans_types;
      toBeRemoved = strcmp(transformation,transformations);
      transformations(toBeRemoved) = [];
      % only commit the change if there is still at least one transformation
      % left
      if ~isempty(transformations)
        self.vocabulary{pfIndex}.(wfType).values.trans_types=transformations;
      end
    end  % method

    % ---------------------------------------------------------------------
    function result=isWFTransformation(self,pfName,wfType,transformation)
      % For per-frame feature pfName, window feature type wfType, returns
      % true iff the given transformation is in the transformation list.
      pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      transformations = self.vocabulary{pfIndex}.(wfType).values.trans_types;  %#ok
      result=ismember(transformation,transformations);
    end  % method

    % ---------------------------------------------------------------------
    function setMaxWindowRadiusForAllWFAmounts(self,newValue)
      % This sets the maximum window radius for all of the preset window
      % feature amounts, for all window feature types.  Note, however, that
      % it does not change the window feature parameters
      % themselves.  So WF params that were previously set using one of the
      % WF amount pre-sets will retain the old values.
      wfAmounts = self.wfAmounts;
      wfTypes = self.wfTypes;
      for i = 1:numel(wfAmounts)
        wfAmount = wfAmounts{i};
        for j = 1:numel(wfTypes)
          wfType=wfTypes{j};
          wfParams=self.wfParamsFromAmount.(wfAmount);
          if isfield(wfParams,wfType),
            self.wfParamsFromAmount.(wfAmount).(wfType).values.max_window_radius = newValue;
          end
        end
      end
    end  % method

    % ---------------------------------------------------------------------
    function setMaxWindowRadiusForAllWFs(self,newValue)
      % This sets the maximum window radius for all per-frame features, for
      % all window-feature types.  Thus for all window features.
      pfNames=self.subdialectPFNames;
      wfTypes = self.wfTypes;
      for i = 1:length(pfNames)
        for j = 1:length(wfTypes)
          wfType=wfTypes{j};
          if isfield(self.vocabulary{i},wfType),
            self.vocabulary{i}.(wfType).values.max_window_radius=newValue;
          end
        end
      end
    end  % method

    % ---------------------------------------------------------------------
    function [thereIsConsensus,consensusValue]=getConsensusMaxWindowRadiusForAllWFs(self)
      % Determines the common value of max_window_radius across all WFs, in
      % all PFs.  thereIsConsensus will be true iff there is a common
      % value.  The common value is returned in consensusValue.  If there
      % is no common value, consensusValue is unspecified.  If there are no
      % WFs at all in the vocabulary, thereIsConsensus will be true, but
      % consensusValue will be empty.
      thereIsConsensus=true;
      consensusValue=[];
      pfNames=self.subdialectPFNames;
      wfTypes = self.wfTypes;
      for i = 1:length(pfNames)
        for j = 1:length(wfTypes)
          wfType=wfTypes{j};
          if isfield(self.vocabulary{i},wfType),
            thisValue=self.vocabulary{i}.(wfType).values.max_window_radius;
            if isempty(consensusValue)
              % If this is the first value we've examined, it becomes the
              % consensus value
              consensusValue=thisValue;
            else
              % if this is not the first value we've examined, compare it
              % to the consensusValue.
              if thisValue~=consensusValue
                % No consensus value---return empty matrix
                thereIsConsensus=false;
                return
              end
            end
          end
        end
      end
    end  % method

    % ---------------------------------------------------------------------
    function copyWFParams(self,pfNdxFrom,wfTypeNdxFrom,pfNdxTo,wfTypeNdxTo)
      wfTypeFrom = self.wfTypes{wfTypeNdxFrom};
      % something to copy from?
      if ~isfield(self.vocabulary{pfNdxFrom},wfTypeFrom),
        return;
      end
      wfTypeTo = self.wfTypes{wfTypeNdxTo};
      self.vocabulary{pfNdxTo}.(wfTypeTo).enabled = self.vocabulary{pfNdxFrom}.(wfTypeFrom).enabled;
      for i = 1:numel(self.wfParamNames),
        wfParamName = self.wfParamNames{i};
        self.vocabulary{pfNdxTo}.(wfTypeTo).values.(wfParamName) = ...
          self.vocabulary{pfNdxFrom}.(wfTypeFrom).values.(wfParamName);
      end
      if ~isempty(self.wfExtraParamNames{wfTypeNdxTo})
        extraParam = self.wfExtraParamNames{wfTypeNdxTo};
        if isfield(self.vocabulary{pfNdxFrom}.(wfTypeFrom).values,extraParam)
          self.vocabulary{pfNdxTo}.(wfTypeTo).values.(extraParam) = ...
            self.vocabulary{pfNdxFrom}.(wfTypeFrom).values.(extraParam);
        else
          self.vocabulary{pfNdxTo}.(wfTypeTo).values.(extraParam) = '';
        end
      end
    end  % method

    % ---------------------------------------------------------------------
    function result=get.wfAmounts(self)
      result=fieldnames(self.wfParamsFromAmount);
    end  % method

    % ---------------------------------------------------------------------
    function result=pfIsInVocabulary(self,pfIndex)
      if ischar(pfIndex)
        % means pfIndex is really a per-frame feature name
        pfName=pfIndex;
        pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      end
      result=self.vocabulary{pfIndex}.enabled;
    end  % method

    % ---------------------------------------------------------------------
    function result=wfTypeIsInVocabulary(self,pfIndex,wfType)
      if ischar(pfIndex)
        % means pfIndex is really a per-frame feature name
        pfName=pfIndex;
        pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      end
      if isnumeric(wfType)
        % means wfType is really a window feature index
        wfIndex=wfType;
        wfType=self.wfTypes{wfIndex};
      end
      result=self.vocabulary{pfIndex}.(wfType).enabled;
    end  % method

    % ---------------------------------------------------------------------
    function windowFeatureParams = getInJLabelDataFormat(self)
      % Converts the feature vocabulary into the format used by JLabelData,
      % returns this.
      % Note that none of the information about window-feature amount
      % presets (e.g. normal, more, less) gets used in computing the
      % output.
      windowFeatureParams = struct;
      for pfIndex = 1:numel(self.subdialectPFNames)
        pfName = self.subdialectPFNames{pfIndex};
        if ~self.vocabulary{pfIndex}.enabled; continue;end
        pfParams = self.vocabulary{pfIndex};
        windowFeatureParams.(pfName).sanitycheck = pfParams.sanitycheck;
        % the default window feature type is always included
        for wfParamsIndex = 1:numel(self.wfParamNames)
          wfParamName = self.wfParamNames{wfParamsIndex};
          windowFeatureParams.(pfName).(wfParamName) = pfParams.default.values.(wfParamName);
        end
        for wfTypeIndex = 2:numel(self.wfTypes)
          wfType = self.wfTypes{wfTypeIndex};
          if pfParams.(wfType).enabled,
            for wfParamsIndex = 1:numel(self.wfParamNames)
              wfParamName = self.wfParamNames{wfParamsIndex};
              windowFeatureParams.(pfName).(wfType).(wfParamName) =...
                  pfParams.(wfType).values.(wfParamName);
            end
            if ~isempty(self.wfExtraParamNames{wfTypeIndex})
              extraParam = self.wfExtraParamNames{wfTypeIndex};
              extraParamVal = pfParams.(wfType).values.(extraParam);
              windowFeatureParams.(pfName).(wfType).(extraParam) = extraParamVal;
            end
          end
        end
      end
    end  % method

    % ---------------------------------------------------------------------
    function pfNames=getPFNamesInCategory(self,pfCategoryName)
      % For the given per-frame feature category, return a cell array of
      % strings containing the names of all per-frame features in that
      % category.
      if isnumeric(pfCategoryName)
        % this means the category name is really a category index
        pfCategoryIndex=pfCategoryName;
        pfCategoryName=self.pfCategoryNames{pfCategoryIndex};
      end
      pfNames=cell(1,0);
      for iPF = 1:numel(self.subdialectPFNames)
        pfName=self.subdialectPFNames{iPF};
        categoriesThisPF=self.pfCategoriesFromName.(pfName);
        if ismember(pfCategoryName,categoriesThisPF)
          % if the selected category contains this per-frame feature,
          % do stuff
          pfNames{1,end+1}=pfName;  %#ok
        end
      end
    end  % method

    % ---------------------------------------------------------------------
    function level=getPFCategoryLevel(self,pfCategoryName)
      % For the given per-frame feature category, calculate what 'level' of
      % per-frame features are in the vocabulary.  Returns one of 'none',
      % 'custom', 'all'.
      pfNamesInCategory=self.getPFNamesInCategory(pfCategoryName);
      nPFsInCategory=length(pfNamesInCategory);
      nPFsInCategoryAndVocab=0;
      for i = 1:nPFsInCategory
        pfName=pfNamesInCategory{i};
        if self.pfIsInVocabulary(pfName),
          nPFsInCategoryAndVocab=nPFsInCategoryAndVocab+1;
        end
      end
      if nPFsInCategoryAndVocab==0
        if nPFsInCategoryAndVocab==nPFsInCategory
          level='n/a';
        else
          level='none';
        end
      elseif nPFsInCategoryAndVocab==nPFsInCategory
        level='all';
      else
        level='custom';
      end
    end  % method

    % ---------------------------------------------------------------------
    function wfAmount=getWFAmountForPFCategory(self,pfCategoryName)
      % Get the current amount for the given per-frame feature category.  Returns
      % one of 'normal', 'more', 'less', 'custom', and 'n/a' (if there are
      % no PFs in the category.
      pfNamesInCategory=self.getPFNamesInCategory(pfCategoryName);
      nPFsInCategory=length(pfNamesInCategory);
      if (nPFsInCategory==0)
        wfAmount='n/a';
        return
      end
      % Get the amount for the first PF
      pfName=pfNamesInCategory{1};
      wfAmountPutative=self.getWFAmountForPF(pfName);
      % If the first one is 'custom', no need to check the rest
      if isequal(wfAmountPutative,'custom')
        wfAmount='custom';
        return
      end
      % If we get here, wfAmountThisPutative is something other than
      % 'custom'
      % The rest have to match wfAmountPutative, or else the amount is
      % 'custom'
      for i = 2:nPFsInCategory
        pfName=pfNamesInCategory{i};
        wfAmountThis=self.getWFAmountForPF(pfName);
        if ~isequal(wfAmountThis,wfAmountPutative)
          wfAmount='custom';
          return
        end
      end
      % If we get here, they all match wfAmountPutative
      wfAmount=wfAmountPutative;
    end

    % ---------------------------------------------------------------------
    function wfAmount=getWFAmountForPF(self,pfIndex)
      % get the current amount for the given per-frame feature.  Returns
      % one of 'normal', 'more', 'less', and 'custom'.
      if ischar(pfIndex)
        % means it's really a name
        pfName=pfIndex;
        pfIndex=find(strcmp(pfName,self.subdialectPFNames));
      else
        % pfIndex is really an index
        pfName=self.subdialectPFNames{pfIndex};
      end
      wfAmount='custom';  % if no match, we default to custom
      defaultPFTransTypes=self.defaultPFTransTypesFromName.(pfName);
      pfParams=self.vocabulary{pfIndex};  % all the WF params for the given PF
      nAmounts=length(self.wfAmounts);
      for j=1:nAmounts
        wfAmountTest=self.wfAmounts{j};
        isMatch=true;
        pfParamsTemplate=self.wfParamsFromAmount.(wfAmountTest);
          % all the WF params for a single PF, from the pre-set amount
        for wfTypeIndex = 1:length(self.wfTypes)
          wfType = self.wfTypes{wfTypeIndex};
          if isfield(pfParams,wfType) && ~isfield(pfParamsTemplate,wfType),
            % If this window-feature type is present in the vocabulary, but not
            % in the template, then skip to next wfType
            continue
          end
          % check that the enablement in the vocab matches the template
          if pfParams.(wfType).enabled ~= pfParamsTemplate.(wfType).enabled ,
            isMatch=false;
            break
          end
          % For all the regular WF parameter names, check that they're
          % equal
          for i = 1:numel(self.wfParamNames),
            wfParamName = self.wfParamNames{i};
            if isequal(wfParamName,'trans_types')
              % Skip transformation types, we do those below
              continue
            end
            if pfParams.(wfType).values.(wfParamName) ~= ...
               pfParamsTemplate.(wfType).values.(wfParamName) ,
              isMatch=false;
              break
            end
          end
          if ~isMatch, break, end
          % If there is an extra WF parameter for this WF type, copy it over
          % also.
          extraWFParamName = self.wfExtraParamNames{wfTypeIndex};
          if ~isempty(extraWFParamName) && ...
             isfield(pfParamsTemplate.(wfType).values,extraWFParamName) && ...
             ~isequal(pfParams.(wfType).values.(extraWFParamName), ...
                      pfParamsTemplate.(wfType).values.(extraWFParamName)) ,
            isMatch=false;
            break
          end
          % Check the transformation types also
          if ~isequal(unique(pfParams.(wfType).values.trans_types) , ...
                      unique(defaultPFTransTypes) ) ,
            isMatch=false;
            break
          end
        end  % for wfTypeIndex = 1:length(self.wfTypes)
        if isMatch,
          wfAmount=wfAmountTest;
          break  % if match, we're done
        else
          continue  % if no match, try next wfAmount
        end
      end  % for wfAmountTest=self.wfAmounts
    end  % method

  end  % public instance methods

  % -----------------------------------------------------------------------
  % private instance methods
  methods (Access=private)
    % ---------------------------------------------------------------------
    function setWindowFeaturesOfSingleType(self,pfNdx,wfTypeNdx,wfParamTemplate)
      % For the per-frame feature with index pfNdx, sets the window features
      % with type given by index wfTypeNdx to the amount given in wfParamTemplate.

      wfType = self.wfTypes{wfTypeNdx};
      %wfParamTemplate=self.wfParamsFromAmount.(wfAmount);
      % Want to do something like:
      %   self.vocabulary{pfNdx}.(wfType)=pfParamTemplate.(wfType)  ,
      % but adapted to the particular PF, I guess ---ALT, Mar 31, 2013
      % something to copy from?
      if isfield(self.vocabulary{pfNdx},wfType) && ~isfield(wfParamTemplate,wfType),
        % If this window-feature type is present in the vocabulary, but not
        % in the template, then disable it in the vocab, and return.
        % (If this ever happens, doesn't it mean the template is broken?
        % And since we control the template, this need never happen.
        % --ALT, Mar 31, 2013)
        self.vocabulary{pfNdx}.(wfType).enabled = false;
        return
      end
      % Set the enablement in the vocab to match the template
      self.vocabulary{pfNdx}.(wfType).enabled = wfParamTemplate.(wfType).enabled;
      % For all the regular WF parameter names, copy them over
      for i = 1:numel(self.wfParamNames),
        wfParamName = self.wfParamNames{i};
        self.vocabulary{pfNdx}.(wfType).values.(wfParamName) = ...
          wfParamTemplate.(wfType).values.(wfParamName);
      end
      % If there is an extra WF parameter for this WF type, copy it over
      % also.
      extraWFParamName = self.wfExtraParamNames{wfTypeNdx};
      if ~isempty(extraWFParamName)
        if isfield(wfParamTemplate.(wfType).values,extraWFParamName) ,
          value=wfParamTemplate.(wfType).values.(extraWFParamName);
        else
          % If no value is specified, set to empty matrix
          value=[];
        end
        self.vocabulary{pfNdx}.(wfType).values.(extraWFParamName) = ...
          value;
      end
    end  % method

  end  % private methods
end  % classdef

author: kristinbranson

repo: JAABA


text:
function [ H] = Hmatrix( Ix, Iy, SizeBig, alfa )

%At each pyramid level, this function generates the Hessian matrix for the
%source image

H = zeros([2 2 size(Ix)-SizeBig]);

for i = 1+SizeBig : size(Ix,1)-SizeBig
    for j = 1+SizeBig : size(Ix,2)-SizeBig

        ix = Ix( i-SizeBig:i+SizeBig, j-SizeBig:j+SizeBig );
        iy = Iy( i-SizeBig:i+SizeBig, j-SizeBig:j+SizeBig );
        H(1,1,i,j) = alfa+sum(sum( ix.^2 ));
        H(2,2,i,j) = alfa+sum(sum( iy.^2 ));
        H(1,2,i,j) = sum(sum( ix .* iy ));
        H(2,1,i,j) = H(1,2,i,j);

    end
end

end

author: Sable

repo: mcbench-benchmarks


text:
function value = i4vec_index ( n, a, aval )

%*****************************************************************************80
%
%% I4VEC_INDEX returns the first location of a given value in an I4VEC
%
%  Discussion:
%
%    An I4VEC is a vector of integer values.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    05 November 2005
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the number of entries in the vector.
%
%    Input, integer A(N), the vector to be searched.
%
%    Input, integer AVAL, the value to be indexed.
%
%    Output, integer VALUE, the first location in A which has the
%    value AVAL, or -1 if no such index exists.
%
  for i = 1 : n
    if ( a(i) == aval )
      value = i;
      return
    end
  end

  value = -1;

  return
end

author: johannesgerer

repo: jburkardt-m


text:
% Usage  [L, inliers] = ransacfitline(XYZ, t, feedback)
%
% Arguments:
%          XYZ - 3xNpts array of xyz coordinates to fit line to.
%          t   - The distance threshold between data point and the line
%                used to decide whether a point is an inlier or not.
%          feedback - Optional flag 0 or 1 to turn on RANSAC feedback
%                     information.
%
% Returns:.
%           V - Line obtained by a simple fitting on the points that
%               are considered inliers.  The line goes through the
%               calculated mean of the inlier points, and is parallel to
%               the principal eigenvector.  The line is scaled by the
%               square root of the largest eigenvalue.
%               This line is a n*2 matrix.  The first column is the
%               beginning point, the second column is the end point of the
%               line.
%           L - The two points in the data set that were found to
%               define a line having the most number of inliers.
%               The two columns of L defining the two points.
%           inliers - The indices of the points that were considered
%                     inliers to the fitted line.
%
% See also:  RANSAC, FITPLANE, RANSACFITPLANE

% Copyright (c) 2003-2006 Peter Kovesi and Felix Duvallet (CMU)
% School of Computer Science & Software Engineering
% The University of Western Australia
% http://www.csse.uwa.edu.au/
%
% Permission is hereby granted, free of charge, to any person obtaining a copy
% of this software and associated documentation files (the "Software"), to deal
% in the Software without restriction, subject to the following conditions:
%
% The above copyright notice and this permission notice shall be included in
% all copies or substantial portions of the Software.
%
% The Software is provided "as is", without warranty of any kind.

% Aug  2006 - created ransacfitline from ransacfitplane
%             author: Felix Duvallet

function [V, L, inliers] = ransacfitparabola(XY, t, feedback)

    if nargin == 2
	feedback = 0;
    end

    [rows, npts] = size(XY);

    if rows ~= 2
        error('data is not 2D');
    end

    if npts < 2
        error('too few points to fit line');
    end

    s = 3;  % Minimum No of points needed to fit a line.

    fittingfn = @defineParabola;
    distfn    = @paraptdist;
    degenfn   = @isdegenerate;

    [L, inliers] = ransac(XY, fittingfn, distfn, degenfn, s, t, feedback);

    V = fitparabola(XY(:, inliers));

%------------------------------------------------------------------------
% Function to define a parabola given 3 data points as required by
% RANSAC.

function L = defineParabola(X)
    L = X;
%------------------------------------------------------------------------
% Function to calculate distances between a line and an array of points.
function [inliers, L] = paraptdist(L, X, t)
    p = polyfit( L(1, :), L(2, :), 2 );
    YEst = polyval(p, X(1, :) );
    d = YEst - X(2, :);
    inliers = find(abs(d) < t);

%------------------------------------------------------------------------
% Function to determine whether a set of 2 points are in a degenerate
% configuration for fitting a line as required by RANSAC.
% In this case two points are degenerate if they are the same point
% or if they are exceedingly close together.

function r = isdegenerate(X)
    %find the norm of the difference of the two points
    % this will be 0 iff the two points are the same (the norm of their
    % difference is zero)
    r0 = norm(X(:,1) - X(:,2)) < eps;
    r1 = norm(X(:,1) - X(:,3)) < eps;
    r = r0 | r1;

author: DrGabor

repo: LiDAR

save_path: github-repos/MATLAB/DrGabor-LiDAR


text:
% BIOSIG2EEGLAB - convert BIOSIG structure to EEGLAB structure
%
% Usage:
%   >> OUTEEG = pop_biosig2eeglab(hdr, data, interval);
%
% Inputs:
%   hdr   - BIOSIG header
%   data  - BIOSIG data array
%
% Optional input:
%   interval - BIOSIG does not remove event which are outside of
%              the data range when importing data range subsets. This
%              parameter helps fix this problem.
%
% Outputs:
%   OUTEEG   - EEGLAB data structure
%
% Author: Arnaud Delorme, SCCN, INC, UCSD, Oct. 29, 2009-
%
% Note: BIOSIG toolbox must be installed. Download BIOSIG at
%       http://biosig.sourceforge.net
%       Contact a.schloegl@ieee.org for troubleshooting using BIOSIG.

% Copyright (C) 2003 Arnaud Delorme, SCCN, INC, UCSD, arno@salk.edu
%
% This file is part of EEGLAB, see http://www.eeglab.org
% for the documentation and details.
%
% Redistribution and use in source and binary forms, with or without
% modification, are permitted provided that the following conditions are met:
%
% 1. Redistributions of source code must retain the above copyright notice,
% this list of conditions and the following disclaimer.
%
% 2. Redistributions in binary form must reproduce the above copyright notice,
% this list of conditions and the following disclaimer in the documentation
% and/or other materials provided with the distribution.
%
% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
% AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
% IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
% ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
% LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
% CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
% SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
% INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
% CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
% ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF
% THE POSSIBILITY OF SUCH DAMAGE.

function EEG = biosig2eeglab(dat, DAT, interval, channels, importevent, importannot)

if nargin < 2
    help biosig2eeglab;
    return;
end
if nargin < 3
    interval = [];
end
if nargin < 4
    channels = [];
end
if nargin < 5
    importevent = 0;
end
if nargin < 6
    importannot = 0;
end

% import data
% -----------
if abs(dat.NS - size(DAT,1)) > 1
    DAT = DAT';
end
EEG = eeg_emptyset;

% convert to seconds for sread
% ----------------------------
if ~isfield(dat, 'InChanSelect') && max(dat.InChanSelect) > size(DAT,1)
    dat.InChanSelect = [1:size(DAT,1)];
end
EEG.nbchan          = length(dat.InChanSelect); %= size(DAT,1);
EEG.srate           = dat.SampleRate(1);
EEG.data            = DAT(dat.InChanSelect,:);  %DAT
clear DAT;
% try  % why would you do the following???????  JO
%     EEG.data            = EEG.data';
% catch,
%     pack;
%     EEG.data            = EEG.data';
% end;
EEG.setname 		= sprintf('%s file', dat.TYPE);
EEG.comments        = [ 'Original file: ' dat.FileName ];
EEG.xmin            = 0;
nepoch              = dat.NRec;
EEG.trials   = nepoch;
EEG.pnts     = size(EEG.data,2)/nepoch;

if isfield(dat,'T0')
    EEG.etc.T0 = dat.T0; % added sjo
end

if isfield(dat, 'Label') && ~isempty(dat.Label)
    if ischar(dat.Label)
        EEG.chanlocs = struct('labels', cellstr(char(dat.Label(dat.InChanSelect)))); % 5/8/2104 insert (dat.InChanSelect) Ramon
    else
        % EEG.chanlocs = struct('labels', dat.Label(1:min(length(dat.Label), size(EEG.data,1))));
        EEG.chanlocs = struct('labels', dat.Label(dat.InChanSelect)); % sjo added 120907 to avoid error below % 5/8/2104 insert (dat.InChanSelect) Ramon
    end
    if length(EEG.chanlocs) > EEG.nbchan, EEG.chanlocs = EEG.chanlocs(1:EEG.nbchan); end

    % check BIOSEMI channel types
    if strcmpi(dat.TYPE, 'BDF')
        bdftype = { 'EXG1'  'EXG2'  'EXG3'  'EXG4'  'EXG5'  'EXG6'  'EXG7'  'EXG8'  'GSR1'  'GSR2'  'Erg1'  'Erg2'  'Resp'  'Plet'  'Temp';
                    'MISC'  'MISC'  'MISC'  'MISC'  'MISC'  'MISC'  'MISC'  'MISC'  'GSR'   'GSR'   'MISC'  'MISC'  'RESP'  'PPG'   'TEMP' };
        for iChan = 1:length(EEG.chanlocs)
            matchChan = strmatch(EEG.chanlocs(iChan).labels, bdftype(1,:), 'exact');
            if ~isempty(matchChan)
                EEG.chanlocs(iChan).type = bdftype{2,matchChan};
            else
                EEG.chanlocs(iChan).type = 'EEG';
            end
        end
    end
end
EEG = eeg_checkset(EEG);

% extract events % this part I totally revamped to work...  JO
% --------------
EEG.event = [];

% startval = mode(EEG.data(end,:)); % my code
% for p = 2:size(EEG.data,2)-1
%     [codeout] = code(EEG.data(end,p));
%     if EEG.data(end,p) > EEG.data(end,p-1) & EEG.data(end,p) >= EEG.data(end,p+1)
%         EEG.event(end+1).latency =  p;
%         EEG.event(end).type = bitand(double(EEG.data(end,p)-startval),255);
%     end
% end

% lastout = mod(EEG.data(end,1),256);newevs = []; % andrey's code 8 bits
% codeout = mod(EEG.data(end,2),256);
% for p = 2:size(EEG.data,2)-1
%     nextcode = mod(EEG.data(end,p+1),256);
%     if codeout > lastout & codeout >= nextcode
%         newevs = [newevs codeout];
%         EEG.event(end+1).latency =  p;
%         EEG.event(end).type = codeout;
%     end
%     lastout = codeout;
%     codeout = nextcode;
% end

%lastout = mod(EEG.data(end,1),256*256);newevs = []; % andrey's code 16 bits
%codeout = mod(EEG.data(end,2),256*256);
%for p = 2:size(EEG.data,2)-1
%    nextcode = mod(EEG.data(end,p+1),256*256);
%    if (codeout > lastout) & (codeout >= nextcode)
%        newevs = [newevs codeout];
%        EEG.event(end+1).latency =  p;
%        EEG.event(end).type = codeout;
%    end
%    lastout = codeout;
%    codeout = nextcode;
%end

% if strcmp(dat.TYPE,'EDF') % sjo added 120907
%     disp('filetype EDF does not support events');
%     importevent = 0;
% end

if importevent
    if isfield(dat, 'BDF')
        if dat.BDF.Status.Channel <= size(EEG.data,1)
            EEG.data(dat.BDF.Status.Channel,:) = [];
        end
        EEG.nbchan = size(EEG.data,1);
        if ~isempty(EEG.chanlocs) && dat.BDF.Status.Channel <= length(EEG.chanlocs)
            EEG.chanlocs(dat.BDF.Status.Channel,:) = [];
        end
    elseif isempty(dat.EVENT.POS)
        disp('Extracting events from last EEG channel...');
        %Modifieded by Andrey (Aug.5,2008) to detect all non-zero codes:
        if length(unique(EEG.data(end, 1:100))) > 20
            disp('Warning: event extraction failure, the last channel contains data');
        elseif length(unique(EEG.data(end, :))) > 1000
            disp('Warning: event extraction failure, the last channel contains data');
        else
            thiscode = 0;
            for p = 1:size(EEG.data,2)*size(EEG.data,3)-1
                prevcode = thiscode;
                thiscode = mod(EEG.data(end,p),256*256);   % andrey's code - 16 bits
                if (thiscode ~= 0) && (thiscode~=prevcode)
                    EEG.event(end+1).latency =  p;
                    EEG.event(end).type = thiscode;
                end
            end
            EEG.data(end,:) = [];
            EEG.chanlocs(end) = [];
        end
        % recreate the epoch field if necessary
        % -------------------------------------
        if EEG.trials > 1
            for i = 1:length(EEG.event)
                EEG.event(i).epoch = ceil(EEG.event(i).latency/EEG.pnts);
            end
        end
        EEG = eeg_checkset(EEG, 'eventconsistency');
    end

    if ~isempty(dat.EVENT.POS)
        if isfield(dat, 'out') % Alois fix for event interval does not work
            if isfield(dat.out, 'EVENT')
                dat.EVENT = dat.out.EVENT;
            end
        end
        EEG.event = biosig2eeglabevent(dat.EVENT, interval, importannot); % Toby's fix

        % recreate the epoch field if necessary
        % -------------------------------------
        if EEG.trials > 1
            for i = 1:length(EEG.event)
                EEG.event(i).epoch = ceil(EEG.event(i).latency/EEG.pnts);
            end
        end

        EEG = eeg_checkset(EEG, 'eventconsistency');
    elseif isempty(EEG.event)
        disp('Warning: no event found. Events might be embedded in a data channel.');
        disp('         To extract events, use menu File > Import Event Info > From data channel');
    end
end

% convert data to single if necessary
% -----------------------------------
EEG = eeg_checkset(EEG,'makeur');   % Make EEG.urevent field
EEG = eeg_checkset(EEG);

author: sccn

repo: eeglab

save_path: github-repos/MATLAB/sccn-eeglab