text:
classdef nnconvt < nntest
  properties (TestParameter)
    depth = {1 2 3}
    numImages = {1 2 3 4}
    numFilters = {1 2 3}
    upx = {1 2 3}
    upy = {1 2 3}
    padx1 = {1 2 3}
    padx2 = {1 2 3}
    pady1 = {1 2 3}
    pady2 = {1 2 3}
    up = {1 2}
    fsx = {1 2}
    crop = {1 2 3 4 5 6 7 8}
    numGroups = {1 2 3}
  end

  methods (Test)
    function basic(test, depth, numImages, numFilters)
      m = depth ;
      n = numImages ;
      k = numFilters;
      x = test.randn(10,12,m,n) ;
      f = test.randn(3,4,k,m) ;
      b = test.randn(1,k) ;
      y = vl_nnconvt(x,f,b) ;
      dzdy = test.randn(size(y)) ;
      [dzdx,dzdf,dzdb] = vl_nnconvt(x,f,b,dzdy) ;
      test.der(@(x) vl_nnconvt(x,f,b), x, dzdy, dzdx, test.range * 1e-2) ;
      test.der(@(f) vl_nnconvt(x,f,b), f, dzdy, dzdf, test.range * 1e-2) ;
      test.der(@(b) vl_nnconvt(x,f,b), b, dzdy, dzdb, test.range) ;
    end

    function upsample_crop(test,upx,upy,padx1,pady1,padx2,pady2)
      m = 3 ; n = 2 ; k = 3;
      opts = {'upsample',[upy upx],'crop',[pady1 pady2 padx1 padx2]} ;
      x = test.randn(5,6,m,n) ;
      f = test.randn(3,4,k,m) ;
      b = test.randn(1,k) ;
      y = vl_nnconvt(x,f,b,opts{:}) ;
      dzdy = test.randn(size(y)) ;
      [dzdx,dzdf,dzdb] = vl_nnconvt(x,f,b,dzdy,opts{:}) ;
      test.der(@(x) vl_nnconvt(x,f,b,opts{:}), x, dzdy, dzdx, test.range * 1e-2) ;
      test.der(@(f) vl_nnconvt(x,f,b,opts{:}), f, dzdy, dzdf, test.range * 1e-2) ;
      test.der(@(b) vl_nnconvt(x,f,b,opts{:}), b, dzdy, dzdb, test.range) ;
    end

    function grouped_filters(test, numGroups, depth, numFilters)
      ng = numGroups ;
      m = depth ;
      k = numFilters ;
      n = 3 ;
      opts = {'numgroups',ng} ;
      x = test.randn(10,12,m*ng,n) ;
      f = test.randn(3,4,k,m*ng) ;
      b = test.randn(1,k*ng) ;
      y = vl_nnconvt(x,f,b,opts{:}) ;
      dzdy = test.randn(size(y)) ;
      [dzdx,dzdf,dzdb] = vl_nnconvt(x,f,b,dzdy,opts{:}) ;
      test.der(@(x) vl_nnconvt(x,f,b,opts{:}), x, dzdy, dzdx, test.range * 1e-2) ;
      test.der(@(f) vl_nnconvt(x,f,b,opts{:}), f, dzdy, dzdf, test.range * 1e-2) ;
      test.der(@(b) vl_nnconvt(x,f,b,opts{:}), b, dzdy, dzdb, test.range) ;
    end

    function one_one_image(test,up,fsx,crop)
      fsx = fsx*up ;
      if crop > fsx-1, return ; end
      m = 3 ;
      n = 4 ;
      k = 3 ;
      fsy = fsx * 3 ;
      x = test.randn(1,1,m,n) ;
      f = test.randn(fsy,fsx,k,m) ;
      b = test.randn(1,k) ;
      croph = floor(crop/2) ;
      opts = {'crop', [croph, crop-croph, croph, crop-croph], 'upsample', [up up]} ;
      y = vl_nnconvt(x,f,b,opts{:}) ;
      dzdy = test.randn(size(y)) ;
      [dzdx,dzdf,dzdb] = vl_nnconvt(x,f,b,dzdy,opts{:}) ;
      test.der(@(x) vl_nnconvt(x,f,b,opts{:}), x, dzdy, dzdx, test.range * 1e-2) ;
      test.der(@(f) vl_nnconvt(x,f,b,opts{:}), f, dzdy, dzdf, test.range * 1e-2) ;
      test.der(@(b) vl_nnconvt(x,f,b,opts{:}), b, dzdy, dzdb, test.range * 1e-1) ;
    end

    function test_gpu_correctnes(test)
      if ~strcmp(test.currentDevice, 'gpu'), return ; end
      opts = {...
        {'crop', [0 0 0 0], 'upsample', [1 1]}, ...
        {'crop', [5 5 8 8], 'upsample', [1 1]}, ...
        {'crop', [5 5 8 8], 'upsample', [3 2]}} ;

      variants = {{'nocudnn'}, ...
                  {'cudnn', 'cudnnworkspacelimit', 0}, ...
                  {'cudnn', 'cudnnworkspacelimit', +inf}} ;

      fh = 11 ;
      fw = 11 ;
      fn = 10 ;
      n = 32 ;
      depth = 32 ;
      x = test.randn(32,32,depth,n) ;
      w = test.randn(fh,fw,fn,depth) ;
      b = test.randn(1,fn) ;

      for o = 1:numel(opts)
        for v = 1:numel(variants)
          %args = horzcat(variants{v}, opts{o}, {'verbose'}) ;
          args = horzcat(variants{v}, opts{o}) ;
          y = vl_nnconvt(x,w,b,args{:}) ;
          dzdy = test.randn(size(y)) ;
          [dzdx,dzdw,dzdb] = vl_nnconvt(x,w,b,dzdy,args{:}) ;

          dzdy_ = gather(dzdy) ;
          y_ = vl_nnconvt(gather(x), gather(w), gather(b), opts{o}{:}) ;
          [dzdx_,dzdw_,dzdb_] = vl_nnconvt(gather(x),gather(w),gather(b), ...
                                           gather(dzdy), opts{o}{:}) ;

          test.eq(y, y_) ;
          test.eq(dzdx, dzdx_) ;
          test.eq(dzdw, dzdw_) ;
          test.eq(dzdb, dzdb_) ;
        end
      end
    end

  end
end

author: willard-yuan

repo: cnn-for-image-retrieval


text:
function [y, dzdg, dzdb, moments] = vl_nnbnorm_wrapper(x, g, b, moments, test, varargin)
%VL_NNBNORM_WRAPPER AutoNN wrapper for MatConvNet's vl_nnbnorm
%   VL_NNBNORM has a non-standard interface (returns a derivative for the
%   moments, even though they are not an input), so we must wrap it.
%   Layer.vl_nnbnorm replaces a standard VL_NNBNORM call with this one.
%
%   This also lets us supports nice features like setting the parameter
%   sizes automatically (e.g. building a net with VL_NNBNORM(X) is valid).

% Copyright (C) 2016-2017 Joao F. Henriques.
% All rights reserved.
%
% This file is part of the VLFeat library and is made available under
% the terms of the BSD license (see the COPYING file).

  % use number of channels in X to extend scalar (i.e. default) params to
  % the correct size. this way the layer can be constructed without
  % knowledge of the number of channels. scalars also permit gradient
  % accumulation with any tensor shape (in CNN_TRAIN_AUTONN).
  if isscalar(g)
    g(1:size(x,3),1) = g ;
  end
  if isscalar(b)
    b(1:size(x,3),1) = b ;
  end
  if isscalar(moments)
    moments(1:size(x,3),1:2) = moments ;
  end

  if isempty(x)
    % vl_nnbnorm throws an error on empty inputs. however, these can happen
    % in some dynamic architectures, so handle this case.
    y = zeros(size(x), 'like', x) ;
    if nargout > 1
      dzdg = zeros(size(g), 'like', g) ;
      dzdb = zeros(size(b), 'like', b) ;
      moments = zeros(size(moments), 'like', moments) ;
    end
    return
  end

  if ~test
    if numel(varargin) >= 1 && isnumeric(varargin{1})
      % training backward mode
      [y, dzdg, dzdb, moments] = vl_nnbnorm(x, g, b, varargin{:}) ;
    else
      % training forward mode
      y = vl_nnbnorm(x, g, b, varargin{:}) ;
    end
  else
    % test mode, pass in moments
    if numel(varargin) >= 1 && isnumeric(varargin{1})
      % test backward mode (to implement e.g. bnorm frozen during training)
      [y, dzdg, dzdb, moments] = vl_nnbnorm(x, g, b, varargin{1}, 'moments', moments, varargin{2:end}) ;
    else
      % test forward mode
      y = vl_nnbnorm(x, g, b, 'moments', moments, varargin{:}) ;
    end
  end

end

repo: MFT


text:
function fails = failedtests(results)

    tr = table(results);

    f = string( tr(tr.Failed,:).Name );

    f = arrayfun( @(x) extractBefore(x, "/"), f);

    fails = unique(f);
end

author: petercorke

repo: robotics-toolbox-matlab


text:
function step = stepsize_alg(iter, options)
% stepsize control algorithm.
%
% Inputs:
%       iter        number of iterations
%       options     options
% Output:
%       step        stepsize
%
% This file is part of SGDLibrary.
%
% Created by H.Kasai on Sep. 25, 2017
% Modified by H.Kasai on Sep. 28, 2017

    % extract options
    if ~isfield(options, 'step_init')
        step_init = 0.1;
    else
        step_init = options.step_init;
    end

    if ~isfield(options, 'step_alg')
        step_alg = 'fix';
    else
        if strcmp(options.step_alg, 'decay')
            step_alg = 'decay';
        elseif strcmp(options.step_alg, 'decay-2')
            step_alg = 'decay-2';
        elseif strcmp(options.step_alg, 'decay-3')
            step_alg = 'decay-3';
        elseif strcmp(options.step_alg, 'fix')
            step_alg = 'fix';
        else
            step_alg = 'decay';
        end
    end

    if ~isfield(options, 'lambda')
        lambda = 0.1;
    else
        lambda = options.lambda;
    end

    % update step-size
    if strcmp(step_alg, 'fix')
        step = step_init;
    elseif strcmp(step_alg, 'decay')
        step = step_init / (1 + step_init * lambda * iter);
    elseif strcmp(step_alg, 'decay-2')
        step = step_init / (1 + iter);
    elseif strcmp(step_alg, 'decay-3')
        step = step_init / (lambda + iter);
    end

end

author: hiroyuki-kasai

repo: SGDLibrary


text:
function exesuff=getexeext()
%
% exesuff=getexeext()
%
% get meshing external tool extension names for the current platform
%
% author: Qianqian Fang (fangq<at> nmr.mgh.harvard.edu)
%
% output:
%     exesuff: file extension for iso2mesh tool binaries
%
% -- this function is part of iso2mesh toolbox (http://iso2mesh.sf.net)
%

exesuff='.exe';
if(isunix)
	exesuff=['.',mexext];
end
if(isoctavemesh)
   if(~ispc)
      if(~ismac)
	   if(isempty(regexp(computer,'86_64')))
	      exesuff='.mexglx';
	   else
              exesuff='.mexa64';
	   end
      else
           if(isempty(regexp(computer,'86_64')))
              exesuff='.mexmaci';
           else
              exesuff='.mexmaci64';
           end
      end
   else
      exesuff='.exe';
   end
end

repo: REST


text:
function aedat = ImportAedatDataVersion1or2(aedat)
%{
Code contributions from Gemma Taverni.

This is a sub-function of importAedat - it process the data where the aedat
file format is determined to be 1 or 2.
The .aedat file format is documented here:
http://inilabs.com/support/software/fileformat/
This function is based on a combination of the loadaerdat function and
sensor-specific address interpretation functions.
There is a single input "info", a structure with the following
fields:
	- beginningOfDataPointer - Points to the byte before the beginning of the
		data
	- fileHandle - handle of the open aedat file in question
	- fileFormat - indicates the version of the aedat file format - should
		be either 1 or 2
	- startTime (optional) - if provided, any data with a timeStamp lower
		than this will not be returned.
	- endTime (optional) - if provided, any data with a timeStamp higher than
		this time will not be returned.
	- startEvent (optional) Any events with a lower count that this will not be returned.
		APS samples, if present, are counted as events.
	- endEvent (optional) Any events with a higher count that this will not be returned.
		APS samples, if present, are counted as events.
	- source - a string containing the name of the chip class from
		which the data came. Options are:
		- dvs128
		- davis240a
		- davis240b
		- davis240c
		- davis128mono
		- davis128rgb
		- davis208rgbw
		- davis208mono
		- davis346rgb
		- davis346mono
		- davis346bsi
		- davis640rgb
		- davis640mono
		- hdavis640mono
		- hdavis640rgbw
		- das1
		'file' and 'network' origins are not acceptable - the chip class is
		necessary for the interpretation of the addresses in address-events.
	- dataTypes (optional) cellarray. If present, only data types specified
		in this cell array are returned. Options are:
		special; polarity; frame; imu6; imu9; sample; ear; config.
	- subtractResetRead - optional, only for aedat2 and DAVIS -

The output is a structure with 2 fields.
	- info - the input structure, embelished with other data
	- data a structure which contains one structure for each type of
		data present. These structures are named according to the type of
		data; the options are:
		- special
		- polarity
		- frame
		- imu6
		- sample
		- ear
		Other data types supported in aedat3.0 are not implemented
		because no chip class currently implements them.
		Within each of these structures, there are typically a set of column
			vectors containing timeStamp, a valid bit and then other data fields,
			where each vector has the same number of elements.
			There are some exceptionss to this.
			In detail the contents of these structures are:
		- special
			- valid (colvector bool)
			- timeStamp (colvector uint32)
			- address (colvector uint32)
		- polarity
			- valid (colvector bool)
			- timeStamp (colvector uint32)
			- x (colvector uint16)
			- y (colvector uint16)
			- polarity (colvector bool)
		- frame
			- valid (bool)
			- frame timeStamp start ???
			- frame timeStamp end ???
			- timeStampExposureStart (uint32)
			- timeStampExposureEnd (uint32)
			- samples (matrix of uint16 r*c, where r is the number of rows and c is
				the number of columns.)
			- xStart (only present if the frame doesn't start from x=0)
			- yStart (only present if the frame doesn't start from y=0)
			- roiId (only present if this frame has an ROI identifier)
			- colChannelId (optional, if its not present, assume a mono array)
		- imu6
			- valid (colvector bool)
			- timeStamp (colvector uint32)
			- accelX (colvector single)
			- accelY (colvector single)
			- accelZ (colvector single)
			- gyroX (colvector single)
			- gyroY (colvector single)
			- gyroZ (colvector single)
			- temperature (colvector single)
		- sample
			- valid (colvector bool)
			- timeStamp (colvector uint32)
			- sampleType (colvector uint8)
			- sample (colvector uint32)
		- ear
			- valid (colvector bool)
			- timeStamp (colvector uint32)
			- position (colvector uint8)
			- channel (colvector uint16)
			- neuron (colvector uint8)
			- filter (colvector uint8)

Implementation: There is an efficient implementation of startEvent and
EndEvent, since the correct file locations to read from can be determined
in advance.
However, the implementation of startTime and endTime is not efficient, since
the file is read and then the timestamps are processed.
It is not possible to do better than this, since a binary search through
the file to find the correct file locations in advance could fail due to
non-monotonic timestamps.
%}

dbstop if error

info = aedat.info;
importParams = aedat.importParams;

% The fileFormat dictates whether there are 6 or 8 bytes per event.
if info.fileFormat == 1
	numBytesPerAddress = 2;
	numBytesPerEvent = 6;
	addrPrecision = 'uint16';
else
	numBytesPerAddress = 4;
	numBytesPerEvent = 8;
	addrPrecision = 'uint32';
end

fileHandle = aedat.importParams.fileHandle;

% Go to the EOF to find out how long it is
fseek(fileHandle, 0, 'eof');

% Calculate the number of events
info.numEventsInFile = floor((ftell(fileHandle) - info.beginningOfDataPointer) / numBytesPerEvent);

% Check the startEvent and endEvent parameters
if isfield(importParams, 'startEvent')
    startEvent = importParams.startEvent;
else
    startEvent = 1;
end
if startEvent > info.numEventsInFile
	error([	'The file contains ' num2str(info.numEventsInFile) ...
			'; the startEvent parameter is ' num2str(startEvent) ]);
end
if isfield(importParams, 'endEvent')
	endEvent = importParams.endEvent;
else
	endEvent = info.numEventsInFile;
end

if endEvent > info.numEventsInFile
	disp([	'The file contains ' num2str(info.numEventsInFile) ...
			'; the endEvent parameter is ' num2str(endEvents) ...
			'; reducing the endEvent parameter accordingly.']);
		endEvent = info.numEventsInFile;
end
if startEvent >= endEvent
	error([	'The startEvent parameter is ' num2str(startEvent) ...
		', but the endEvent parameter is ' num2str(endEvent) ]);
end

if isfield(importParams, 'startPacket')
	error('The startPacket parameter is set, but range by packets is not available for .aedat version < 3 files')
end
if isfield(importParams, 'endPacket')
	error('The endPacket parameter is set, but range by events is not available for .aedat version < 3 files')
end

numEventsToRead = endEvent - startEvent + 1;

% Read addresses
disp('Reading addresses ...')
fseek(fileHandle, info.beginningOfDataPointer + numBytesPerEvent * startEvent, 'bof');
allAddr = uint32(fread(fileHandle, numEventsToRead, addrPrecision, 4, 'b'));

% Read timestamps
disp('Reading timestamps ...')
fseek(fileHandle, info.beginningOfDataPointer + numBytesPerEvent * startEvent + numBytesPerAddress, 'bof');
allTs = uint32(fread(fileHandle, numEventsToRead, addrPrecision, numBytesPerAddress, 'b'));

% Trim events outside time window
% This is an inefficent implementation, which allows for
% non-monotonic timestamps.

if isfield(importParams, 'startTime')
    disp('Trimming to start time ...')
	tempIndex = allTs >= startTime * 1e6;
	allAddr = allAddr(tempIndex);
	allTs	= allTs(tempIndex);
end

if isfield(importParams, 'endTime')
    disp('Trimming to end time ...')
	tempIndex = allTs <= endTime * 1e6;
	allAddr = allAddr(tempIndex);
	allTs	= allTs(tempIndex);
end

% Interpret the addresses
%{
- Split between DVS/DAVIS and DAS.
	For DAS1:
		- Special events - external injected events has never been
		implemented for DAS
		- Split between Address events and ADC samples
		- Intepret address events
		- Interpret ADC samples
	For DVS128:
		- Special events - external injected events are on bit 15 = 1;
		there is a more general label for special events which is bit 31 =
		1, but this has ambiguous interpretations; it is also overloaded
		for the stereo pair encoding - ignore this.
		- Intepret address events
	For DAVIS:
		- Special events
			- Interpret IMU events from special events
		- Interpret DVS events according to chip class
		- Interpret APS events according to chip class
%}

% Declare function for finding specific event types in eventTypes cell array
cellFind = @(string)(@(cellContents)(strcmpi(string, cellContents)));

% Create structure to put all the data in
data = struct;

%% DAS1

if strcmp(info.source, 'Das1')
	% DAS1
	sampleMask = hex2dec('1000');
	sampleLogical = bitand(allAddr, sampleMask);
	earLogical = ~sampleLogical;
	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('ear'), importParams.dataTypes))) && any(sampleLogical)
		% ADC Samples
		data.sample.timeStamp = allTs(sampleLogical);
		% Sample type
		sampleTypeMask = hex2dec('1c00'); % take ADC scanner sync and ADC channel together for this value - kludge - the alternative would be to introduce a special event type to represent the scanner wrapping around
		sampleTypeShiftBits = 10;
		data.sample.sampleType = uint8(bitshift(bitand(allAddr(sampleLogical), sampleTypeMask), -sampleTypeShiftBits));
		% Sample data
		sampleDataMask = hex2dec('3FF'); % take ADC scanner sync and ADC channel together for this value - kludge - the alternative would be to introduce a special event type to represent the scanner wrapping around
		data.sample.sample = uint32(bitand(allAddr(sampleLogical), sampleTypeMask));
	end
	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('ear'), importParams.dataTypes))) && any(earLogical)
		% EAR events
		data.ear.timeStamp = allTs(earLogical);
		% Filter (0 = BPF, 1 = SOS)
		filterMask     = hex2dec('0001');
		data.ear.filter = uint8(bitand(allAddr, filterMask));
		% Position (0 = left; 1 = right)
		positionMask   = hex2dec('0002');
		positionShiftBits = 1;
		data.ear.position = uint8(bitshift(bitand(allAddr, positionMask), -positionShiftBits));
		% Channel (0 (high freq) to 63 (low freq))
		channelMask = hex2dec('00FC');
		channelShiftBits = 2;
		data.ear.channel = uint16(bitshift(bitand(allAddr, channelMask), -channelShiftBits));
		% Neuron (in the range 0-3)
		neuronMask  = hex2dec('0300');
		neuronShiftBits = 8;
		data.ear.neuron = uint8(bitshift(bitand(allAddr, neuronMask), -neuronShiftBits));
    end

%% DVS128

elseif strcmp(info.source, 'Dvs128')
	% DVS128
	specialMask = hex2dec ('8000');
	specialLogical = boolean(bitand(allAddr, specialMask));
	polarityLogical = ~specialLogical;
	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('special'), importParams.dataTypes))) && any(specialLogical)
		% Special events
		data.special.timeStamp = allTs(specialLogical);
		% No need to create address field, since there is only one type of special event
	end
	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('polarity'), importParams.dataTypes))) && any(polarityLogical)
		% Polarity events
		data.polarity.timeStamp = allTs(polarityLogical); % Use the negation of the special mask for polarity events
		% Y addresses
		yMask = hex2dec('7F00');
		yShiftBits = 8;
		data.polarity.y = uint16(bitshift(bitand(allAddr(polarityLogical), yMask), -yShiftBits));
		% X addresses
		xMask = hex2dec('fE');
		xShiftBits = 1;
		data.polarity.x = uint16(bitshift(bitand(allAddr(polarityLogical), xMask), -xShiftBits));
		% Polarity bit
		polBit = 1;
		data.polarity.polarity = bitget(allAddr(polarityLogical), polBit) == 0;
    end

%% Davis

elseif (~isempty(strfind(info.source, 'Davis')) ...
        && strfind(info.source, 'Davis')) ...
       || (~isempty(strfind(info.source, 'SecDvs')) ...
           &&strfind(info.source, 'SecDvs'))
	% DAVIS
	% In the 32-bit address:
	% bit 32 (1-based) being 1 indicates an APS sample
	% bit 11 (1-based) being 1 indicates a special event
	% bits 11 and 32 (1-based) both being zero signals a polarity event
    disp('Constructing logical indices for event types ...')
	apsOrImuMask = hex2dec ('80000000');
	apsOrImuLogical = bitand(allAddr, apsOrImuMask);
	ImuOrPolarityMask = hex2dec ('800');
	ImuOrPolarityLogical = bitand(allAddr, ImuOrPolarityMask);
	signalOrSpecialMask = hex2dec ('400');
	signalOrSpecialLogical = bitand(allAddr, signalOrSpecialMask);
	frameLogical = apsOrImuLogical & ~ImuOrPolarityLogical;
	imuLogical = apsOrImuLogical & ImuOrPolarityLogical;
	polarityLogical = ~apsOrImuLogical & ~signalOrSpecialLogical;
	specialLogical = ~apsOrImuLogical & signalOrSpecialLogical;

	% These masks are used for both frames and polarity events, so are
	% defined outside of the following if statement
	yMask = hex2dec('7FC00000');
	yShiftBits = 22;
	xMask = hex2dec('003FF000');
	xShiftBits = 12;

%% Davis special events

	% Special events
	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('special'), importParams.dataTypes))) && any(specialLogical)
		disp('Importing special events ...')
        data.special.timeStamp = allTs(specialLogical);
		% No need to create address field, since there is only one type of special event
    end

%% Davis polarity events

	% Polarity (DVS) events
	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('polarity'), importParams.dataTypes))) && any(polarityLogical)
		disp('Importing polarity events ...')
		data.polarity.timeStamp = allTs(polarityLogical);
		% Y addresses
		data.polarity.y = uint16(bitshift(bitand(allAddr(polarityLogical), yMask), -yShiftBits));
		% X addresses
		data.polarity.x = uint16(bitshift(bitand(allAddr(polarityLogical), xMask), -xShiftBits));
		% Polarity bit
		polBit = 12;
		data.polarity.polarity = bitget(allAddr(polarityLogical), polBit) == 1;
	end

%% Davis frame events

	% NOTE This code currently only handles global shutter readout ...

	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('frame'), importParams.dataTypes))) && any(frameLogical)
		disp('Importing frames ...')
		% These two are defined in the format, but not actually necessary to establish the frame boundaries
		% frameLastEventMask = hex2dec ('FFFFFC00');
		% frameLastEvent = hex2dec     ('80000000'); %starts with biggest address

		frameSampleMask = bin2dec('1111111111');

		frameData = allAddr(frameLogical);
		frameTs = allTs(frameLogical);

		frameX = uint16(bitshift(bitand(frameData, xMask),-xShiftBits));
		frameY = uint16(bitshift(bitand(frameData, yMask),-yShiftBits));
		frameSignal = logical(bitand(frameData, signalOrSpecialMask));
		frameSample = uint16(bitand(frameData, frameSampleMask));

		% In general the ramp of address values could be in either
		% direction and either x or y could be the outer(inner) loop
		% Search for a discontinuity in both x and y simultaneously
		frameXDouble = double(frameX);
		frameYDouble = double(frameY);
		frameXDiscont = abs(frameXDouble(2 : end) - frameXDouble(1 : end - 1)) > 1;
		frameYDiscont = abs(frameYDouble(2 : end) - frameYDouble(1 : end - 1)) > 1;
		frameStarts = [1; find(frameXDiscont & frameYDiscont) + 1; length(frameData) + 1];
		% Now we have the indices of the first sample in each frame, plus
		% an additional index just beyond the end of the array
		numFrames = length(frameStarts) - 1;

		data.frame.reset = false(numFrames, 1);
		data.frame.timeStampStart = zeros(numFrames, 1, 'uint32');
		data.frame.timeStampEnd = zeros(numFrames, 1, 'uint32');
		data.frame.samples = cell(numFrames, 1);
		data.frame.xLength = zeros(numFrames, 1, 'uint16');
		data.frame.yLength = zeros(numFrames, 1, 'uint16');
		data.frame.xPosition = zeros(numFrames, 1, 'uint16');
		data.frame.yPosition = zeros(numFrames, 1, 'uint16');

		for frameIndex = 1 : numFrames
			disp(['Processing frame ' num2str(frameIndex)])
			% All within a frame should be either reset or signal. I could
			% implement a check here to see that that's true, but I haven't
			% done so; rather I just take the first value
			data.frame.reset(frameIndex) = ~frameSignal(frameStarts(frameIndex));

			% in aedat 2 format we don't have the four timestamps of aedat 3 format
			% We expect to find all the same timestamps;
			% nevertheless search for lowest and highest
			data.frame.timeStampStart(frameIndex) = min(frameTs(frameStarts(frameIndex) : frameStarts(frameIndex + 1) - 1));
			data.frame.timeStampEnd(frameIndex) = max(frameTs(frameStarts(frameIndex) : frameStarts(frameIndex + 1) - 1));

			tempXPosition = min(frameX(frameStarts(frameIndex) : frameStarts(frameIndex + 1) - 1));
			data.frame.xPosition(frameIndex) = tempXPosition;
			tempYPosition = min(frameY(frameStarts(frameIndex) : frameStarts(frameIndex + 1) - 1));
			data.frame.yPosition(frameIndex) = tempYPosition;
			data.frame.xLength(frameIndex) = max(frameX(frameStarts(frameIndex) : frameStarts(frameIndex + 1) - 1)) - data.frame.xPosition(frameIndex) + 1;
			data.frame.yLength(frameIndex) = max(frameY(frameStarts(frameIndex) : frameStarts(frameIndex + 1) - 1)) - data.frame.yPosition(frameIndex) + 1;
			% If we worked out which way the data is ramping in each
			% direction, and if we could exclude data loss, then we could
			% do some nice clean matrix transformations; but I'm just going
			% to iterate through the samples, putting them in the right
			% place in the array according to their address

			% first create a temporary array - there is no concept of
			% colour channels in aedat2
            %{
            Code before I learned about 'accumarray':
			tempSamples = zeros(data.frame.yLength(frameIndex), data.frame.xLength(frameIndex), 'uint16');
			for sampleIndex = frameStarts(frameIndex) : frameStarts(frameIndex + 1) - 1
				tempSamples(frameY(sampleIndex) - data.frame.yPosition(frameIndex) + 1, ...
							frameX(sampleIndex) - data.frame.xPosition(frameIndex) + 1) ...
							= frameSample(sampleIndex);
			end
			data.frame.samples{frameIndex} = tempSamples;
            %}
            sampleIndexRange = frameStarts(frameIndex) : frameStarts(frameIndex + 1) - 1;
            data.frame.samples{frameIndex} = accumarray(...
                [frameY(sampleIndexRange) - data.frame.yPosition(frameIndex) + 1 ...
                 frameX(sampleIndexRange) - data.frame.xPosition(frameIndex) + 1], ...
                 frameSample(sampleIndexRange), ...
                 [data.frame.yLength(frameIndex) data.frame.xLength(frameIndex)]);
		end
		% By default, subtract the reset read
		if isfield(importParams, 'subtractResetRead')
            subtractResetRead = importParams.subtractResetRead;
        else
            subtractResetRead = true;
        end
		if subtractResetRead && isfield(data.frame, 'reset')
            disp('Performing frame subtraction ...')
			% Make a second pass through the frames, subtracting reset
			% reads from signal reads
			frameCount = 0;
			for frameIndex = 1 : numFrames
				if data.frame.reset(frameIndex)
					resetFrame = data.frame.samples{frameIndex};
					resetXPosition = data.frame.xPosition(frameIndex);
					resetYPosition = data.frame.yPosition(frameIndex);
					resetXLength = data.frame.xLength(frameIndex);
					resetYLength = data.frame.yLength(frameIndex);
				else
					frameCount = frameCount + 1;
					% If a resetFrame has not yet been found,
					% push through the signal frame as is
					if ~exist('resetFrame', 'var')
						data.frame.samples{frameCount} ...
							= data.frame.samples{frameIndex};
					else
						% If the resetFrame and signalFrame are not the same size,
						% don't attempt subtraction
						% (there is probably a cleaner solution than this - could be improved)
						if resetXPosition ~= data.frame.xPosition(frameIndex) ...
							|| resetYPosition ~= data.frame.yPosition(frameIndex) ...
							|| resetXLength ~= data.frame.xLength(frameIndex) ...
							|| resetYLength ~= data.frame.yLength(frameIndex)
							data.frame.samples{frameCount} ...
								= data.frame.samples{frameIndex};
						else
							% Do the subtraction
							data.frame.samples{frameCount} ...
								= resetFrame - data.frame.samples{frameIndex};
                                			% This operation was on unsigned integers, set negatives to zero
                                			data.frame.samples{frameCount}(data.frame.samples{frameCount} > 32767) = 0;

						end
						% Copy over the rest of the info
						data.frame.xPosition(frameCount) = data.frame.xPosition(frameIndex);
						data.frame.yPosition(frameCount) = data.frame.yPosition(frameIndex);
						data.frame.xLength(frameCount) = data.frame.xLength(frameIndex);
						data.frame.yLength(frameCount) = data.frame.yLength(frameIndex);
						data.frame.timeStampStart(frameCount) = data.frame.timeStampStart(frameIndex);
						data.frame.timeStampEnd(frameCount) = data.frame.timeStampEnd(frameIndex);
					end
				end
			end
			% Clip the arrays
			data.frame.xPosition = data.frame.xPosition(1 : frameCount);
			data.frame.yPosition = data.frame.yPosition(1 : frameCount);
			data.frame.xLength = data.frame.xLength(1 : frameCount);
			data.frame.yLength = data.frame.yLength(1 : frameCount);
			data.frame.timeStampStart = data.frame.timeStampStart(1 : frameCount);
			data.frame.timeStampEnd = data.frame.timeStampEnd(1 : frameCount);
			data.frame.samples = data.frame.samples(1 : frameCount);
			data.frame = rmfield(data.frame, 'reset'); % reset is no longer needed
		end
	end
%% Davis IMU events

		% These come in blocks of 7, for the 7 different values produced in
		% a single sample; the following code recomposes these
		% 7 words are sent in series, these being 3 axes for accel, temperature, and 3 axes for gyro
	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('imu6'), importParams.dataTypes))) && any(imuLogical)
		disp('Importing imu6 events ...')
		if mod(nnz(imuLogical), 7) > 0
			%error('The number of IMU samples is not divisible by 7, so IMU samples are not interpretable')
            % There's a problem here, but chop off the last words and hope
            % for the best ...
		end
		data.imu6.timeStamp = allTs(imuLogical);
		data.imu6.timeStamp = data.imu6.timeStamp(1 : 7 : end - mod(nnz(imuLogical), 7));

		%Conversion factors
		accelScale = 1/8192;
		gyroScale = 1/65.5;
		temperatureScale = 1/340;
		temperatureOffset=35;

		imuDataMask = hex2dec('0FFFF000');
		imuDataShiftBits = 12;
        rawData = bitshift(bitand(allAddr(imuLogical), imuDataMask), -imuDataShiftBits);
        % Now RawData is uint32, but the data is a signed int in the ls 16 bits.
		rawData = uint16(rawData); % Now any negative indicators are in the MSB
        rawData = typecast(rawData, 'int16'); % Now the same numbers are interpretted as negatives
        rawData = double(rawData);
        % Now the data is floating point, ready for conversion to physical units

		data.imu6.accelX			= rawData(1 : 7 : end - mod(nnz(imuLogical), 7)) * accelScale;
		data.imu6.accelY			= rawData(2 : 7 : end - mod(nnz(imuLogical), 7)) * accelScale;
		data.imu6.accelZ			= rawData(3 : 7 : end - mod(nnz(imuLogical), 7)) * accelScale;
		data.imu6.temperature	= rawData(4 : 7 : end - mod(nnz(imuLogical), 7)) * temperatureScale + temperatureOffset;
		data.imu6.gyroX			= rawData(5 : 7 : end - mod(nnz(imuLogical), 7)) * gyroScale;
		data.imu6.gyroY			= rawData(6 : 7 : end - mod(nnz(imuLogical), 7)) * gyroScale;
		data.imu6.gyroZ			= rawData(7 : 7 : end - mod(nnz(imuLogical), 7)) * gyroScale;

	end

	% If you want to do chip-specific address shifts or subtractions,
	% this would be the place to do it.

elseif strfind(info.source, 'SecDvs')
	if (~isfield(importParams, 'dataTypes') || any(cellfun(cellFind('polarity'), importParams.dataTypes)))

		data.polarity.timeStamp = allTs;
		% Y addresses
		yMask = hex2dec('ff800');
		yShiftBits = 11;
		data.polarity.y = uint16(bitshift(bitand(allAddr, yMask), -yShiftBits));
		% X addresses
		xMask = hex2dec('7fe');
		xShiftBits = 1;
		data.polarity.x = uint16(bitshift(bitand(allAddr, xMask), -xShiftBits));
		% Polarity bit
		polBit = 1;
		data.polarity.polarity = bitget(allAddr, polBit) == 1;
    end
end

%% Pack data

% aedat.importParams is already there and should be unchanged
aedat.info = info;
aedat.data = data;

%% Find first and last time stamps

aedat = FindFirstAndLastTimeStamps(aedat);

%% Add NumEvents field for each data type

aedat = NumEventsByType(aedat);

disp('Import finished')

author: panpanfei


text:
function s = zeroSingFun()
%ZEROSINGFUN   Constructs the zero SINGFUN. The output SINGFUN object has a
%   zero smooth part with no singularities at any end points.
%
% See also SINGFUN.

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

% Create an empty SINGFUN object:
s = singfun();

% Create a zero smooth part:
s.smoothPart = singfun.constructSmoothPart(0, [], []);

% No singularities at any end points:
s.exponents = [0, 0];

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
function q = divine(im)
% DIIVINE Software release.
%
%
% ========================================================================
%
% -----------COPYRIGHT NOTICE STARTS WITH THIS LINE------------
% Copyright (c) 2010 The University of Texas at Austin
% All rights reserved.
%
% Permission is hereby granted, without written agreement and without license or royalty fees, to use, copy,
% modify, and distribute this code (the source files) and its documentation for
% any purpose, provided that the copyright notice in its entirety appear in all copies of this code, and the
% original source of this code, Laboratory for Image and Video Engineering (LIVE, http://live.ece.utexas.edu)
% and Center for Perceptual Systems (CPS, http://www.cps.utexas.edu) at the University of Texas at Austin (UT Austin,
% http://www.utexas.edu), is acknowledged in any publication that reports research using this code. The research
% is to be cited in the bibliography as:
%
% 1. A. K. Moorthy and A. C. Bovik, "Blind Image Quality Assessment: From Natural
% Scene Statistics to Perceptual Quality", IEEE Transactions on Image Processing, to appear (2011).
%
% 2. A. K. Moorthy and A. C. Bovik, "DIVINE Software Release",
% URL: http://live.ece.utexas.edu/research/quality/DIIVINE_release.zip, 2010
%
% IN NO EVENT SHALL THE UNIVERSITY OF TEXAS AT AUSTIN BE LIABLE TO ANY PARTY FOR DIRECT, INDIRECT, SPECIAL, INCIDENTAL,
% OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OF THIS DATABASE AND ITS DOCUMENTATION, EVEN IF THE UNIVERSITY OF TEXAS
% AT AUSTIN HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
%
% THE UNIVERSITY OF TEXAS AT AUSTIN SPECIFICALLY DISCLAIMS ANY WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
% WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE DATABASE PROVIDED HEREUNDER IS ON AN "AS IS" BASIS,
% AND THE UNIVERSITY OF TEXAS AT AUSTIN HAS NO OBLIGATION TO PROVIDE MAINTENANCE, SUPPORT, UPDATES, ENHANCEMENTS, OR MODIFICATIONS.
%
% -----------COPYRIGHT NOTICE ENDS WITH THIS LINE------------%
%
% Author  : Anush Krishna Moorthy
% Version : 1.1
%
% The authors are with the Laboratory for Image and Video Engineering
% (LIVE), Department of Electrical and Computer Engineering, The
% University of Texas at Austin, Austin, TX.
%
% Kindly report any suggestions or corrections to anushmoorthy@gmail.com
%
% ========================================================================
%
% This is a demonstration of the Distortion Identification based image Verity and INtegrity Evaluation (DIVINE) index.
% It is an implementation of the BIQI in the reference.
% The algorithm is described in:
% A. K. Moorthy and A. C. Bovik, "Blind Image Quality Assessment: From Natural
% Scene Statistics to Perceptual Quality",  IEEE Transactions on Image Processing, to appear (2011).
%
% You can change this program as you like and use it anywhere, but please
% refer to its original source (cite our paper and our web page at
% http://live.ece.utexas.edu/research/quality/DIIVINE_release.zip).
%
% Input : A test 8bits/pixel grayscale image loaded in a 2-D array
% Output: A quality score of the image. The score typically has a value
%        between 0 and 100 (0 represents the best quality, 100 the worst).
%
% Usage:
%
% 1. Load the image, for example
%
%   image = rgb2gray(imread('testimage.jpg'));
%
% 2. Call this function to calculate the quality score:
%
%    quality = divine(image)
%
% Dependencies:
% Steerable Pyramid Toolbox, Download from: http://www.cns.nyu.edu/~eero/steerpyr/
% LibSVM package for MATLAB, Download from: http://www.csie.ntu.edu.tw/~cjlin/libsvm/
% You may need the MATLAB Image Processing Toolbox
%
% Dependencies--
%
% MATLAB files:  ssim_index_new.m, norm_sender_normalized.m, find_spatial_hist_fast.m, divine_overall_quality.m
%               divine_feature_extract.m,  map_matrix_to_closest_vec.m (provided with release)
%
% Data files: data_live_trained.mat (provided with release)
%
% This code has been tested on Windows and Mac OSX (Snow Leopard)
%
% ========================================================================%
% Note on training:
% This release version of BIQI was trained on the entire LIVE database.
%
%

import divine.*

f = divine_feature_extract(im);
q = divine_overall_quality(f);

author: dsoellinger


text:
function [ grid_weight, grid_point ] = sparse_grid_own ( dim_num, ...
  level_max, rule, point_num )

%*****************************************************************************80
%
%% SPARSE_GRID_OWN computes a sparse grid based on an OWN 1D rule.
%
%  Discussion:
%
%    The 1D quadrature rule is assumed to be Open Weakly Nested.
%    Such rules include Gauss Hermite and Gauss Legendre rules.
%
%    A Smolyak construction is used to create a multidimensional sparse grid.
%
%    The user specifies:
%    * the spatial dimension of the quadrature region,
%    * the level that defines the Smolyak grid,
%    * the rule;
%    * the number of points.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    02 July 2008
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    Fabio Nobile, Raul Tempone, Clayton Webster,
%    A Sparse Grid Stochastic Collocation Method for Partial Differential
%    Equations with Random Input Data,
%    SIAM Journal on Numerical Analysis,
%    Volume 46, Number 5, 2008, pages 2309-2345.
%
%  Parameters:
%
%    Input, integer DIM_NUM, the spatial dimension.
%
%    Input, integer LEVEL_MAX, controls the size of the
%    sparse grid.
%
%    Input, integer RULE, the index of the rule.
%    1, "CC", Clenshaw Curtis Closed Fully Nested rule.
%    2, "F1", Fejer 1 Open Fully Nested rule.
%    3, "F2", Fejer 2 Open Fully Nested rule.
%    4, "GP", Gauss Patterson Open Fully Nested rule.
%    5, "GL", Gauss Legendre Open Weakly Nested rule.
%    6, "GH", Gauss Hermite Open Weakly Nested rule.
%    7, "LG", Gauss Laguerre Open Non Nested rule.
%
%    Input, integer POINT_NUM, the number of points in the grid,
%    as determined by LEVELS_INDEX_SIZE_OWN.
%
%    Output, real GRID_WEIGHT(POINT_NUM), the weights.
%
%    Output, real GRID_POINT(DIM_NUM,POINT_NUM), the points.
%
  grid_weight(1:point_num) = 0.0;
  grid_point = zeros ( dim_num, point_num );
%
%  The outer loop generates LEVELs from LEVEL_MIN to LEVEL_MAX.
%
  point_num2 = 0;

  level_min = max ( 0, level_max + 1 - dim_num );

  if ( dim_num == 1 )
    level_min2 = level_min;
  else
    level_min2 = 0;
  end

  for level = level_min2 : level_max
%
%  The middle loop generates the next partition LEVEL_1D(1:DIM_NUM)
%  that adds up to LEVEL.
%
    level_1d = [];
    more = 0;
    h = 0;
    t = 0;

    while ( 1 )

      [ level_1d, more, h, t ] = comp_next ( level, dim_num, level_1d, more, h, t );
%
%  Transform each 1D level to a corresponding 1D order.
%  The relationship is the same as for other OPEN rules.
%  The GL rule differs from the other OPEN rules only in the nesting behavior.
%
      order_1d = level_to_order_open ( dim_num, level_1d );

      grid_base2(1:dim_num) = round ( ( order_1d(1:dim_num) - 1 ) / 2 );
%
%  The product of the 1D orders gives us the number of points in this subgrid.
%
      order_nd = prod ( order_1d(1:dim_num) );
%
%  Compute the weights for this product grid.
%
      grid_weight2 = product_weights ( dim_num, order_1d, order_nd, rule );
%
%  Now determine the coefficient of the weight.
%
      coeff = r8_mop ( level_max - level ) ...
        * r8_choose ( dim_num - 1, level_max - level );
%
%  The inner (hidden) loop generates all points corresponding to given grid.
%  The grid indices will be between -M to +M, where 2*M + 1 = ORDER_1D(DIM).
%
      grid_index2 = multigrid_index_own ( dim_num, order_1d, order_nd );
%
%  Determine the first level of appearance of each of the points.
%  This allows us to flag certain points as being repeats of points
%  generated on a grid of lower level.
%
%  This is SLIGHTLY tricky.
%
      grid_level = index_level_own ( level, level_max, dim_num, order_nd, ...
        grid_index2, grid_base2 );
%
%  Only keep those points which first appear on this level.
%
      for point = 1 : order_nd
%
%  Either a "new" point (increase count, create point, create weight)
%
        if ( grid_level(point) == level )

          point_num2 = point_num2 + 1;

          if ( point_num < point_num2 )
            fprintf ( 1, '\n' );
            fprintf ( 1, 'SPARSE_GRID_OWN - Fatal error!\n' );
            fprintf ( 1, ...
            '  Exceeding maximum point index POINT_NUM = %d\n', point_num );
            error ( 'SPARSE_GRID_OWN - Fatal error!' );
          end

          if ( rule == 5 )
            grid_point(1:dim_num,point_num2) = gl_abscissa ( dim_num, 1, ...
              grid_index2(1:dim_num,point), grid_base2(1:dim_num) );
          elseif ( rule == 6 )
            grid_point(1:dim_num,point_num2) = gh_abscissa ( dim_num, 1, ...
              grid_index2(1:dim_num,point), grid_base2(1:dim_num) );
          else
            fprintf ( 1, '\n' );
            fprintf ( 1, 'SPARSE_GRID_OWN - Fatal error!\n' );
            fprintf ( 1, '  Unrecognized rule number = %d\n', rule );
            error ( 'SPARSE_GRID_OWN - Fatal error!' );
          end

          if ( level_min <= level )
            grid_weight(point_num2) = coeff * grid_weight2(point);
          end
%
%  or an already existing point (create point temporarily, find match,
%  add weight to matched point's weight).
%
        else

          if ( level_min <= level )

            if ( rule == 5 )
              grid_point_temp(1:dim_num) = gl_abscissa ( dim_num, 1, ...
                grid_index2(1:dim_num,point), grid_base2(1:dim_num) );
            elseif ( rule == 6 )
              grid_point_temp(1:dim_num) = gh_abscissa ( dim_num, 1, ...
                grid_index2(1:dim_num,point), grid_base2(1:dim_num) );
            else
              fprintf ( 1, '\n' );
              fprintf ( 1, 'SPARSE_GRID_OWN - Fatal error!\n' );
              fprintf ( 1, '  Unrecognized rule number = %d\n', rule );
              error ( 'SPARSE_GRID_OWN - Fatal error!' );
            end

            point3 = -1;

            for point2 = 1 : point_num2
              if ( all ( grid_point(1:dim_num,point2) == grid_point_temp(1:dim_num)' ) )
                point3 = point2;
                break
              end
            end

            if ( point3 == -1 )
              fprintf ( 1, '\n' );
              fprintf ( 1, 'SPARSE_GRID_OWN - Fatal error!\n' );
              fprintf ( 1, '  Could not match point.\n' );
              error ( 'SPARSE_GRID_OWN - Fatal error!' );
            end

            grid_weight(point3) = grid_weight(point3) ...
              + coeff * grid_weight2(point);

          end

        end

      end

      if ( ~more )
        break
      end

    end

  end

  if ( point_num2 < point_num )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'SPARSE_GRID_OWN - Fatal error!\n' );
    fprintf ( 1,'  Set fewer points than POINT_NUM = %d\n', point_num );
    error ( 'SPARSE_GRID_OWN - Fatal error!' );
  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function disp(t,name)
%DISP Command window display of a ttensor.
%
%   DISP(T) displays a ttensor with no name.
%
%   DISP(T,NAME) display a ttensor with the given name.
%
%   See also TTENSOR.
%
%MATLAB Tensor Toolbox.
%Copyright 2010, Sandia Corporation.

% This is the MATLAB Tensor Toolbox by Brett Bader and Tamara Kolda.
% http://csmr.ca.sandia.gov/~tgkolda/TensorToolbox.
% Copyright (2010) Sandia Corporation. Under the terms of Contract
% DE-AC04-94AL85000, there is a non-exclusive license for use of this
% work by or on behalf of the U.S. Government. Export of this data may
% require a license from the United States Government.
% The full license terms can be found in tensor_toolbox/LICENSE.txt
% $Id: disp.m,v 1.14 2010/03/19 23:46:31 tgkolda Exp $

if ~exist('name','var')
    name = 'ans';
end

fprintf(1,'%s is a ttensor of size %s\n', name, tt_size2str(size(t)));
disp(t.core, sprintf('\t%s.core',name));

for j = 1 : ndims(t)
    fprintf('\t%s.U{%d} = \n', name, j);
    output = tt_matrix2cellstr(t.u{j});
    fprintf('\t\t%s\n',output{:});
end

author: andrewssobral

repo: mtt

save_path: github-repos/MATLAB/andrewssobral-mtt


text:
% STD_SELCOMP - Helper function for STD_ERPPLOT, STD_SPECPLOT
%                       and STD_ERSPPLOT to select specific
%                       components prior to plotting.
% Usage:
%  >> std_selcomp( STUDY, data, cluster, setinds, compinds, comps)
%
% Inputs:
%  STUDY -  EEGLAB STUDY structure.
%  data  -  [cell array] mean data for each subject group and/or data
%           condition. For example, to compute mean ERPs statistics from a
%           STUDY for epochs of 800 frames in two conditions from three
%           groups of 12 subjects,
%           >> data = { [800x12] [800x12] [800x12];... % 3 groups, cond 1
%                       [800x12] [800x12] [800x12] };  % 3 groups, cond 2
% cluster - [integer] cluster index
% setinds - [cell array] set indices for each of the last dimension of the
%           data cell array.
%           >> setinds = { [12] [12] [12];... % 3 groups, cond 1
%                          [12] [12] [12] };  % 3 groups, cond 2
% compinds - [cell array] component indices for each of the last dimension
%           of the data cell array.
%           >> compinds = { [12] [12] [12];... % 3 groups, cond 1
%                           [12] [12] [12] };  % 3 groups, cond 2
% comps    - [integer] find and select specific component index in array
%
% Output:
%  data       - [cell array] data array with the subject or component selected
%  subject    - [string] subject name (for component selection)
%  comp_names - [cell array] component names (for component selection)
%
% Author: Arnaud Delorme, CERCO, CNRS, 2006-
%
% See also: STD_ERPPLOT, STD_SPECPLOT and STD_ERSPPLOT

% Copyright (C) 2006 Arnaud Delorme
%
% 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 [data, subject, comp_names] = std_selcomp(STUDY, data, clust, setinds, compinds, compsel)

if nargin < 2
    help std_selcomp;
    return;
end

optndims = ndims(data{1});
comp_names = {};
subject    = '';

% find and select group
% ---------------------
if isempty(compsel), return; end
sets   = STUDY.cluster(clust).sets(:,compsel);
comps  = STUDY.cluster(clust).comps(compsel);
%grp    = STUDY.datasetinfo(sets(1)).group;
%grpind = strmatch( grp, STUDY.group );
%if isempty(grpind), grpind = 1; end
%data = data(:,grpind);

% find component
% --------------
for c = 1:length(data(:))
    rminds = 1:size(data{c},optndims);

    for ind = length(compinds{c}):-1:1
        setindex = STUDY.design(STUDY.currentdesign).cell(setinds{c}(ind)).dataset;
        if compinds{c}(ind) == comps && any(setindex == sets)
            rminds(ind) = [];
        end
    end

    if optndims == 2
        data{c}(:,rminds) = []; %2-D
    elseif optndims == 3
        data{c}(:,:,rminds) = []; %3-D
    else
        data{c}(:,:,:,rminds) = []; %3-D
    end
    comp_names{c,1} = comps;
end
% for c = 1:size(data,1)
%     for ind = 1:length(compinds{1,grpind})
%         if compinds{1,grpind}(ind) == comps && any(setinds{1,grpind}(ind) == sets)
%             if optndims == 2
%                  data{c} = data{c}(:,ind);
%             else data{c} = data{c}(:,:,ind);
%             end
%             comp_names{c,1} = comps;
%         end
%     end
% end
subject = STUDY.datasetinfo(sets(1)).subject;

author: sccn

repo: eeglab

save_path: github-repos/MATLAB/sccn-eeglab


text:
function [GB] = PB2GB(PB)
% Convert computery things from petabytes to gigabytes.
% Chad A. Greene 2012
GB = PB*1048576;

author: Sable

repo: mcbench-benchmarks


text:
function [cal] = kcal2cal(kcal)
% Convert energy or work from kilocalories to calories.
% Note: kilocalories are the same thing as nutritional Calories.
% Chad A. Greene 2012
cal = kcal*1000;


text:
function [Nm] = ftlb2Nm(ftlb)
% Convert energy or work from foot-pounds to newtons-meters.
% Chad A. Greene 2012
Nm = ftlb*1.3558179;


text:
%% pro_Create: Create a new project structure
%
% Usage:
%   pro = pro_Create()
%
% Inputs:
%
% Output:
%     pro                project structure
%                          .Inputs.pdfs: cell-array of the model inputs with the pdf handles
%                          .Inputs.Names: cell-array with the input names
%                          .N: scalar, number of samples of crude Monte Carlo
%                          .Model.handle: handle to the model function
%                          .Model.Name: string name of the model
%
% ------------------------------------------------------------------------
% See also
%
% Author : Flavio Cannavo'
% e-mail: flavio(dot)cannavo(at)gmail(dot)com
% Release: 1.0
% Date   : 28-01-2011
%
% History:
% 1.0  28-01-2011  First release.
%%

function pro = pro_Create()

pro.Inputs.pdfs = {};
pro.Inputs.Names = {};
pro.N = 10000;
pro.Model.handle = [];
pro.Model.Name = [];


text:
function varargout = Sim_Sensitivity_Analysis_GUI(varargin)
% SIM_SENSITIVITY_ANALYSIS_GUI MATLAB code for Sim_Sensitivity_Analysis_GUI.fig
%      SIM_SENSITIVITY_ANALYSIS_GUI, by itself, creates a new SIM_SENSITIVITY_ANALYSIS_GUI or raises the existing
%      singleton*.
%
%      H = SIM_SENSITIVITY_ANALYSIS_GUI returns the handle to a new SIM_SENSITIVITY_ANALYSIS_GUI or the handle to
%      the existing singleton*.
%
%      SIM_SENSITIVITY_ANALYSIS_GUI('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in SIM_SENSITIVITY_ANALYSIS_GUI.M with the given input arguments.
%
%      SIM_SENSITIVITY_ANALYSIS_GUI('Property','Value',...) creates a new SIM_SENSITIVITY_ANALYSIS_GUI or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before Sim_Sensitivity_Analysis_GUI_OpeningFcn gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to Sim_Sensitivity_Analysis_GUI_OpeningFcn via varargin.
%
%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one
%      instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help Sim_Sensitivity_Analysis_GUI

% Last Modified by GUIDE v2.5 27-Jul-2017 17:23:32
if moxunit_util_platform_is_octave, warndlg('Graphical user interface not available on octave... use command lines instead'); return; end

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
    'gui_Singleton',  gui_Singleton, ...
    'gui_OpeningFcn', @Sim_Sensitivity_Analysis_GUI_OpeningFcn, ...
    'gui_OutputFcn',  @Sim_Sensitivity_Analysis_GUI_OutputFcn, ...
    'gui_LayoutFcn',  [] , ...
    'gui_Callback',   []);
if nargin && ischar(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT

% --- Executes just before Sim_Sensitivity_Analysis_GUI is made visible.
function Sim_Sensitivity_Analysis_GUI_OpeningFcn(hObject, eventdata, handles, varargin)
handles.output = hObject;
handles.Model = varargin{1};
if ~isfield(handles,'opened')
    % clear axe
    axes(handles.SimVaryAxe)

    % fill table
    Nparam=length(handles.Model.xnames);
    FitOptTable(:,1)=handles.Model.xnames(:);
    if isprop(handles.Model,'fx') && ~isempty(handles.Model.fx),    FitOptTable(:,2)=mat2cell(~logical(handles.Model.fx(:)),ones(Nparam,1)); end

    if isprop(handles.Model,'ub') && ~isempty(handles.Model.ub)
        FitOptTable(:,4)=mat2cell(handles.Model.lb(:),ones(Nparam,1));
        FitOptTable(:,5)=mat2cell(handles.Model.ub(:),ones(Nparam,1));
    end
    if isprop(handles.Model,'st') && ~isempty(handles.Model.st)
        FitOptTable(:,3)=mat2cell(handles.Model.st(:),ones(Nparam,1));
    elseif size(FitOptTable,2)==5
        FitOptTable(:,3) = mat2cell(mean(cat(2,handles.Model.lb(:),handles.Model.ub(:)),2),ones(Nparam,1));
    else
        FitOptTable(:,3) = mat2cell(ones(Nparam,1),ones(Nparam,1));
    end

    if size(FitOptTable,2)<5
        FitOptTable(:,4)=mat2cell(cell2mat(FitOptTable(:,3))/2,ones(Nparam,1));
        FitOptTable(:,5)=mat2cell(cell2mat(FitOptTable(:,3))*2,ones(Nparam,1));
    end
    set(handles.SimVaryOptTable,'Data',FitOptTable)
    % fill parameters
    set(handles.SimVaryPlotX,'String',handles.Model.xnames')
    set(handles.SimVaryPlotY,'String',handles.Model.xnames')

    % Options
    if isprop(handles.Model,'Sim_Single_Voxel_Curve_buttons')
        opts = handles.Model.Sim_Single_Voxel_Curve_buttons;
    else
        opts = {'SNR',50};
    end
    if isprop(handles.Model,'Sim_Sensitivity_Analysis_buttons'), opts = cat(2,opts,handles.Model.Sim_Sensitivity_Analysis_buttons);
    else
        opts = cat(2,opts,{'# of run',20});
    end

    handles.options = GenerateButtonsWithPanels(opts,handles.OptionsPanel);
    handles.opened = 1;

    % Create CALLBACK for buttons
    ff = fieldnames(handles.options);
    for ii=1:length(ff)
        %set(handles.OptionsPanel_handle.(ff{ii}),'Callback',@(src,event) ModelOptions_Callback(handles));
        switch get(handles.options.(ff{ii}),'Style')
            case 'togglebutton'
                set(handles.options.(ff{ii}),'Callback',@(src,event) ModelSimOptions_Callback(handles));
        end
    end

end
% Update handles structure
guidata(hObject, handles);

% UIWAIT makes Sim_Sensitivity_Analysis_GUI wait for user response (see UIRESUME)
% uiwait(handles.Simu);

% --- Outputs from this function are returned to the command line.
function varargout = Sim_Sensitivity_Analysis_GUI_OutputFcn(hObject, eventdata, handles)
% varargout  cell array for returning output args (see VARARGOUT);
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure
varargout{1} = handles.output;

% --------------------------------------------------------------------
function FileMenu_Callback(hObject, eventdata, handles)
% hObject    handle to FileMenu (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% --------------------------------------------------------------------
function OpenMenuItem_Callback(hObject, eventdata, handles)
% hObject    handle to OpenMenuItem (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
file = uigetfile('*.fig');
if ~isequal(file, 0)
    open(file);
end

% --------------------------------------------------------------------
function PrintMenuItem_Callback(hObject, eventdata, handles)
% hObject    handle to PrintMenuItem (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
printdlg(handles.Simu)

% --------------------------------------------------------------------
function CloseMenuItem_Callback(hObject, eventdata, handles)
% hObject    handle to CloseMenuItem (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
selection = questdlg(['Close ' get(handles.Simu,'Name') '?'],...
    ['Close ' get(handles.Simu,'Name') '...'],...
    'Yes','No','Yes');
if strcmp(selection,'No')
    return;
end

delete(handles.Simu)

% --- Executes on selection change in popupmenu1.
function popupmenu1_Callback(hObject, eventdata, handles)
% hObject    handle to popupmenu1 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: contents = get(hObject,'String') returns popupmenu1 contents as cell array
%        contents{get(hObject,'Value')} returns selected item from popupmenu1

% --- Executes during object creation, after setting all properties.
function popupmenu1_CreateFcn(hObject, eventdata, handles)
% hObject    handle to popupmenu1 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: popupmenu controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

set(hObject, 'String', {'plot(rand(5))', 'plot(sin(1:0.01:25))', 'bar(1:.5:10)', 'plot(membrane)', 'surf(peaks)'});

% --- Executes on selection change in SimVaryPlotX.
function SimVaryPlotX_Callback(hObject, eventdata, handles)
SimVaryPlotResults(handles)

% --- Executes during object creation, after setting all properties.
function SimVaryPlotX_CreateFcn(hObject, eventdata, handles)
% hObject    handle to SimVaryPlotX (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: popupmenu controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

% --- Executes on selection change in SimVaryPlotY.
function SimVaryPlotY_Callback(hObject, eventdata, handles)
SimVaryPlotResults(handles)

% --- Executes during object creation, after setting all properties.
function SimVaryPlotY_CreateFcn(hObject, eventdata, handles)
% hObject    handle to SimVaryPlotY (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: popupmenu controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

% --- Executes on button press in SimVaryUpdate.
function SimVaryUpdate_Callback(hObject, eventdata, handles)
Model_new = getappdata(0,'Model');
if ~isempty(Model_new) && strcmp(class(Model_new),class(handles.Model))
    handles.Model = Model_new;
end
FitOptTable = get(handles.SimVaryOptTable,'Data'); FitOptTable(:,2)=mat2cell(~[FitOptTable{:,2}]',ones(size(FitOptTable,1),1), 1);
FitOptTable = cell2struct(FitOptTable,{'xnames','fx','st','lb','ub'},2);
Opts = button_handle2opts(handles.options);
handles.SimVaryResults = handles.Model.Sim_Sensitivity_Analysis(FitOptTable,Opts);
SetSimVaryResults(handles)
guidata(hObject, handles);

% --- Executes on button press in Save.
function Save_Callback(hObject, eventdata, handles)
if isfield(handles,'SimVaryResults')
    Method = class(handles.Model);
    [FileName,PathName] = uiputfile([Method '_SimResults.mat']);
    if PathName == 0, return; end
    SimVaryResults = handles.SimVaryResults;
    save(fullfile(PathName,FileName),'SimVaryResults')
end

% --- Executes on button press in Options panel.
function ModelSimOptions_Callback(handles)
xtable = get(handles.SimVaryOptTable,'Data');
x=cell2mat(xtable(~cellfun(@isempty,xtable(:,3)),3))';
xnew = SimOpt(handles.Model,x,button_handle2opts(handles.options));
if ~isempty(xnew) % update the ParamTable in the GUI
    Nparam = length(handles.Model.xnames);
    xtable(1:Nparam,3) = mat2cell(xnew',ones(Nparam,1));
    set(handles.SimVaryOptTable,'Data',xtable);
end

function Load_Callback(hObject, eventdata, handles)
Method = class(handles.Model);
[FileName,PathName] = uigetfile([Method '_SimResults.mat']);
if PathName == 0, return; end
load(fullfile(PathName,FileName));
handles.SimVaryResults = SimVaryResults;
SetSimVaryResults(handles)
guidata(hObject, handles);

function SetSimVaryResults(handles)
ff=fieldnames(handles.SimVaryResults);
set(handles.SimVaryPlotX,'String',ff);
ff=fieldnames(handles.SimVaryResults.(ff{1}));
set(handles.SimVaryPlotY,'String',ff(~ismember(ff,{'x','fit'})));
SimVaryPlotResults(handles)

function SimVaryPlotResults(handles)
if isfield(handles,'SimVaryResults')
    Xaxis = get(handles.SimVaryPlotX,'String'); Xaxis = Xaxis{get(handles.SimVaryPlotX,'Value')};
    Yaxis = get(handles.SimVaryPlotY,'String'); Yaxis = Yaxis{get(handles.SimVaryPlotY,'Value')};
    SimVaryPlot(handles.SimVaryResults,Xaxis,Yaxis)
end

function SimVaryOptTable_CellEditCallback(hObject, eventdata, handles)

author: qMRLab

repo: qMRLab

save_path: github-repos/MATLAB/qMRLab-qMRLab


text:
%--------------------------------------------------------------------------
%this function is to transfer image to data suitable for transmission
%
%Chen Zhifeng
%2007-05-19
%zhifeng@ecel.ufl.edu
%%I left this function as is since I'm clueless and a neophyte in image
%processing. JC 7/16/08
%--------------------------------------------------------------------------

function [data, row_im, col_im, third_im] = image2data(im, M)
% im = imread('photo.bmp');
% M =4;
[row_im, col_im, third_im] =size(im);
V_im = zeros(1, row_im*col_im*third_im);
for i=1:third_im,
    for j=1:col_im,
        V_im((i-1)*row_im*col_im+(j-1)*row_im +1 : (i-1)*row_im*col_im+j*row_im) = im(:,j,i);
    end
end

%I use dec2base here, then if M>8, for example M=16, the result string may
%include characters, this is not appropriate to use str2num below. However,
%due to the time limit, I will not add function here to deal with this.
%Actually, this may be done by dec2bin function.
Tm = dec2base(V_im,M);
[row, col] = size(Tm);
% data = [];
% for i = 1:col,
%     data = [data; Tm(:,i)];
% end
%
%data=[];
for i =1:col,
    data(row*(i-1)+1:row*i)=str2num(Tm(:,i));   %very important to add str2num here
end
%Tm = str2num(Tm_char);

%data = str2num(data);


text:
function [textures] = getGLSZMtextures_temp(GLSZM)
% -------------------------------------------------------------------------
% function [textures] = getGLSZMtextures_temp(GLSZM)
% -------------------------------------------------------------------------
% DESCRIPTION:
% This function computes texture features from an input Gray-Level Size
% Zone Matrix (GLSZM).
% -------------------------------------------------------------------------
% REFERENCES:
% [1] Galloway, M. M. (1975). Texture analysis using gray level run lengths.
%     Computer Graphics and Image Processing, 4(2), 172–179.
% [2] Chu, A., Sehgal, C. M., & Greenleaf, J. F. (1990). Use of gray value
%     distribution of run lengths for texture analysis. Pattern Recognition
%     Letters, 11(6), 415-419.
% [3] Dasarathy, B. V., & Holder, E. B. (1991). Image characterizations
%     based on joint gray level-run length distributions. Pattern
%     Recognition Letters, 12(8), 497-502.
% [4] Thibault, G., Fertil, B., Navarro, C., Pereira, S., Cau, P., Levy,
%     N., Mari, J.-L. (2009). Texture Indexes and Gray Level Size Zone
%     Matrix. Application to Cell Nuclei Classification. In Pattern
%     Recognition and Information Processing (PRIP) (pp. 140–145).
% -------------------------------------------------------------------------
% INPUTS:
% - GLSZM: Gray-Level Size Zone Matrix.
%
% ** 'GLSZM' should be the output from 'getGLSZM.m' **
% -------------------------------------------------------------------------
% OUTPUTS:
% - textures: Struture specifying the values of different GLSZM texture
%             features as defined below.
% -------------------------------------------------------------------------
% AUTHOR(S): Martin Vallieres <mart.vallieres@gmail.com>
% -------------------------------------------------------------------------
% HISTORY:
% - Creation: January 2013
% - Revision: May 2015
% -------------------------------------------------------------------------
% STATEMENT:
% This file is part of <https://github.com/mvallieres/radiomics/>,
% a package providing MATLAB programming tools for radiomics analysis.
% --> Copyright (C) 2015  Martin Vallieres
%
%    This package is free software: you can redistribute it and/or modify
%    it under the terms of the GNU General Public License as published by
%    the Free Software Foundation, either version 3 of the License, or
%    (at your option) any later version.
%
%    This package is distributed in the hope that it will be useful,
%    but WITHOUT ANY WARRANTY; without even the implied warranty of
%    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%    GNU General Public License for more details.
%
%    You should have received a copy of the GNU General Public License
%    along with this package.  If not, see <http://www.gnu.org/licenses/>.
% -------------------------------------------------------------------------

% USEFUL MATRICES, VECTORS AND QUANTITIES
sz = size(GLSZM); % Size of GLSZM
nRuns = sum(GLSZM(:));
cVect = 1:sz(2); rVect = 1:sz(1);% Row and column vectors
[cMat,rMat] = meshgrid(cVect,rVect); % Column and row indicators for each entry of the GLSZM
pg = sum(GLSZM,2)'; % Gray-Level Run-Number Vector
pr = sum(GLSZM); % Run-Length Run-Number Vector

% COMPUTATION OF TEXTURE FEATURES
% 1. Small Zone Emphasis (SZE), Ref.[1,4]
textures.SZE = (pr*(cVect.^(-2))')/nRuns;

% 2. Large Zone Emphasis (LZE), Ref.[1,4]
textures.LZE = (pr*(cVect.^2)')/nRuns;

% 3. Gray-Level Nonuniformity (GLN), adapted from Ref.[1,4]
textures.GLN = sum(pg.^2)/nRuns;

% 4. Zone-Size Nonuniformity (ZSN), adapted from Ref.[1,4]
textures.ZSN = sum(pr.^2)/nRuns;

% 5. Zone Percentage (ZP), adapted from Ref.[1,4]
textures.ZP = nRuns/(pr*cVect');

% 6. Low Gray-Level Zone Emphasis (LGZE), Ref.[2,4]
textures.LGZE = (pg*(rVect.^(-2))')/nRuns;

% 7. High Gray-Level Zone Emphasis (HGZE), Ref.[2,4]
textures.HGZE = (pg*(rVect.^2)')/nRuns;

% 8. Small Zone Low Gray-Level Emphasis (SZLGE), Ref.[3,4]
textures.SZLGE = sum(sum(GLSZM.*(rMat.^(-2)).*(cMat.^(-2))))/nRuns;

% 9. Small Zone High Gray-Level Emphasis (SZHGE), Ref.[3,4]
textures.SZHGE = sum(sum(GLSZM.*(rMat.^2).*(cMat.^(-2))))/nRuns;

% 10. Large Zone Low Gray-Level Emphasis (LZLGE), Ref.[3,4]
textures.LZLGE = sum(sum(GLSZM.*(rMat.^(-2)).*(cMat.^2)))/nRuns;

% 11. Large Zone High Gray-Level Emphasis (LZHGE), Ref.[3,4]
textures.LZHGE = sum(sum(GLSZM.*(rMat.^2).*(cMat.^2)))/nRuns;

% New features according to Ref.[4]
% GLSZM = GLSZM./nRuns; % In the future, this operation will be applied at the beginning of the function
% pg=sum(GLSZM,2)'; pr=sum(GLSZM);
ug = (pg*rVect')/(sz(1)*sz(2));
ur = (pr*cVect')/(sz(1)*sz(2));

% 12. Gray-Level Variance (GLV), adapted from Ref.[4]
GLV = 0;
for g = 1:sz(1)
    for r = 1:sz(2)
        GLV = GLV + (GLSZM(g,r)*g-ug)^2;
    end
end
textures.GLV = GLV/(sz(1)*sz(2));

% 13. Zone-Size Variance (ZSV), adapted from Ref.[4]
ZSV = 0;
for g = 1:sz(1)
    for r = 1:sz(2)
        ZSV = ZSV + (GLSZM(g,r)*r-ur)^2;
    end
end
textures.ZSV = ZSV/(sz(1)*sz(2));

end

author: mvallieres

repo: radiomics


text:
function A = old_get_gain(fid,ndx)
% GET_GAIN: Get a set of source foward fields from a .bin gain matrix file.
%
% USAGE:  A = old_get_gain(fid,ndx);
%
% INPUT:
%    - fid : a file identifier following a fopen call on a valid binary gain matrix file
%    - ndx : vector array of source indices for which the forward fields are requested
% OUTPUT:
%    - A   : array of requested forward fields
%
% NOTES:
%    - The user is responsible for opening the file in the proper machine format.
%    - For head model matrices, this routine should only be accessed from read_gain,
%      so that a consistent machine format is used to read the file.
%
% SEE ALSO: READ_GAIN, LOAD_RAW, SAVE_RAW

% @=============================================================================
% This function is part of the Brainstorm software:
% https://neuroimage.usc.edu/brainstorm
%
% Copyright (c) University of Southern California & McGill University
% This software is distributed under the terms of the GNU General Public License
% as published by the Free Software Foundation. Further details on the GPLv3
% license can be found at http://www.gnu.org/copyleft/gpl.html.
%
% FOR RESEARCH PURPOSES ONLY. THE SOFTWARE IS PROVIDED "AS IS," AND THE
% UNIVERSITY OF SOUTHERN CALIFORNIA AND ITS COLLABORATORS DO NOT MAKE ANY
% WARRANTY, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO WARRANTIES OF
% MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, NOR DO THEY ASSUME ANY
% LIABILITY OR RESPONSIBILITY FOR THE USE OF THIS SOFTWARE.
%
% For more information type "brainstorm license" at command prompt.
% =============================================================================@
%
% Authors: Sylvain Baillet, October 2002

frewind(fid);

rows = fread(fid,1,'uint32');

if nargin == 2
   cols = length(ndx);

   A = zeros(rows,cols);

   for i = 1:cols,
      % 4 bytes per element, find starting point
      offset = 4 + (ndx(i)-1)*rows*4;
      status = fseek(fid,offset,'bof');
      if(status == -1),
         error('Error reading file at column %.0f',i);
      end

      A(:,i) = fread(fid,[rows,1],'float32');

   end

   return

else % ndx is not specified: read the whole matrix

   offset = 4;
   fseek(fid,offset,'bof');
   fbegin = ftell(fid);
   fseek(fid,0,'eof');
   fend = ftell(fid);

   frewind(fid);
   offset = 4;
   fseek(fid,offset,'bof');

   cols = (fend - fbegin + 1) / (4* rows); % Number of sources / columns

   A = fread(fid,[rows, cols],'float32');

end

author: brainstorm-tools


text:
function test_nint_test03 ( )

%*****************************************************************************80
%
%% TEST03 applies a composite midpoint rule to box regions.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    05 June 2007
%
%  Author:
%
%    John Burkardt
%
  clear

  fprintf ( 1, '\n' );
  fprintf ( 1, 'TEST03\n' );
  fprintf ( 1, '  Use a simple product rule on box regions.\n' );
  fprintf ( 1, '  Use a fixed spatial dimension.\n' );

  problem_num = get_problem_num ( );

  fprintf ( 1, '\n' );
  fprintf ( 1, ...
    '  Prob   Dim  Subs       Approx          Exact          Error\n' );
  fprintf ( 1, '\n' );

  for problem = 1 : problem_num

    dim_num = 3;
%
%  Set problem data to default values.
%
    p00_default ( problem, dim_num );
%
%  Get the region type.
%
    region = p00_region ( problem );

    if ( strncmpi ( region, 'box', 3 ) )

      for sub_num = 1 : 2 : 5

        result = p00_box_gl05 ( problem, dim_num, sub_num );

        exact = p00_exact ( problem, dim_num );

        if ( exact == r8_huge ( ) )

          fprintf ( 1, '  %4d  %4d  %4d  %14f  %s  %s\n', ...
            problem, dim_num, sub_num, result, ...
            '--------------', '--------------' );

        else

          error = abs ( result - exact );

          fprintf ( 1, '  %4d  %4d  %4d  %14f  %14f  %14f\n', ...
            problem, dim_num, sub_num, result, exact, error );

        end

      end

      fprintf ( 1, '\n' );

    end

  end

  return
end


text:
function newval=meshremap(fromval,elemid,elembary,toelem,nodeto)
%
% newval=meshremap(fromval,elemid,elembary,toelem,nodeto)
%
% Redistribute nodal values from the source mesh to the target mesh so that
% the sum of each property on each mesh is the same
%
% author: Qianqian Fang (q.fang at neu.edu)
%
% input:
%	 fromval: values defined at the source mesh nodes, the row or column
%	          number must be the same as the source mesh node number, which
%	          is the same as the elemid length
%	 elemid: the IDs of the target mesh element that encloses the nodes of
%            the source mesh nodes; a vector of length of src mesh node
%            count; elemid and elembary can be generated by calling
%
%           [elemid,elembary]=tsearchn(node_target, elem_target, node_src);
%
%           note that the mapping here is inverse to that in meshinterp()
%
%	 elembary: the bary-centric coordinates of each source mesh nodes
%	         within the target mesh elements, sum of each row is 1, expect
%	         3 or 4 columns (or can be N-D)
%    toelem: the element list of the target mesh
%    nodeto: the total number of target mesh nodes
%
%
% output:
%	 newval: a 2D array with rows equal to the target mesh nodes (nodeto),
%            and columns equals to the value numbers defined at each source
%            mesh node
% example:
%
%    [n1,f1,e1]=meshabox([0 0 0],[10 20 5],1); % src mesh
%    [n2,f2,e2]=meshabox([0 0 0],[10 20 5],2); % target mesh
%    [id, ww]=tsearchn(n2,e2,n1);              % project src to target mesh
%    value_src=n1(:,[2 3 1]);             % create dummy values at src mesh
%    newval=meshremap(value_src,id,ww,e2,size(n2,1)); % map to target
%
% -- this function is part of iso2mesh toolbox (http://iso2mesh.sf.net)
%

if(size(fromval,1)==1)
    fromval=fromval(:);
end

if(size(fromval,2)==length(elemid))
    fromval=fromval.';
end

newval=zeros(nodeto,size(fromval,2));

idx=~isnan(elemid);
fromval=fromval(idx,:);
elembary=elembary(idx,:);
idx=elemid(idx);

nodeval=repmat(fromval,[1,1,size(elembary,2)]).*repmat(permute(elembary,[1,3,2]),[1,size(fromval,2),1]);

for i=1:size(elembary,2)
    [ix,iy]=meshgrid(toelem(idx,i),1:size(fromval,2));
    nval=nodeval(:,:,i).';
    newval=newval + accumarray([ix(:),iy(:)],nval(:), size(newval));
end

author: fangq


text:
function g = ice(varargin)
%ICE Launch Interactive Color Editor.
%
%   OUT = ICE('Name', 'Value', ...) transforms an
%   image's color components based on interactively specified mapping
%   functions. Inputs are Name-Value pairs:
%
%     Name            Value
%     ------------    -----------------------------------------------
%     'image'         An RGB or monochrome input image to be
%                     transformed by interactively specified
%                     mappings.
%     'space'         The color space of the components to be
%                     modified. Possible values are 'rgb', 'cmy',
%                     'hsi', 'hsv', 'ntsc' (or 'yiq'), 'ycbcr'. When
%                     omitted, the RGB color space is assumed.
%     'wait'          If 'on' (the default), OUT is the mapped input
%                     image and ICE returns to the calling function
%                     or workspace when closed. If 'off', OUT is the
%                     handle of the mapped input image and ICE
%                     returns immediately.
%
%   EXAMPLES:
%     ice OR ice('wait', 'off')            % Demo user interface
%     ice('image', f)                      % Map RGB or mono image
%     ice('image', f, 'space', 'hsv')      % Map HSV of RGB image
%     g = ice('image', f)                  % Return mapped image
%     g = ice('image', f, 'wait', 'off');  % Return its handle
%
%   ICE displays one popup menu selectable mapping function at a
%   time. Each image component is mapped by a dedicated curve (e.g.,
%   R, G, or B) and then by an all-component curve (e.g., RGB). Each
%   curve's control points are depicted as circles that can be moved,
%   added, or deleted with a two- or three-button mouse:
%
%     Mouse Button    Editing Operation
%     ------------    -----------------------------------------------
%     Left            Move control point by pressing and dragging.
%     Middle          Add and position a control point by pressing
%                     and dragging. (Optionally Shift-Left)
%     Right           Delete a control point. (Optionally
%                     Control-Left)
%
%   Checkboxes determine how mapping functions are computed, whether
%   the input image and reference pseudo- and full-color bars are
%   mapped, and the displayed reference curve information (e.g.,
%   PDF):
%
%     Checkbox        Function
%     ------------    -----------------------------------------------
%     Smooth          Checked for cubic spline (smooth curve)
%                     interpolation. If unchecked, piecewise linear.
%     Clamp Ends      Checked to force the starting and ending curve
%                     slopes in cubic spline interpolation to 0. No
%                     effect on piecewise linear.
%     Show PDF        Display probability density function(s) [i.e.,
%                     histogram(s)] of the image components affected
%                     by the mapping function.
%     Show CDF        Display cumulative distributions function(s)
%                     instead of PDFs.
%                     <Note: Show PDF/CDF are mutually exclusive.>
%     Map Image       If checked, image mapping is enabled; else
%                     not.
%     Map Bars        If checked, pseudo- and full-color bar mapping
%                     is enabled; else display the unmapped bars (a
%                     gray wedge and hue wedge, respectively).
%
%   Mapping functions can be initialized via pushbuttons:
%
%     Button          Function
%     ------------    -----------------------------------------------
%     Reset           Init the currently displayed mapping function
%                     and uncheck all curve parameters.
%     Reset All       Initialize all mapping functions.
%
%   Copyright 2002-2020 Gatesmark
%
%   This function, and other functions in the DIPUM Toolbox, are based
%   on the theoretical and practical foundations established in the
%   book Digital Image Processing Using MATLAB, 3rd ed., Gatesmark
%   Press, 2020.
%
%   Book website: http://www.imageprocessingplace.com
%   License: https://github.com/dipum/dipum-toolbox/blob/master/LICENSE.txt

parser = inputParser;
addParameter(parser,'image',imread('peppers.png'));
addParameter(parser,'space','rgb');
addParameter(parser,'wait','on');
parse(parser,varargin{:});

app = iceControlPanel(parser.Results.image,parser.Results.space);
output_fig = app.OutputImageFigure;

if strcmp(parser.Results.wait,'on')
   uiwait(app.Figure);
   im_handle = findobj(output_fig,'type','image');
   g = im_handle.CData;
else
   g = app.OutputImageFigure;
end

author: dipum

repo: dipum-toolbox

save_path: github-repos/MATLAB/dipum-dipum-toolbox


text:
% This file is part of the project NILM-Eval (https://github.com/beckel/nilm-eval).
% Licence: GPL 2.0 (http://www.gnu.org/licenses/gpl-2.0.html)
% Copyright: ETH Zurich, 2014
% Author: Romano Cicchetti

function [usedParameterValues] = updateEventMetrics(usedParameterValues, summary)

    % update event metrics in 'usedParametervalues'

    for event_metric = fieldnames(summary.events)'
       if (~ismember(event_metric{1}, usedParameterValues.event_metrics))
            metrics_cell = usedParameterValues.event_metrics;
            metrics_cell{length(metrics_cell)+1} = event_metric{1};
            usedParameterValues.event_metrics = metrics_cell;
       end
    end

end

author: beckel

repo: nilm-eval

save_path: github-repos/MATLAB/beckel-nilm-eval


text:
function sparse_grid_composite_test01 ( dim_min, dim_max, level_max_min, ...
  level_max_max )

%*****************************************************************************80
%
%% SPARSE_GRID_COMPOSITE_TEST01 tests SPARSE_GRID_COMPOSITE_SIZE.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    11 August 2009
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer DIM_MIN, the minimum spatial dimension to consider.
%
%    Input, integer DIM_MAX, the maximum spatial dimension to consider.
%
%    Input, integer LEVEL_MAX_MIN, the minimum value of LEVEL_MAX to consider.
%
%    Input, integer LEVEL_MAX_MAX, the maximum value of LEVEL_MAX to consider.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'SPARSE_GRID_COMPOSITE_TEST01\n' );
  fprintf ( 1, '  SPARSE_GRID_COMPOSITE_SIZE returns the number of distinct\n' );
  fprintf ( 1, '  points in a composite sparse grid.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Each sparse grid is of spatial dimension DIM,\n' );
  fprintf ( 1, '  and is made up of all product grids of levels up to LEVEL_MAX.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '   DIM: ' );

  for dim_num = dim_min : dim_max
    fprintf ( 1, '  %8d', dim_num );
  end

  fprintf ( 1, '\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '   LEVEL_MAX\n' );
  fprintf ( 1, '\n' );

  for level_max = level_max_min : level_max_max
    fprintf ( 1, '    %4d', level_max );
    for dim_num = dim_min : dim_max
      point_num = sparse_grid_composite_size ( dim_num, level_max );
      fprintf ( 1, '  %8d', point_num );
    end
    fprintf ( 1, '\n' );
  end

  return
end


text:
function matches=associate2(first_list, second_list, offset, max_difference)
% first_list first column are timestamps, so are second_list, both
% first_list and second_list are of same columns, offset is to be applied
% to second_list when doing association,
% matches first col the index in the first_list and second col, index in
% second_list
% author Jianzhu Huai 2014

    firstArray=[first_list(:,1), ones(size(first_list,1),1), (1:size(first_list,1))'; inf, 1, -1];
    secondArray=[second_list(:,1)+offset, ones(size(second_list,1),1)*2, (1:size(second_list,1))'; inf, 2, -1];
    resultArray=merge_sorted_arrays(firstArray, secondArray);
    secondUsed=zeros(size(second_list,1),1);% whether second data was matched
    matches=zeros(size(firstArray,1),2);
    matchCount=0;
    for i=1:size(resultArray, 1)
        if(resultArray(i, 2)==1) %from first_list
            % check both sides
            leftdiff=inf;
            rightdiff=inf;

            if(i>1 && resultArray(i-1,2)==2 && secondUsed(resultArray(i-1,3))==0)
                   leftdiff= abs(firstArray(resultArray(i,3), 1)-secondArray(resultArray(i-1,3), 1));
            end
            if(i<size(resultArray,1)&& resultArray(i+1,2)==2 && secondUsed(resultArray(i+1,3))==0)
                   rightdiff= abs(firstArray(resultArray(i,3), 1)-secondArray(resultArray(i+1,3), 1));
            end
            if(leftdiff<max_difference || rightdiff< max_difference)
                if(leftdiff<rightdiff)
                    assert(leftdiff<max_difference);
                    matchCount=matchCount+1;
                    matches(matchCount, :)=[resultArray(i,3), resultArray(i-1,3)];
                    secondUsed(resultArray(i-1,3))= resultArray(i,3);
                else
                    matchCount=matchCount+1;
                    matches(matchCount, :)=[resultArray(i,3), resultArray(i+1, 3)];
                    assert(rightdiff<max_difference);
                    secondUsed(resultArray(i+1,3))= resultArray(i,3);
                end
            end
        end
    end
    matches=matches(1:matchCount,:);

    function A=merge_sorted_arrays(L, R)
        % L mxn, R pxn, ascendingly sorted arrays,
        % n fields, m-1 and p-1 observations, last row is [inf, ...] as sentinel
        % see introduction to algorithms 2nd edition p29
        r=size(L,1)+size(R,1)-2;
        iota=1;
        j=1;
        A=zeros(r,size(L,2));
        for k=1:r
            if(L(iota,1)<=R(j,1))
                A(k,:)=L(iota,:);
                iota=iota+1;
            else
                A(k, :)=R(j,:);
                j=j+1;
            end
        end
    end

end

repo: ekfmonoslam


text:
function  x = SoftTh(y,thld)
%    x = SoftTh(y,thld);
%
%    SOFTTH soft thresholds the input signal y with the threshold value
%    thld.
%
%    Input:
%       y    : 1D or 2D signal to be thresholded
%       thld : Threshold value
%
%    Output:
%       x : Soft thresholded output (x = sign(y)(|y|-thld)_+)
%
%  HERE'S AN EASY WAY TO RUN THE EXAMPLES:
%  Cut-and-paste the example you want to run to a new file
%  called ex.m, for example. Delete out the % at the beginning
%  of each line in ex.m (Can use search-and-replace in your editor
%  to replace it with a space). Type 'ex' in matlab and hit return.
%
%
%    Example:
%       y = makesig('Doppler',8);
%       thld = 0.2;
%       x = SoftTh(y,thld)
%       x = 0 0 0 -0.0703 0 0.2001 0.0483 0
%
%    See also: HardTh
%
%    Reference:
%       "De-noising via Soft-Thresholding" Tech. Rept. Statistics,
%       Stanford, 1992. D.L. Donoho.
%

%File Name: SoftTh.m
%Last Modification Date: 8/15/95	17:49:48
%Current Version: SoftTh.m	2.4
%File Creation Date: Mon Mar  7 10:38:45 1994
%Author: Haitao Guo  <harry@jazz.rice.edu>
%
%Copyright (c) 2000 RICE UNIVERSITY. All rights reserved.
%Created by Haitao Guo, Department of ECE, Rice University.
%
%This software is distributed and licensed to you on a non-exclusive
%basis, free-of-charge. Redistribution and use in source and binary forms,
%with or without modification, are permitted provided that the following
%conditions are met:
%
%1. Redistribution of source code must retain the above copyright notice,
%   this list of conditions and the following disclaimer.
%2. Redistribution 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.
%3. All advertising materials mentioning features or use of this software
%   must display the following acknowledgment: This product includes
%   software developed by Rice University, Houston, Texas and its contributors.
%4. Neither the name of the University nor the names of its contributors
%   may be used to endorse or promote products derived from this software
%   without specific prior written permission.
%
%THIS SOFTWARE IS PROVIDED BY WILLIAM MARSH RICE UNIVERSITY, HOUSTON, TEXAS,
%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 RICE UNIVERSITY
%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 INTERRUPTIONS) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
%WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
%OTHERWISE), PRODUCT LIABILITY, OR OTHERWISE ARISING IN ANY WAY OUT OF THE
%USE OF THIS SOFTWARE,  EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
%
%For information on commercial licenses, contact Rice University's Office of
%Technology Transfer at techtran@rice.edu or (713) 348-6173
%
%Change History:
%

x = abs(y);
x = sign(y).*(x >= thld).*(x - thld);

author: thomaskuestner


text:
function[imdb] = setupImdbDetection(trainName, testName, net)
% [imdb] = setupImdbDetection(trainName, testName, net)

global DATAopts;

%%% Setup the Imdb
% Get and test images
trainIms = textread(sprintf(DATAopts.imgsetpath, trainName), '%s'); %#ok<DTXTRD>
testIms = textread(sprintf(DATAopts.imgsetpath, testName), '%s'); %#ok<DTXTRD>

% Make train, val, and test set. For Pascal, I illegally use part of the test images
% as validation set. This is to match Girshick performance while still having
% meaningful graphs for the validation set.
% Note: allIms are just all images. datasetIdx determines how these are divided over
% train, val, and test.
allIms = cat(1, trainIms, testIms);
datasetIdx = cell(3, 1);
datasetIdx{1} = (1:length(trainIms))';  % Jasper: Use all training images. Only for comparison Pascal Girshick
datasetIdx{2} = (length(trainIms)+1:length(trainIms)+501)'; % Use part of the test images for validation. Not entirely legal, but otherwise it will take much longer to get where we want.
datasetIdx{3} = (length(trainIms)+1:length(allIms))';

imdb = ImdbDetectionFullSupervision(DATAopts.imgpath(1:end-6), ...        % path
    DATAopts.imgpath(end-3:end), ...      % image extension
    DATAopts.gStructPath, ...             % gStruct path
    allIms, ...                           % all images
    datasetIdx, ...                       % division into train/val/test
    net.meta.normalization.averageImage);      % average image used to pretrain network

% Usually instance weighting gives better performance. But not Girshick style
% imdbPascal.SetInstanceWeighting(true);

% Store lists for use in eval
imdb.misc.trainIms  = trainIms;
imdb.misc.testIms   = testIms;

author: nightrome

repo: matconvnet-calvin


text:
% labelmePrepareData
DO_DB = 0;
DO_LABELS = 0;
DO_TRANSFER_BMAPS = 0;
GET_TRAINING_DATA = 1;

addpath('~/src/labelme');

set = 'train'; % train, test

basedir = '~/data/occlusion/labelme';
imdir = fullfile(basedir, set, 'Images');
anndir = fullfile(basedir, set, 'Annotations');
labeldir = fullfile(basedir, 'labels2');

if DO_DB
    outdir = basedir;
    db = LMdatabase(anndir);
    save(fullfile(outdir, [set '_db']), 'db');
end

if DO_LABELS
    filestr = '*.xml';
    maxsize = 1024;
    outdir = fullfile(basedir, 'labels3');
    ignoreparts = true;
    processDirectory(anndir, filestr, outdir, '_labels.mat', ...
        @processIm2LabelmeLabels, imdir, maxsize, ignoreparts);
%     load(fullfile(outdir, [set '_db']));
%     parfor k = 1:numel(db)
%         [lim{k}, objnames{k}] = annotation2labels(db(k).annotation, imdir, 1024);
%     end
%     save(fullfile(outdir, [set '_labels']), 'lim', 'objnames');
end

if DO_TRANSFER_BMAPS
    setOcclusionDirectories;
    %load(fullfile(basedir, [set '_db']));
    isFullyLabeled = false(size(db));

    for k = 1:numel(db)

        if mod(k, 100)==0
            disp(num2str(k));
        end
        try
            folder = db(k).annotation.folder;
            bn = strtok(db(k).annotation.filename, '.');
            occfn = fullfile(occdir, folder, [bn '_occlusion']);
            labelfn = fullfile(labeldir, folder, [bn '_labels']);
            load(occfn, 'bndinfo_all');
            load(labelfn, 'lim');

            if mean(lim(:)>0)>0.97 % fully labeled

                isFullyLabeled(k) = true;

                % resize object segmentation image to be same size as
                [imh, imw] = size(bndinfo_all{1}.wseg);
                lim = imresize(lim, [imh imw], 'nearest');

                % get distance of pixels from ground truth boundary
                gtmap = seg2bmap(lim, imw, imh);
                gtdist = bwdist(gtmap);
                maxdist = 0.01*sqrt(imh.^2+imw.^2);

                ne = bndinfo_all{1}.ne;
                labels{k} = false(ne, 1);
                labels_near{k} = false(ne, 1);

                [lab, labim, err] = transferRegionLabels(lim, double(bndinfo_all{1}.wseg));
                for k2 = 1:ne
                    LR = bndinfo_all{1}.edges.spLR(k2, :);
                    labels{k}(k2) = lab(LR(1))~=lab(LR(2));
                    labels_near{k}(k2) = labels{k}(k2) && ...
                        (mean(gtdist(bndinfo_all{1}.edges.indices{k2}))<maxdist);
                end
            end
        catch
            disp(lasterr);
        end
    end
    save(fullfile(traindir, 'tranferredBoundaryLabels.mat'), 'labels', 'labels_near', 'isFullyLabeled');

end

if GET_TRAINING_DATA

    %N = 1000; % number of pixels per image to sample

    %load(fullfile(basedir, [set '_db']));
    load(fullfile(traindir, 'tranferredBoundaryLabels.mat'));
    x = cell(numel(labels_near), 1);
    y = cell(numel(labels_near), 1);
    for k = 1:numel(labels_near)

        if mod(k, 100)==0
            disp(num2str(k));
        end

        if isFullyLabeled(k)

            folder = db(k).annotation.folder;
            bn = strtok(db(k).annotation.filename, '.');

            occfn = fullfile(occdir, folder, [bn '_occlusion']);
            load(occfn, 'bndinfo_all');
            y{k} = labels_near{k};
            y{k} = 2*y{k}-1; % convert to {-1, 1}

            % set label of pixels near border to 0
            brd = 10;
            ind  = getBoundaryCenterIndices(bndinfo_all{1});
            [yk, xk] = ind2sub(ind, occ.bndinfo.imsize(1:2));
            valid = (yk>brd) & (yk <= occ.bndinfo_all{1}.imsize(1)-brd) & ...
                (xk > brd) & (xk <=occ.bndinfo_all{1}.imsize(2)-brd);
            y{k}(~valid) = 0;

            occ.po_all = getOcclusionMaps(bndinfo_all);
            x{k} = getBoundaryEdgletFeatures(bndinfo_all, occ);
        end

    end
    x = cat(1, x{:});
    y = cat(1, y{:});
    featDescription = {'pocc1', 'pocc2', 'pocc3', 'pocc4', 'pocc_max'};
    save(fullfile(traindir, 'trainData.mat'), 'x', 'y', 'featDescription');
end

author: Cloud-CV

repo: object-proposals


text:
%-------------------------------------------------------------------------------------------------------------------%
%
% IB2d is an Immersed Boundary Code (IB) for solving fully coupled non-linear
% 	fluid-structure interaction models. This version of the code is based off of
%	Peskin's Immersed Boundary Method Paper in Acta Numerica, 2002.
%
% Author: Nicholas A. Battista
% Email:  battistn@tcnj.edu
% Date Modified: April 2021
% Current Institution: TCNJ
%
% IB2d Date Created: May 27th, 2015
% Institution Created: UNC-CH
%
% This code is capable of creating Lagrangian Structures using:
% 	1. Springs
% 	2. Beams (torsional springs or non-invariant beams)
% 	3. Target Points
%	4. Muscle-Model (combined Force-Length-Velocity model, "HIll+(Length-Tension)")
%
% One is able to update those Lagrangian Structure parameters, e.g., spring constants, resting lengths, etc
%
% There are a number of built in Examples, mostly used for teaching purposes.
%
% If you would like us to add a specific muscle model, please let Nick (battistn@tcnj.edu) know.
%
%--------------------------------------------------------------------------------------------------------------------%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% FUNCTION: updates the target point positions
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function targets = update_Target_Point_Positions(dt,current_time,targets)

%IDs = targets(:,1);                 % Stores Lag-Pt IDs in col vector
%xPts= targets(:,2);                % Original x-Values of x-Target Pts.
%yPts= targets(:,3);                 % Original y-Values of y-Target Pts.
%kStiffs = targets(:,4);             % Stores Target Stiffnesses

%-----------------------------------------------------
% Geometric Translations to Quadrant 1
%    since stored .mat geo data centered at (0,0)
%-----------------------------------------------------
xoffset = 1;   % To put geometry into QUADRANT-1
yoffset = 0.5; % To put geometry into QUADRANT-1

%-----------------------------------------------
% Note: (1) ds for stored data is 0.6/(2*1024)
%       (2) 'ratio' is comparing 1024:desired resolution
%-----------------------------------------------
Lx=2;          % Horizontal Eulerian grid length
Nx=64;        % Eulerian grid resolution
dx=Lx/Nx;      % Eulerian grid spacing
ds=dx/2;       % Lagrangian point spacing

%-----------------------------------------------
% Load prescribed position data for tentacles
%-----------------------------------------------
load('total_coeffs.mat')
load('coral_coeff_30.mat')

%-----------------------------------------------
% Arclength
%-----------------------------------------------
s=(0:ds:total_meanL)/total_meanL;

%---------------------------------------------------------
% Get index correponding to current time in
%   simulation to determine how interpolation occurs
%---------------------------------------------------------
indx=ceil(current_time/dt)+1;

%-----------------------------------------------
% Load geometry state data
%-----------------------------------------------
load('cval.mat')

%-------------------------------------------------
% Get interpolation polynomial coefficients and
%   then determine geometry of LEFT tentacle
%-------------------------------------------------
C1=c1vals(indx,:);
C2=c2vals(indx,:);

XbL_1=(C1(1)*s.^3+C1(2)*s.^2+C1(3)*s+C1(4));
YbL_1=(C2(1)*s.^3+C2(2)*s.^2+C2(3)*s+C2(4));
L_ten=sum(sqrt((XbL_1(2:end)-XbL_1(1:end-1)).^2 +(YbL_1(2:end)-YbL_1(1:end-1)).^2 ));

XbL=(C1(1)*s.^3+C1(2)*s.^2+C1(3)*s+C1(4))*total_meanL/L_ten+total_offset;
YbL=(C2(1)*s.^3+C2(2)*s.^2+C2(3)*s+C2(4))*total_meanL/L_ten;

%-------------------------------------------------------
% Set up symmetric RIGHT tentacle interpolation states
%-------------------------------------------------------
XbR=-XbL;
YbR=YbL;

%-------------------------------------------------------
% Get coral polyp STEM geometry
%-------------------------------------------------------
YbStem=(YbL(1))+0*XbStem;

%-------------------------------------------------------
% Combine geometry into one vector and translate
%-------------------------------------------------------
x=[flip(XbR) XbStem XbL];
y=[flip(YbR) YbStem YbL];
%
x = x+xoffset;
y = y+yoffset;

%-------------------------------------------------------
% Update target point positions!
%-------------------------------------------------------
targets(:,2) = x;   % Store new xVals
targets(:,3) = y;   % Store new yVals

author: nickabattista