text:
function [imageData, alpha] = export_fig(varargin)
%EXPORT_FIG  Exports figures in a publication-quality format
%
% Examples:
%   imageData = export_fig
%   [imageData, alpha] = export_fig
%   export_fig filename
%   export_fig filename -format1 -format2
%   export_fig ... -nocrop
%   export_fig ... -transparent
%   export_fig ... -native
%   export_fig ... -m<val>
%   export_fig ... -r<val>
%   export_fig ... -a<val>
%   export_fig ... -q<val>
%   export_fig ... -p<val>
%   export_fig ... -d<gs_option>
%   export_fig ... -depsc
%   export_fig ... -<renderer>
%   export_fig ... -<colorspace>
%   export_fig ... -append
%   export_fig ... -bookmark
%   export_fig ... -clipboard
%   export_fig ... -update
%   export_fig ... -nofontswap
%   export_fig(..., handle)
%
% This function saves a figure or single axes to one or more vector and/or
% bitmap file formats, and/or outputs a rasterized version to the workspace,
% with the following properties:
%   - Figure/axes reproduced as it appears on screen
%   - Cropped borders (optional)
%   - Embedded fonts (vector formats)
%   - Improved line and grid line styles
%   - Anti-aliased graphics (bitmap formats)
%   - Render images at native resolution (optional for bitmap formats)
%   - Transparent background supported (pdf, eps, png)
%   - Semi-transparent patch objects supported (png only)
%   - RGB, CMYK or grayscale output (CMYK only with pdf, eps, tiff)
%   - Variable image compression, including lossless (pdf, eps, jpg)
%   - Optionally append to file (pdf, tiff)
%   - Vector formats: pdf, eps
%   - Bitmap formats: png, tiff, jpg, bmp, export to workspace
%
% This function is especially suited to exporting figures for use in
% publications and presentations, because of the high quality and
% portability of media produced.
%
% Note that the background color and figure dimensions are reproduced
% (the latter approximately, and ignoring cropping & magnification) in the
% output file. For transparent background (and semi-transparent patch
% objects), use the -transparent option or set the figure 'Color' property
% to 'none'. To make axes transparent set the axes 'Color' property to
% 'none'. PDF, EPS and PNG are the only formats that support a transparent
% background, while only PNG format supports transparency of patch objects.
%
% The choice of renderer (opengl, zbuffer or painters) has a large impact
% on the quality of output. The default value (opengl for bitmaps, painters
% for vector formats) generally gives good results, but if you aren't
% satisfied then try another renderer.  Notes: 1) For vector formats (EPS,
% PDF), only painters generates vector graphics. 2) For bitmaps, only
% opengl can render transparent patch objects correctly. 3) For bitmaps,
% only painters will correctly scale line dash and dot lengths when
% magnifying or anti-aliasing. 4) Fonts may be substitued with Courier when
% using painters.
%
% When exporting to vector format (PDF & EPS) and bitmap format using the
% painters renderer, this function requires that ghostscript is installed
% on your system. You can download this from:
%   http://www.ghostscript.com
% When exporting to eps it additionally requires pdftops, from the Xpdf
% suite of functions. You can download this from:
%   http://www.foolabs.com/xpdf
%
% Inputs:
%   filename - string containing the name (optionally including full or
%              relative path) of the file the figure is to be saved as. If
%              a path is not specified, the figure is saved in the current
%              directory. If no name and no output arguments are specified,
%              the default name, 'export_fig_out', is used. If neither a
%              file extension nor a format are specified, a ".png" is added
%              and the figure saved in that format.
%   -format1, -format2, etc. - strings containing the extensions of the
%                              file formats the figure is to be saved as.
%                              Valid options are: '-pdf', '-eps', '-png',
%                              '-tif', '-jpg' and '-bmp'. All combinations
%                              of formats are valid.
%   -nocrop - option indicating that the borders of the output are not to
%             be cropped.
%   -transparent - option indicating that the figure background is to be
%                  made transparent (png, pdf and eps output only).
%   -m<val> - option where val indicates the factor to magnify the
%             on-screen figure pixel dimensions by when generating bitmap
%             outputs (does not affect vector formats). Default: '-m1'.
%   -r<val> - option val indicates the resolution (in pixels per inch) to
%             export bitmap and vector outputs at, keeping the dimensions
%             of the on-screen figure. Default: '-r864' (for vector output
%             only). Note that the -m option overides the -r option for
%             bitmap outputs only.
%   -native - option indicating that the output resolution (when outputting
%             a bitmap format) should be such that the vertical resolution
%             of the first suitable image found in the figure is at the
%             native resolution of that image. To specify a particular
%             image to use, give it the tag 'export_fig_native'. Notes:
%             This overrides any value set with the -m and -r options. It
%             also assumes that the image is displayed front-to-parallel
%             with the screen. The output resolution is approximate and
%             should not be relied upon. Anti-aliasing can have adverse
%             effects on image quality (disable with the -a1 option).
%   -a1, -a2, -a3, -a4 - option indicating the amount of anti-aliasing to
%                        use for bitmap outputs. '-a1' means no anti-
%                        aliasing; '-a4' is the maximum amount (default).
%   -<renderer> - option to force a particular renderer (painters, opengl or
%                 zbuffer). Default value: opengl for bitmap formats or
%                 figures with patches and/or transparent annotations;
%                 painters for vector formats without patches/transparencies.
%   -<colorspace> - option indicating which colorspace color figures should
%                   be saved in: RGB (default), CMYK or gray. CMYK is only
%                   supported in pdf, eps and tiff output.
%   -q<val> - option to vary bitmap image quality (in pdf, eps and jpg
%             files only).  Larger val, in the range 0-100, gives higher
%             quality/lower compression. val > 100 gives lossless
%             compression. Default: '-q95' for jpg, ghostscript prepress
%             default for pdf & eps. Note: lossless compression can
%             sometimes give a smaller file size than the default lossy
%             compression, depending on the type of images.
%   -p<val> - option to pad a border of width val to exported files, where
%             val is either a relative size with respect to cropped image
%             size (i.e. p=0.01 adds a 1% border). For EPS & PDF formats,
%             val can also be integer in units of 1/72" points (abs(val)>1).
%             val can be positive (padding) or negative (extra cropping).
%             If used, the -nocrop flag will be ignored, i.e. the image will
%             always be cropped and then padded. Default: 0 (i.e. no padding).
%   -append - option indicating that if the file (pdfs only) already
%             exists, the figure is to be appended as a new page, instead
%             of being overwritten (default).
%   -bookmark - option to indicate that a bookmark with the name of the
%               figure is to be created in the output file (pdf only).
%   -clipboard - option to save output as an image on the system clipboard.
%                Note: background transparency is not preserved in clipboard
%   -d<gs_option> - option to indicate a ghostscript setting. For example,
%                   -dMaxBitmap=0 or -dNoOutputFonts (Ghostscript 9.15+).
%   -depsc -  option to use EPS level-3 rather than the default level-2 print
%             device. This solves some bugs with Matlab's default -depsc2 device
%             such as discolored subplot lines on images (vector formats only).
%   -update - option to download and install the latest version of export_fig
%   -nofontswap - option to avoid font swapping. Font swapping is automatically
%             done in vector formats (only): 11 standard Matlab fonts are
%             replaced by the original figure fonts. This option prevents this.
%   handle -  The handle of the figure, axes or uipanels (can be an array of
%             handles, but the objects must be in the same figure) to be
%             saved. Default: gcf.
%
% Outputs:
%   imageData - MxNxC uint8 image array of the exported image.
%   alpha     - MxN single array of alphamatte values in the range [0,1],
%               for the case when the background is transparent.
%
%   Some helpful examples and tips can be found at:
%      https://github.com/altmany/export_fig
%
%   See also PRINT, SAVEAS, ScreenCapture (on the Matlab File Exchange)

%{
% Copyright (C) Oliver Woodford 2008-2014, Yair Altman 2015-

% The idea of using ghostscript is inspired by Peder Axensten's SAVEFIG
% (fex id: 10889) which is itself inspired by EPS2PDF (fex id: 5782).
% The idea for using pdftops came from the MATLAB newsgroup (id: 168171).
% The idea of editing the EPS file to change line styles comes from Jiro
% Doke's FIXPSLINESTYLE (fex id: 17928).
% The idea of changing dash length with line width came from comments on
% fex id: 5743, but the implementation is mine :)
% The idea of anti-aliasing bitmaps came from Anders Brun's MYAA (fex id:
% 20979).
% The idea of appending figures in pdfs came from Matt C in comments on the
% FEX (id: 23629)

% Thanks to Roland Martin for pointing out the colour MATLAB
% bug/feature with colorbar axes and transparent backgrounds.
% Thanks also to Andrew Matthews for describing a bug to do with the figure
% size changing in -nodisplay mode. I couldn't reproduce it, but included a
% fix anyway.
% Thanks to Tammy Threadgill for reporting a bug where an axes is not
% isolated from gui objects.
%}
%{
% 23/02/12: Ensure that axes limits don't change during printing
% 14/03/12: Fix bug in fixing the axes limits (thanks to Tobias Lamour for reporting it).
% 02/05/12: Incorporate patch of Petr Nechaev (many thanks), enabling bookmarking of figures in pdf files.
% 09/05/12: Incorporate patch of Arcelia Arrieta (many thanks), to keep tick marks fixed.
% 12/12/12: Add support for isolating uipanels. Thanks to michael for suggesting it.
% 25/09/13: Add support for changing resolution in vector formats. Thanks to Jan Jaap Meijer for suggesting it.
% 07/05/14: Add support for '~' at start of path. Thanks to Sally Warner for suggesting it.
% 24/02/15: Fix Matlab R2014b bug (issue #34): plot markers are not displayed when ZLimMode='manual'
% 25/02/15: Fix issue #4 (using HG2 on R2014a and earlier)
% 25/02/15: Fix issue #21 (bold TeX axes labels/titles in R2014b)
% 26/02/15: If temp dir is not writable, use the user-specified folder for temporary EPS/PDF files (Javier Paredes)
% 27/02/15: Modified repository URL from github.com/ojwoodford to /altmany
%           Indented main function
%           Added top-level try-catch block to display useful workarounds
% 28/02/15: Enable users to specify optional ghostscript options (issue #36)
% 06/03/15: Improved image padding & cropping thanks to Oscar Hartogensis
% 26/03/15: Fixed issue #49 (bug with transparent grayscale images); fixed out-of-memory issue
% 26/03/15: Fixed issue #42: non-normalized annotations on HG1
% 26/03/15: Fixed issue #46: Ghostscript crash if figure units <> pixels
% 27/03/15: Fixed issue #39: bad export of transparent annotations/patches
% 28/03/15: Fixed issue #50: error on some Matlab versions with the fix for issue #42
% 29/03/15: Fixed issue #33: bugs in Matlab's print() function with -cmyk
% 29/03/15: Improved processing of input args (accept space between param name & value, related to issue #51)
% 30/03/15: When exporting *.fig files, then saveas *.fig if figure is open, otherwise export the specified fig file
% 30/03/15: Fixed edge case bug introduced yesterday (commit #ae1755bd2e11dc4e99b95a7681f6e211b3fa9358)
% 09/04/15: Consolidated header comment sections; initialize output vars only if requested (nargout>0)
% 14/04/15: Workaround for issue #45: lines in image subplots are exported in invalid color
% 15/04/15: Fixed edge-case in parsing input parameters; fixed help section to show the -depsc option (issue #45)
% 21/04/15: Bug fix: Ghostscript croaks on % chars in output PDF file (reported by Sven on FEX page, 15-Jul-2014)
% 22/04/15: Bug fix: Pdftops croaks on relative paths (reported by Tintin Milou on FEX page, 19-Jan-2015)
% 04/05/15: Merged fix #63 (Kevin Mattheus Moerman): prevent tick-label changes during export
% 07/05/15: Partial fix for issue #65: PDF export used painters rather than opengl renderer (thanks Nguyenr)
% 08/05/15: Fixed issue #65: bad PDF append since commit #e9f3cdf 21/04/15 (thanks Robert Nguyen)
% 12/05/15: Fixed issue #67: exponent labels cropped in export, since fix #63 (04/05/15)
% 28/05/15: Fixed issue #69: set non-bold label font only if the string contains symbols (\beta etc.), followup to issue #21
% 29/05/15: Added informative error message in case user requested SVG output (issue #72)
% 09/06/15: Fixed issue #58: -transparent removed anti-aliasing when exporting to PNG
% 19/06/15: Added -update option to download and install the latest version of export_fig
% 07/07/15: Added -nofontswap option to avoid font-swapping in EPS/PDF
% 16/07/15: Fixed problem with anti-aliasing on old Matlab releases
%}

    if nargout
        [imageData, alpha] = deal([]);
    end
    hadError = false;
    displaySuggestedWorkarounds = true;

    % Ensure the figure is rendered correctly _now_ so that properties like axes limits are up-to-date
    drawnow;
    pause(0.05);  % this solves timing issues with Java Swing's EDT (http://undocumentedmatlab.com/blog/solving-a-matlab-hang-problem)

    % Parse the input arguments
    fig = get(0, 'CurrentFigure');
    [fig, options] = parse_args(nargout, fig, varargin{:});

    % Ensure that we have a figure handle
    if isequal(fig,-1)
        return;  % silent bail-out
    elseif isempty(fig)
        error('No figure found');
    end

    % Isolate the subplot, if it is one
    cls = all(ismember(get(fig, 'Type'), {'axes', 'uipanel'}));
    if cls
        % Given handles of one or more axes, so isolate them from the rest
        fig = isolate_axes(fig);
    else
        % Check we have a figure
        if ~isequal(get(fig, 'Type'), 'figure');
            error('Handle must be that of a figure, axes or uipanel');
        end
        % Get the old InvertHardcopy mode
        old_mode = get(fig, 'InvertHardcopy');
    end

    % Hack the font units where necessary (due to a font rendering bug in print?).
    % This may not work perfectly in all cases.
    % Also it can change the figure layout if reverted, so use a copy.
    magnify = options.magnify * options.aa_factor;
    if isbitmap(options) && magnify ~= 1
        fontu = findall(fig, 'FontUnits', 'normalized');
        if ~isempty(fontu)
            % Some normalized font units found
            if ~cls
                fig = copyfig(fig);
                set(fig, 'Visible', 'off');
                fontu = findall(fig, 'FontUnits', 'normalized');
                cls = true;
            end
            set(fontu, 'FontUnits', 'points');
        end
    end

    try
        % MATLAB "feature": axes limits and tick marks can change when printing
        Hlims = findall(fig, 'Type', 'axes');
        if ~cls
            % Record the old axes limit and tick modes
            Xlims = make_cell(get(Hlims, 'XLimMode'));
            Ylims = make_cell(get(Hlims, 'YLimMode'));
            Zlims = make_cell(get(Hlims, 'ZLimMode'));
            Xtick = make_cell(get(Hlims, 'XTickMode'));
            Ytick = make_cell(get(Hlims, 'YTickMode'));
            Ztick = make_cell(get(Hlims, 'ZTickMode'));
            Xlabel = make_cell(get(Hlims, 'XTickLabelMode'));
            Ylabel = make_cell(get(Hlims, 'YTickLabelMode'));
            Zlabel = make_cell(get(Hlims, 'ZTickLabelMode'));
        end

        % Set all axes limit and tick modes to manual, so the limits and ticks can't change
        % Fix Matlab R2014b bug (issue #34): plot markers are not displayed when ZLimMode='manual'
        set(Hlims, 'XLimMode', 'manual', 'YLimMode', 'manual');
        set_tick_mode(Hlims, 'X');
        set_tick_mode(Hlims, 'Y');
        if ~using_hg2(fig)
            set(Hlims,'ZLimMode', 'manual');
            set_tick_mode(Hlims, 'Z');
        end
    catch
        % ignore - fix issue #4 (using HG2 on R2014a and earlier)
    end

    % Fix issue #21 (bold TeX axes labels/titles in R2014b when exporting to EPS/PDF)
    try
        if using_hg2(fig) && isvector(options)
            % Set the FontWeight of axes labels/titles to 'normal'
            % Fix issue #69: set non-bold font only if the string contains symbols (\beta etc.)
            texLabels = findall(fig, 'type','text', 'FontWeight','bold');
            symbolIdx = ~cellfun('isempty',strfind({texLabels.String},'\'));
            set(texLabels(symbolIdx), 'FontWeight','normal');
        end
    catch
        % ignore
    end

    % Fix issue #42: non-normalized annotations on HG1 (internal Matlab bug)
    annotationHandles = [];
    try
        if ~using_hg2(fig)
            annotationHandles = findall(fig,'Type','hggroup','-and','-property','Units','-and','-not','Units','norm');
            originalUnits = get(annotationHandles,'Units');
            set(annotationHandles,'Units','norm');
        end
    catch
        % should never happen, but ignore in any case - issue #50
    end

    % Fix issue #46: Ghostscript crash if figure units <> pixels
    oldFigUnits = get(fig,'Units');
    set(fig,'Units','pixels');

    % Set to print exactly what is there
    set(fig, 'InvertHardcopy', 'off');
    % Set the renderer
    switch options.renderer
        case 1
            renderer = '-opengl';
        case 2
            renderer = '-zbuffer';
        case 3
            renderer = '-painters';
        otherwise
            renderer = '-opengl'; % Default for bitmaps
    end

    try
        % Do the bitmap formats first
        if isbitmap(options)
            if abs(options.bb_padding) > 1
                displaySuggestedWorkarounds = false;
                error('For bitmap output (png,jpg,tif,bmp) the padding value (-p) must be between -1<p<1')
            end
            % Get the background colour
            if options.transparent && (options.png || options.alpha)
                % Get out an alpha channel
                % MATLAB "feature": black colorbar axes can change to white and vice versa!
                hCB = findall(fig, 'Type','axes', 'Tag','Colorbar');
                if isempty(hCB)
                    yCol = [];
                    xCol = [];
                else
                    yCol = get(hCB, 'YColor');
                    xCol = get(hCB, 'XColor');
                    if iscell(yCol)
                        yCol = cell2mat(yCol);
                        xCol = cell2mat(xCol);
                    end
                    yCol = sum(yCol, 2);
                    xCol = sum(xCol, 2);
                end
                % MATLAB "feature": apparently figure size can change when changing
                % colour in -nodisplay mode
                pos = get(fig, 'Position');
                % Set the background colour to black, and set size in case it was
                % changed internally
                tcol = get(fig, 'Color');
                set(fig, 'Color', 'k', 'Position', pos);
                % Correct the colorbar axes colours
                set(hCB(yCol==0), 'YColor', [0 0 0]);
                set(hCB(xCol==0), 'XColor', [0 0 0]);

                % The following code might cause out-of-memory errors
                try
                    % Print large version to array
                    B = print2array(fig, magnify, renderer);
                    % Downscale the image
                    B = downsize(single(B), options.aa_factor);
                catch
                    % This is more conservative in memory, but kills transparency (issue #58)
                    B = single(print2array(fig, magnify/options.aa_factor, renderer));
                end

                % Set background to white (and set size)
                set(fig, 'Color', 'w', 'Position', pos);
                % Correct the colorbar axes colours
                set(hCB(yCol==3), 'YColor', [1 1 1]);
                set(hCB(xCol==3), 'XColor', [1 1 1]);

                % The following code might cause out-of-memory errors
                try
                    % Print large version to array
                    A = print2array(fig, magnify, renderer);
                    % Downscale the image
                    A = downsize(single(A), options.aa_factor);
                catch
                    % This is more conservative in memory, but kills transparency (issue #58)
                    A = single(print2array(fig, magnify/options.aa_factor, renderer));
                end

                % Set the background colour (and size) back to normal
                set(fig, 'Color', tcol, 'Position', pos);
                % Compute the alpha map
                alpha = round(sum(B - A, 3)) / (255 * 3) + 1;
                A = alpha;
                A(A==0) = 1;
                A = B ./ A(:,:,[1 1 1]);
                clear B
                % Convert to greyscale
                if options.colourspace == 2
                    A = rgb2grey(A);
                end
                A = uint8(A);
                % Crop the background
                if options.crop
                    %[alpha, v] = crop_borders(alpha, 0, 1);
                    %A = A(v(1):v(2),v(3):v(4),:);
                    [alpha, vA, vB] = crop_borders(alpha, 0, options.bb_padding);
                    if ~any(isnan(vB)) % positive padding
                        B = repmat(uint8(zeros(1,1,size(A,3))),size(alpha));
                        B(vB(1):vB(2), vB(3):vB(4), :) = A(vA(1):vA(2), vA(3):vA(4), :); % ADDED BY OH
                        A = B;
                    else  % negative padding
                        A = A(vA(1):vA(2), vA(3):vA(4), :);
                    end
                end
                if options.png
                    % Compute the resolution
                    res = options.magnify * get(0, 'ScreenPixelsPerInch') / 25.4e-3;
                    % Save the png
                    imwrite(A, [options.name '.png'], 'Alpha', double(alpha), 'ResolutionUnit', 'meter', 'XResolution', res, 'YResolution', res);
                    % Clear the png bit
                    options.png = false;
                end
                % Return only one channel for greyscale
                if isbitmap(options)
                    A = check_greyscale(A);
                end
                if options.alpha
                    % Store the image
                    imageData = A;
                    % Clear the alpha bit
                    options.alpha = false;
                end
                % Get the non-alpha image
                if isbitmap(options)
                    alph = alpha(:,:,ones(1, size(A, 3)));
                    A = uint8(single(A) .* alph + 255 * (1 - alph));
                    clear alph
                end
                if options.im
                    % Store the new image
                    imageData = A;
                end
            else
                % Print large version to array
                if options.transparent
                    % MATLAB "feature": apparently figure size can change when changing
                    % colour in -nodisplay mode
                    pos = get(fig, 'Position');
                    tcol = get(fig, 'Color');
                    set(fig, 'Color', 'w', 'Position', pos);
                    A = print2array(fig, magnify, renderer);
                    set(fig, 'Color', tcol, 'Position', pos);
                    tcol = 255;
                else
                    [A, tcol] = print2array(fig, magnify, renderer);
                end
                % Crop the background
                if options.crop
                    A = crop_borders(A, tcol, options.bb_padding);
                end
                % Downscale the image
                A = downsize(A, options.aa_factor);
                if options.colourspace == 2
                    % Convert to greyscale
                    A = rgb2grey(A);
                else
                    % Return only one channel for greyscale
                    A = check_greyscale(A);
                end
                % Outputs
                if options.im
                    imageData = A;
                end
                if options.alpha
                    imageData = A;
                    alpha = zeros(size(A, 1), size(A, 2), 'single');
                end
            end
            % Save the images
            if options.png
                res = options.magnify * get(0, 'ScreenPixelsPerInch') / 25.4e-3;
                imwrite(A, [options.name '.png'], 'ResolutionUnit', 'meter', 'XResolution', res, 'YResolution', res);
            end
            if options.bmp
                imwrite(A, [options.name '.bmp']);
            end
            % Save jpeg with given quality
            if options.jpg
                quality = options.quality;
                if isempty(quality)
                    quality = 95;
                end
                if quality > 100
                    imwrite(A, [options.name '.jpg'], 'Mode', 'lossless');
                else
                    imwrite(A, [options.name '.jpg'], 'Quality', quality);
                end
            end
            % Save tif images in cmyk if wanted (and possible)
            if options.tif
                if options.colourspace == 1 && size(A, 3) == 3
                    A = double(255 - A);
                    K = min(A, [], 3);
                    K_ = 255 ./ max(255 - K, 1);
                    C = (A(:,:,1) - K) .* K_;
                    M = (A(:,:,2) - K) .* K_;
                    Y = (A(:,:,3) - K) .* K_;
                    A = uint8(cat(3, C, M, Y, K));
                    clear C M Y K K_
                end
                append_mode = {'overwrite', 'append'};
                imwrite(A, [options.name '.tif'], 'Resolution', options.magnify*get(0, 'ScreenPixelsPerInch'), 'WriteMode', append_mode{options.append+1});
            end
        end

        % Now do the vector formats
        if isvector(options)
            % Set the default renderer to painters
            if ~options.renderer
                if isempty(findall(fig,'-property','FaceAlpha','-and','-not','FaceAlpha',1)) && ...
                        isempty(findall(fig,'type','patch'))
                    renderer = '-painters';
                else
                    % This is *MUCH* slower, but more accurate for patches and transparent annotations (issue #39)
                    renderer = '-opengl';
                end
            end
            % Generate some filenames
            tmp_nam = [tempname '.eps'];
            try
                % Ensure that the temp dir is writable (Javier Paredes 30/1/15)
                fid = fopen(tmp_nam,'w');
                fwrite(fid,1);
                fclose(fid);
                delete(tmp_nam);
                isTempDirOk = true;
            catch
                % Temp dir is not writable, so use the user-specified folder
                [dummy,fname,fext] = fileparts(tmp_nam); %#ok<ASGLU>
                fpath = fileparts(options.name);
                tmp_nam = fullfile(fpath,[fname fext]);
                isTempDirOk = false;
            end
            if isTempDirOk
                pdf_nam_tmp = [tempname '.pdf'];
            else
                pdf_nam_tmp = fullfile(fpath,[fname '.pdf']);
            end
            if options.pdf
                pdf_nam = [options.name '.pdf'];
                try copyfile(pdf_nam, pdf_nam_tmp, 'f'); catch, end  % fix for issue #65
            else
                pdf_nam = pdf_nam_tmp;
            end
            % Generate the options for print
            p2eArgs = {renderer, sprintf('-r%d', options.resolution)};
            if options.colourspace == 1  % CMYK
                % Issue #33: due to internal bugs in Matlab's print() function, we can't use its -cmyk option
                %p2eArgs{end+1} = '-cmyk';
            end
            if ~options.crop
                % Issue #56: due to internal bugs in Matlab's print() function, we can't use its internal cropping mechanism,
                % therefore we always use '-loose' (in print2eps.m) and do our own cropping (in crop_borders)
                %p2eArgs{end+1} = '-loose';
            end
            if any(strcmpi(varargin,'-depsc'))
                % Issue #45: lines in image subplots are exported in invalid color.
                % The workaround is to use the -depsc parameter instead of the default -depsc2
                p2eArgs{end+1} = '-depsc';
            end
            try
                % Generate an eps
                print2eps(tmp_nam, fig, [options.bb_padding, options.crop, options.fontswap], p2eArgs{:});
                % Remove the background, if desired
                if options.transparent && ~isequal(get(fig, 'Color'), 'none')
                    eps_remove_background(tmp_nam, 1 + using_hg2(fig));
                end
                % Fix colorspace to CMYK, if requested (workaround for issue #33)
                if options.colourspace == 1  % CMYK
                    % Issue #33: due to internal bugs in Matlab's print() function, we can't use its -cmyk option
                    change_rgb_to_cmyk(tmp_nam);
                end
                % Add a bookmark to the PDF if desired
                if options.bookmark
                    fig_nam = get(fig, 'Name');
                    if isempty(fig_nam)
                        warning('export_fig:EmptyBookmark', 'Bookmark requested for figure with no name. Bookmark will be empty.');
                    end
                    add_bookmark(tmp_nam, fig_nam);
                end
                % Generate a pdf
                eps2pdf(tmp_nam, pdf_nam_tmp, 1, options.append, options.colourspace==2, options.quality, options.gs_options);
                % Ghostscript croaks on % chars in the output PDF file, so use tempname and then rename the file
                try movefile(pdf_nam_tmp, pdf_nam, 'f'); catch, end
            catch ex
                % Delete the eps
                delete(tmp_nam);
                rethrow(ex);
            end
            % Delete the eps
            delete(tmp_nam);
            if options.eps
                try
                    % Generate an eps from the pdf
                    % since pdftops can't handle relative paths (e.g., '..\'), use a temp file
                    eps_nam_tmp = strrep(pdf_nam_tmp,'.pdf','.eps');
                    pdf2eps(pdf_nam, eps_nam_tmp);
                    movefile(eps_nam_tmp,  [options.name '.eps'], 'f');
                catch ex
                    if ~options.pdf
                        % Delete the pdf
                        delete(pdf_nam);
                    end
                    try delete(eps_nam_tmp); catch, end
                    rethrow(ex);
                end
                if ~options.pdf
                    % Delete the pdf
                    delete(pdf_nam);
                end
            end
        end

        % Revert the figure or close it (if requested)
        if cls || options.closeFig
            % Close the created figure
            close(fig);
        else
            % Reset the hardcopy mode
            set(fig, 'InvertHardcopy', old_mode);
            % Reset the axes limit and tick modes
            for a = 1:numel(Hlims)
                try
                    set(Hlims(a), 'XLimMode', Xlims{a}, 'YLimMode', Ylims{a}, 'ZLimMode', Zlims{a},...
                                  'XTickMode', Xtick{a}, 'YTickMode', Ytick{a}, 'ZTickMode', Ztick{a},...
                                  'XTickLabelMode', Xlabel{a}, 'YTickLabelMode', Ylabel{a}, 'ZTickLabelMode', Zlabel{a});
                catch
                    % ignore - fix issue #4 (using HG2 on R2014a and earlier)
                end
            end
            % Revert the tex-labels font weights
            try set(texLabels, 'FontWeight','bold'); catch, end
            % Revert annotation units
            for handleIdx = 1 : numel(annotationHandles)
                try
                    oldUnits = originalUnits{handleIdx};
                catch
                    oldUnits = originalUnits;
                end
                try set(annotationHandles(handleIdx),'Units',oldUnits); catch, end
            end
            % Revert figure units
            set(fig,'Units',oldFigUnits);
        end

        % Output to clipboard (if requested)
        if options.clipboard
            % Delete the output file if unchanged from the default name ('export_fig_out.png')
            if strcmpi(options.name,'export_fig_out')
                try
                    fileInfo = dir('export_fig_out.png');
                    if ~isempty(fileInfo)
                        timediff = now - fileInfo.datenum;
                        ONE_SEC = 1/24/60/60;
                        if timediff < ONE_SEC
                            delete('export_fig_out.png');
                        end
                    end
                catch
                    % never mind...
                end
            end

            % Save the image in the system clipboard
            % credit: Jiro Doke's IMCLIPBOARD: http://www.mathworks.com/matlabcentral/fileexchange/28708-imclipboard
            try
                error(javachk('awt', 'export_fig -clipboard output'));
            catch
                warning('export_fig -clipboard output failed: requires Java to work');
                return;
            end
            try
                % Import necessary Java classes
                import java.awt.Toolkit.*
                import java.awt.image.BufferedImage
                import java.awt.datatransfer.DataFlavor

                % Get System Clipboard object (java.awt.Toolkit)
                cb = getDefaultToolkit.getSystemClipboard();

                % Add java class (ImageSelection) to the path
                if ~exist('ImageSelection', 'class')
                    javaaddpath(fileparts(which(mfilename)), '-end');
                end

                % Get image size
                ht = size(imageData, 1);
                wd = size(imageData, 2);

                % Convert to Blue-Green-Red format
                try
                    imageData2 = imageData(:, :, [3 2 1]);
                catch
                    % Probably gray-scaled image (2D, without the 3rd [RGB] dimension)
                    imageData2 = imageData(:, :, [1 1 1]);
                end

                % Convert to 3xWxH format
                imageData2 = permute(imageData2, [3, 2, 1]);

                % Append Alpha data (unused - transparency is not supported in clipboard copy)
                alphaData2 = uint8(permute(255*alpha,[3,2,1])); %=255*ones(1,wd,ht,'uint8')
                imageData2 = cat(1, imageData2, alphaData2);

                % Create image buffer
                imBuffer = BufferedImage(wd, ht, BufferedImage.TYPE_INT_RGB);
                imBuffer.setRGB(0, 0, wd, ht, typecast(imageData2(:), 'int32'), 0, wd);

                % Create ImageSelection object from the image buffer
                imSelection = ImageSelection(imBuffer);

                % Set clipboard content to the image
                cb.setContents(imSelection, []);
            catch
                warning('export_fig -clipboard output failed: %s', lasterr); %#ok<LERR>
            end
        end

        % Don't output the data to console unless requested
        if ~nargout
            clear imageData alpha
        end
    catch err
        % Display possible workarounds before the error message
        if displaySuggestedWorkarounds
            if ~hadError,  fprintf(2, 'export_fig error. ');  end
            fprintf(2, 'Please ensure:\n');
            fprintf(2, '  that you are using the <a href="https://github.com/altmany/export_fig/archive/master.zip">latest version</a> of export_fig\n');
            fprintf(2, '  and that you have <a href="http://www.ghostscript.com">Ghostscript</a> installed\n');
            try
                if options.eps
                    fprintf(2, '  and that you have <a href="http://www.foolabs.com/xpdf">pdftops</a> installed\n');
                end
            catch
                % ignore - probably an error in parse_args
            end
            fprintf(2, '  and that you do not have <a href="matlab:which export_fig -all">multiple versions</a> of export_fig installed by mistake\n');
            fprintf(2, '  and that you did not made a mistake in the <a href="matlab:help export_fig">expected input arguments</a>\n');
            fprintf(2, '\nIf the problem persists, then please <a href="https://github.com/altmany/export_fig/issues">report a new issue</a>.\n\n');
        end
        rethrow(err)
    end
end

function [fig, options] = parse_args(nout, fig, varargin)
    % Parse the input arguments
    % Set the defaults
    options = struct(...
        'name', 'export_fig_out', ...
        'crop', true, ...
        'transparent', false, ...
        'renderer', 0, ... % 0: default, 1: OpenGL, 2: ZBuffer, 3: Painters
        'pdf', false, ...
        'eps', false, ...
        'png', false, ...
        'tif', false, ...
        'jpg', false, ...
        'bmp', false, ...
        'clipboard', false, ...
        'colourspace', 0, ... % 0: RGB/gray, 1: CMYK, 2: gray
        'append', false, ...
        'im',    nout == 1, ...
        'alpha', nout == 2, ...
        'aa_factor', 0, ...
        'bb_padding', 0, ...
        'magnify', [], ...
        'resolution', [], ...
        'bookmark', false, ...
        'closeFig', false, ...
        'quality', [], ...
        'update', false, ...
        'fontswap', true, ...
        'gs_options', {{}});
    native = false; % Set resolution to native of an image

    % Go through the other arguments
    skipNext = false;
    for a = 1:nargin-2
        if skipNext
            skipNext = false;
            continue;
        end
        if all(ishandle(varargin{a}))
            fig = varargin{a};
        elseif ischar(varargin{a}) && ~isempty(varargin{a})
            if varargin{a}(1) == '-'
                switch lower(varargin{a}(2:end))
                    case 'nocrop'
                        options.crop = false;
                    case {'trans', 'transparent'}
                        options.transparent = true;
                    case 'opengl'
                        options.renderer = 1;
                    case 'zbuffer'
                        options.renderer = 2;
                    case 'painters'
                        options.renderer = 3;
                    case 'pdf'
                        options.pdf = true;
                    case 'eps'
                        options.eps = true;
                    case 'png'
                        options.png = true;
                    case {'tif', 'tiff'}
                        options.tif = true;
                    case {'jpg', 'jpeg'}
                        options.jpg = true;
                    case 'bmp'
                        options.bmp = true;
                    case 'rgb'
                        options.colourspace = 0;
                    case 'cmyk'
                        options.colourspace = 1;
                    case {'gray', 'grey'}
                        options.colourspace = 2;
                    case {'a1', 'a2', 'a3', 'a4'}
                        options.aa_factor = str2double(varargin{a}(3));
                    case 'append'
                        options.append = true;
                    case 'bookmark'
                        options.bookmark = true;
                    case 'native'
                        native = true;
                    case 'clipboard'
                        options.clipboard = true;
                        options.im = true;
                        options.alpha = true;
                    case 'svg'
                        msg = ['SVG output is not supported by export_fig. Use one of the following alternatives:\n' ...
                               '  1. saveas(gcf,''filename.svg'')\n' ...
                               '  2. plot2svg utility: http://github.com/jschwizer99/plot2svg\n' ...
                               '  3. export_fig to EPS/PDF, then convert to SVG using generic (non-Matlab) tools\n'];
                        error(sprintf(msg)); %#ok<SPERR>
                    case 'update'
                        % Download the latest version of export_fig into the export_fig folder
                        try
                            zipFileName = 'https://github.com/altmany/export_fig/archive/master.zip';
                            folderName = fileparts(which(mfilename('fullpath')));
                            targetFileName = fullfile(folderName, datestr(now,'yyyy-mm-dd.zip'));
                            urlwrite(zipFileName,targetFileName);
                        catch
                            error('Could not download %s into %s\n',zipFileName,targetFileName);
                        end

                        % Unzip the downloaded zip file in the export_fig folder
                        try
                            unzip(targetFileName,folderName);
                        catch
                            error('Could not unzip %s\n',targetFileName);
                        end
                    case 'nofontswap'
                        options.fontswap = false;
                    otherwise
                        try
                            if strcmpi(varargin{a}(1:2),'-d')
                                varargin{a}(2) = 'd';  % ensure lowercase 'd'
                                options.gs_options{end+1} = varargin{a};
                            else
                                val = str2double(regexp(varargin{a}, '(?<=-(m|M|r|R|q|Q|p|P))-?\d*.?\d+', 'match'));
                                if isempty(val) || isnan(val)
                                    % Issue #51: improved processing of input args (accept space between param name & value)
                                    val = str2double(varargin{a+1});
                                    if isscalar(val) && ~isnan(val)
                                        skipNext = true;
                                    end
                                end
                                if ~isscalar(val) || isnan(val)
                                    error('option %s is not recognised or cannot be parsed', varargin{a});
                                end
                                switch lower(varargin{a}(2))
                                    case 'm'
                                        options.magnify = val;
                                    case 'r'
                                        options.resolution = val;
                                    case 'q'
                                        options.quality = max(val, 0);
                                    case 'p'
                                        options.bb_padding = val;
                                end
                            end
                        catch
                            error(['Unrecognized export_fig input option: ''' varargin{a} '''']);
                        end
                end
            else
                [p, options.name, ext] = fileparts(varargin{a});
                if ~isempty(p)
                    options.name = [p filesep options.name];
                end
                switch lower(ext)
                    case {'.tif', '.tiff'}
                        options.tif = true;
                    case {'.jpg', '.jpeg'}
                        options.jpg = true;
                    case '.png'
                        options.png = true;
                    case '.bmp'
                        options.bmp = true;
                    case '.eps'
                        options.eps = true;
                    case '.pdf'
                        options.pdf = true;
                    case '.fig'
                        % If no open figure, then load the specified .fig file and continue
                        if isempty(fig)
                            fig = openfig(varargin{a},'invisible');
                            varargin{a} = fig;
                            options.closeFig = true;
                        else
                            % save the current figure as the specified .fig file and exit
                            saveas(fig(1),varargin{a});
                            fig = -1;
                            return
                        end
                    case '.svg'
                        msg = ['SVG output is not supported by export_fig. Use one of the following alternatives:\n' ...
                               '  1. saveas(gcf,''filename.svg'')\n' ...
                               '  2. plot2svg utility: http://github.com/jschwizer99/plot2svg\n' ...
                               '  3. export_fig to EPS/PDF, then convert to SVG using generic (non-Matlab) tools\n'];
                        error(sprintf(msg)); %#ok<SPERR>
                    otherwise
                        options.name = varargin{a};
                end
            end
        end
    end

    % Quick bail-out if no figure found
    if isempty(fig),  return;  end

    % Do border padding with repsect to a cropped image
    if options.bb_padding
        options.crop = true;
    end

    % Set default anti-aliasing now we know the renderer
    if options.aa_factor == 0
        try isAA = strcmp(get(ancestor(fig, 'figure'), 'GraphicsSmoothing'), 'on'); catch, isAA = false; end
        options.aa_factor = 1 + 2 * (~(using_hg2(fig) && isAA) | (options.renderer == 3));
    end

    % Convert user dir '~' to full path
    if numel(options.name) > 2 && options.name(1) == '~' && (options.name(2) == '/' || options.name(2) == '\')
        options.name = fullfile(char(java.lang.System.getProperty('user.home')), options.name(2:end));
    end

    % Compute the magnification and resolution
    if isempty(options.magnify)
        if isempty(options.resolution)
            options.magnify = 1;
            options.resolution = 864;
        else
            options.magnify = options.resolution ./ get(0, 'ScreenPixelsPerInch');
        end
    elseif isempty(options.resolution)
        options.resolution = 864;
    end

    % Set the default format
    if ~isvector(options) && ~isbitmap(options)
        options.png = true;
    end

    % Check whether transparent background is wanted (old way)
    if isequal(get(ancestor(fig(1), 'figure'), 'Color'), 'none')
        options.transparent = true;
    end

    % If requested, set the resolution to the native vertical resolution of the
    % first suitable image found
    if native && isbitmap(options)
        % Find a suitable image
        list = findall(fig, 'Type','image', 'Tag','export_fig_native');
        if isempty(list)
            list = findall(fig, 'Type','image', 'Visible','on');
        end
        for hIm = list(:)'
            % Check height is >= 2
            height = size(get(hIm, 'CData'), 1);
            if height < 2
                continue
            end
            % Account for the image filling only part of the axes, or vice
            % versa
            yl = get(hIm, 'YData');
            if isscalar(yl)
                yl = [yl(1)-0.5 yl(1)+height+0.5];
            else
                if ~diff(yl)
                    continue
                end
                yl = yl + [-0.5 0.5] * (diff(yl) / (height - 1));
            end
            hAx = get(hIm, 'Parent');
            yl2 = get(hAx, 'YLim');
            % Find the pixel height of the axes
            oldUnits = get(hAx, 'Units');
            set(hAx, 'Units', 'pixels');
            pos = get(hAx, 'Position');
            set(hAx, 'Units', oldUnits);
            if ~pos(4)
                continue
            end
            % Found a suitable image
            % Account for stretch-to-fill being disabled
            pbar = get(hAx, 'PlotBoxAspectRatio');
            pos = min(pos(4), pbar(2)*pos(3)/pbar(1));
            % Set the magnification to give native resolution
            options.magnify = (height * diff(yl2)) / (pos * diff(yl));
            break
        end
    end
end

function A = downsize(A, factor)
    % Downsample an image
    if factor == 1
        % Nothing to do
        return
    end
    try
        % Faster, but requires image processing toolbox
        A = imresize(A, 1/factor, 'bilinear');
    catch
        % No image processing toolbox - resize manually
        % Lowpass filter - use Gaussian as is separable, so faster
        % Compute the 1d Gaussian filter
        filt = (-factor-1:factor+1) / (factor * 0.6);
        filt = exp(-filt .* filt);
        % Normalize the filter
        filt = single(filt / sum(filt));
        % Filter the image
        padding = floor(numel(filt) / 2);
        for a = 1:size(A, 3)
            A(:,:,a) = conv2(filt, filt', single(A([ones(1, padding) 1:end repmat(end, 1, padding)],[ones(1, padding) 1:end repmat(end, 1, padding)],a)), 'valid');
        end
        % Subsample
        A = A(1+floor(mod(end-1, factor)/2):factor:end,1+floor(mod(end-1, factor)/2):factor:end,:);
    end
end

function A = rgb2grey(A)
    A = cast(reshape(reshape(single(A), [], 3) * single([0.299; 0.587; 0.114]), size(A, 1), size(A, 2)), class(A)); %#ok<ZEROLIKE>
end

function A = check_greyscale(A)
    % Check if the image is greyscale
    if size(A, 3) == 3 && ...
            all(reshape(A(:,:,1) == A(:,:,2), [], 1)) && ...
            all(reshape(A(:,:,2) == A(:,:,3), [], 1))
        A = A(:,:,1); % Save only one channel for 8-bit output
    end
end

function eps_remove_background(fname, count)
    % Remove the background of an eps file
    % Open the file
    fh = fopen(fname, 'r+');
    if fh == -1
        error('Not able to open file %s.', fname);
    end
    % Read the file line by line
    while count
        % Get the next line
        l = fgets(fh);
        if isequal(l, -1)
            break; % Quit, no rectangle found
        end
        % Check if the line contains the background rectangle
        if isequal(regexp(l, ' *0 +0 +\d+ +\d+ +r[fe] *[\n\r]+', 'start'), 1)
            % Set the line to whitespace and quit
            l(1:regexp(l, '[\n\r]', 'start', 'once')-1) = ' ';
            fseek(fh, -numel(l), 0);
            fprintf(fh, l);
            % Reduce the count
            count = count - 1;
        end
    end
    % Close the file
    fclose(fh);
end

function b = isvector(options)
    b = options.pdf || options.eps;
end

function b = isbitmap(options)
    b = options.png || options.tif || options.jpg || options.bmp || options.im || options.alpha;
end

% Helper function
function A = make_cell(A)
    if ~iscell(A)
        A = {A};
    end
end

function add_bookmark(fname, bookmark_text)
    % Adds a bookmark to the temporary EPS file after %%EndPageSetup
    % Read in the file
    fh = fopen(fname, 'r');
    if fh == -1
        error('File %s not found.', fname);
    end
    try
        fstrm = fread(fh, '*char')';
    catch ex
        fclose(fh);
        rethrow(ex);
    end
    fclose(fh);

    % Include standard pdfmark prolog to maximize compatibility
    fstrm = strrep(fstrm, '%%BeginProlog', sprintf('%%%%BeginProlog\n/pdfmark where {pop} {userdict /pdfmark /cleartomark load put} ifelse'));
    % Add page bookmark
    fstrm = strrep(fstrm, '%%EndPageSetup', sprintf('%%%%EndPageSetup\n[ /Title (%s) /OUT pdfmark',bookmark_text));

    % Write out the updated file
    fh = fopen(fname, 'w');
    if fh == -1
        error('Unable to open %s for writing.', fname);
    end
    try
        fwrite(fh, fstrm, 'char*1');
    catch ex
        fclose(fh);
        rethrow(ex);
    end
    fclose(fh);
end

function set_tick_mode(Hlims, ax)
    % Set the tick mode of linear axes to manual
    % Leave log axes alone as these are tricky
    M = get(Hlims, [ax 'Scale']);
    if ~iscell(M)
        M = {M};
    end
    M = cellfun(@(c) strcmp(c, 'linear'), M);
    set(Hlims(M), [ax 'TickMode'], 'manual');
    %set(Hlims(M), [ax 'TickLabelMode'], 'manual');  % this hides exponent label in HG2!
end

function change_rgb_to_cmyk(fname)  % convert RGB => CMYK within an EPS file
    % Do post-processing on the eps file
    try
        % Read the EPS file into memory
        fstrm = read_write_entire_textfile(fname);

        % Replace all gray-scale colors
        fstrm = regexprep(fstrm, '\n([\d.]+) +GC\n', '\n0 0 0 ${num2str(1-str2num($1))} CC\n');

        % Replace all RGB colors
        fstrm = regexprep(fstrm, '\n[0.]+ +[0.]+ +[0.]+ +RC\n', '\n0 0 0 1 CC\n');  % pure black
        fstrm = regexprep(fstrm, '\n([\d.]+) +([\d.]+) +([\d.]+) +RC\n', '\n${sprintf(''%.4g '',[1-[str2num($1),str2num($2),str2num($3)]/max([str2num($1),str2num($2),str2num($3)]),1-max([str2num($1),str2num($2),str2num($3)])])} CC\n');

        % Overwrite the file with the modified contents
        read_write_entire_textfile(fname, fstrm);
    catch
        % never mind - leave as is...
    end
end

author: mtex-toolbox

repo: mtex

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


text:
function [shifts,Mf] = register_frame(Yt,fftTempMat,fftTempPatches,patches,options)

minY = min(Yt(:));
maxY = max(Yt(:));
dim = size(Yt);
Yc = split_frame(Yt,patches);
fftY = fft(fft(fft(Yc,[],1),[],2),[],3);
nd = ndims(Yt);
n_patches = [length(unique(patches(:,1))),length(unique(patches(:,3))),length(unique(patches(:,5)))];
shifts = zeros([n_patches,nd]);
N_patches = size(patches,1);

if N_patches > 1
    if nd == 2; out_rig = dftregistration_min_max(fftTempMat,fftn(Yt),options.us_fac,-options.max_shift,options.max_shift,options.phase_flag); lb = out_rig(3:4); ub = out_rig(3:4); end
    if nd == 3; out_rig = dftregistration_min_max_3d(fftTempMat,fftn(Yt),options.us_fac,-options.max_shift,options.max_shift,options.phase_flag); lb = out_rig(3:5); ub = out_rig(3:5); end
    max_dev = options.max_dev;
else
    lb = -options.max_shift(1,nd);
    ub = options.max_shift(1,nd);
    max_dev = 0*options.max_dev;
end

for i = 1:N_patches
    [ix,iy,iz] = ind2sub(n_patches,i);
    if nd == 2
        output = dftregistration_min_max(fftTempPatches(:,:,:,i),fftY(:,:,:,i),options.us_fac,lb-max_dev(1:2),ub+max_dev(1:2),options.phase_flag);
    elseif nd == 3
        output = dftregistration_min_max_3d(fftTempPatches(:,:,:,i),fftY(:,:,:,i),options.us_fac,lb-max_dev,ub+max_dev,options.phase_flag);
        shifts(ix,iy,iz,3) = output(5);
    end

    shifts(ix,iy,iz,1) = output(3);
    shifts(ix,iy,iz,2) = output(4);
end

if nd == 3
    shifts_up = zeros([dim,3]);
    do = size(shifts_up)./size(shifts);
    ds = size(shifts);
    if (0)
        [Xq,Yq,Zq] = meshgrid(linspace((2+0*1/do(2))/2,ds(2)+0*(0-1/do(2))/2,dim(2)),linspace((2+0*1/do(1))/2,ds(1)+0*(0-1/do(1))/2,dim(1)),linspace((2+0*1/do(3))/2,ds(3)+0*(0-1/do(3))/2,dim(3)));
    else
        [Xq,Yq,Zq] = meshgrid(linspace((1+1/do(2))/2,ds(2)+(1-1/do(2))/2,dim(2)),linspace((1+1/do(1))/2,ds(1)+(1-1/do(1))/2,dim(1)),linspace((1+1/do(3))/2,ds(3)+(1-1/do(3))/2,dim(3)));
        Xq(Xq<1)=1; Xq(Xq>dim(2))=dim(2);
        Yq(Yq<1)=1; Yq(Yq>dim(1))=dim(1);
        Zq(Zq<1)=1; Zq(Zq>dim(3))=dim(3);
    end
    if numel(shifts) > 3
%       tform = affine3d(diag([options.mot_uf(:);1]));
%       for dm = 1:3; shifts_up(:,:,:,dm) = imwarp(shifts(:,:,:,dm),tform,'OutputView',imref3d([options.d1,options.d2,options.d3])); end
        for dm = 1:3; shifts_up(:,:,:,dm) = interp3(shifts(:,:,:,dm),Xq,Yq,Zq,'makima'); end
    else
        for dm = 1:3; shifts_up(:,:,:,dm) = shifts(dm); end
    end
    shifts_up(2:2:end,:,:,2) = shifts_up(2:2:end,:,:,2) + options.col_shift;
    Mf = imwarp(Yt,-cat(4,shifts_up(:,:,:,2),shifts_up(:,:,:,1),shifts_up(:,:,:,3)),'bicubic','FillValues',options.add_value);
else
    shifts_up = imresize(shifts,[options.d1,options.d2]);
    shifts_up(2:2:end,:,2) = shifts_up(2:2:end,:,2) + options.col_shift;
    Mf = imwarp(Yt,-cat(3,shifts_up(:,:,2),shifts_up(:,:,1)),'bicubic','FillValues',options.add_value);
end
Mf(Mf<minY)=minY;
Mf(Mf>maxY)=maxY;

end

author: flatironinstitute

repo: NoRMCorre


text:
function out=prox_max(x,alpha)
%PROX_MAX computes the proximal operator of the function alpha*max(x(:))
%
%  Usage:
%  out = PROX_MAX(x,alpha)
%  ===========================================
%  INPUT:
%  x - point to be projected (vector/matrix)
%  alpha - positive scalar
%  ===========================================
%  Output:
%  out - proximal operator at x

% This file is part of the FOM package - a collection of first order methods for solving convex optimization problems
% Copyright (C) 2017 Amir and Nili Beck
%
% This program 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 program 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 program.  If not, see <http://www.gnu.org/licenses/>.

if (nargin < 2)
    error ('usage: prox_max(x,alpha)') ;
end

if (alpha < 0)
    error('usage: prox_max(x,alpha) - alpha should be positive')
end

out = x - alpha * proj_simplex (x/alpha,1,'eq') ;

end

author: hiroyuki-kasai

repo: SGDLibrary


text:
function varargout = ctrInit(varargin)
% CTRINIT M-file for ctrInit.fig
%      CTRINIT, by itself, creates a new CTRINIT or raises the existing
%      singleton*.
%
%      H = CTRINIT returns the handle to a new CTRINIT or the handle to
%      the existing singleton*.
%
%      CTRINIT('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in CTRINIT.M with the given input arguments.
%
%      CTRINIT('Property','Value',...) creates a new CTRINIT or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before ctrInit_OpeningFunction gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to ctrInit_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 ctrInit

% Last Modified by GUIDE v2.5 17-Sep-2008 16:05:29

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @ctrInit_OpeningFcn, ...
                   'gui_OutputFcn',  @ctrInit_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

return;

% --- Executes just before ctrInit is made visible.
function ctrInit_OpeningFcn(hObject, eventdata, handles, varargin)
% This function has no output args, see OutputFcn.
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
% varargin   command line arguments to ctrInit (see VARARGIN)

% Choose default command line output for ctrInit
handles.output = hObject;

% Update handles structure
guidata(hObject, handles);

% UIWAIT makes ctrInit wait for user response (see UIRESUME)
% uiwait(handles.figure1);
btnTimeStamp_Callback(hObject, eventdata, handles);

menuEditRefresh_Callback(hObject, eventdata, handles)

return;

% --- Outputs from this function are returned to the command line.
function varargout = ctrInit_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;

return;

% --- Executes on button press in btnCreateContrackFiles.
function btnCreateContrackFiles_Callback(hObject, eventdata, handles)
% Button:  Generate Fibers

params = getParams(hObject);

% This does ALMOST EVERYTHING
% 1. Creates wmprob.nii.gz
% 2. Creates pdf.nii.gz
% 3. Create ROI mask.nii.gz
% 4. Create the ctrSampler_timestamp.txt file
params = ctrInitParamsFile(params);

% One last thing: make the ctrScript_timestamp.sh file
ctrScript(params);

return;

% --- Executes on button press in btnDT6file.
function btnDT6File_Callback(hObject, eventdata, handles)
% Button:  Browse for DT6 file
hndl = guihandles(hObject);
fullName = mrvSelectFile('r','*','DT6');
if isempty(fullName), return; end
set(hndl.editFilenameDti6,'String',fullName);
menuEditRefresh_Callback(gcbf);
return;

% --- Executes on button press in btnROI1file.
function btnROI1File_Callback(hObject, eventdata, handles)
% Button:  Browse for ROI1 file
hndl = guihandles(hObject);
startDir = fileparts(get(hndl.editFilenameDti6,'String'));
fullName = mrvSelectFile('r','*','ROI 1',startDir);
if isempty(fullName), return; end
set(hndl.editFilenameROI1,'String',fullName);
menuEditRefresh_Callback(gcbf);
return;

% --- Executes on button press in btnROI2file.
function btnROI2File_Callback(hObject, eventdata, handles)
% Button:  Browse for ROI1 file
hndl = guihandles(hObject);
startDir = fileparts(get(hndl.editFilenameROI1,'String'));
if(isempty(startDir))
    startDir = fileparts(get(hndl.editFilenameDti6,'String'));
end
fullName = mrvSelectFile('r','*','ROI 2',startDir);
if isempty(fullName), return; end
set(hndl.editFilenameROI2,'String',fullName);
menuEditRefresh_Callback(gcbf);
return;

function editFilenameDti6_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return;

% --- Executes during object creation, after setting all properties.
function editFilenameDti6_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end
return;

function editFilenameROI1_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return;

% --- Executes during object creation, after setting all properties.
function editFilenameROI1_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end
return;

% --- Executes during object creation, after setting all properties.
function editFilenameROI2_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end
return;

function editFilenameROI2_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return;

% --------------------------------------------------------------------
function menuEdit_Callback(hObject, eventdata, handles)
return;

% --------------------------------------------------------------------
function menuFile_Callback(hObject, eventdata, handles)
return;

% --------------------------------------------------------------------
function menuFileLoad_Callback(hObject, eventdata, handles)
% File | Load params
ctr = ctrLoad;
setParams(hObject,ctr);
return;
% --------------------------------------------------------------------
function menuFileSave_Callback(hObject, eventdata, handles)
% File | Save params
params = getParams(hObject);
ctrInitParamsFile(params);
return;

% --------------------------------------------------------------------
function menuFileClose_Callback(hObject, eventdata, handles)
% File | Close
closereq;
return;

% --------------------------------------------------------------------
function menuFileSaveClose_Callback(hObject, eventdata, handles)
% File | Save and close
params = getParams(hObject);
ctrInitParamsFile(params);
closereq;
return;

% --------------------------------------------------------------------
function menuHelpVistaWiki_Callback(hObject, eventdata, handles)
% Help | Vista Wiki ConTrack
web('http://white.stanford.edu/newlm/index.php/DTI#ConTrack','-browser');
return;

function editDesiredSamples_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return

% --- Executes during object creation, after setting all properties.
function editDesiredSamples_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end
return;

function editMaxNodes_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return;
% --- Executes during object creation, after setting all properties.
function editMaxNodes_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end
return;

function editMinNodes_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return;

% --- Executes during object creation, after setting all properties.
function editMinNodes_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end
return;

function editStepSize_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return;

% --- Executes during object creation, after setting all properties.
function editStepSize_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end
return;

% --------------------------------------------------------------------
function menuHelp_Callback(hObject, eventdata, handles)
return;

% --------------------------------------------------------------------
function menuEditClearFilenames_Callback(hObject, eventdata, handles)
% Edit | Clear file names

gd = guidata(hObject);

set(gd.editFilenameDti6,'String','');
set(gd.editFilenameROI1,'String','');
set(gd.editFilenameROI2,'String','');
set(gd.editWaypointFilename,'String','');

return;

% --------------------------------------------------------------------
function menuEditResetDefaultParameters_Callback(hObject, eventdata, handles)
disp('Todo: Reset default parameters');
return

% --- Executes on button press in checkboxROI1.
function checkboxROI1_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return

% --- Executes on button press in checkboxROI2.
function checkboxROI2_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return

% --- Executes on button press in checkboxPDDPDF.
function checkboxPDDPDF_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return

% --- Executes on button press in checkboxWMFile.
function checkboxWMFile_Callback(hObject, eventdata, handles)
menuEditRefresh_Callback(gcbf);
return

% --------------------------------------------------------------------
function menuEditFullParameters_Callback(hObject, eventdata, handles)
disp('ToDO:  Bring up window with full parameter list and store the list');

params = getParams(hObject);
[params]=ctrSetParamsDlg(hObject,params);
menuEditRefresh_Callback(gcbf);
%guidata(hObject);
return;

% --------------------------------------------------------------------
function menuEditRefresh_Callback(hObject, eventdata, handles)
% Edit | Refresh
% Updates changes to the information panel

gd = guidata(hObject);

% File name section - Show user computed file names
str = get(gd.editFilenameDti6,'String');
if isempty(str), return; end

% Otherwise, we fill in some of the fields: Subject directory
subjDir = fileparts(str);
if length(subjDir) > 30, str = subjDir((end-30):end); else str = subjDir; end
txt = sprintf('Subject''s DT6 dir: ... %s\n\n',str);

% Processed ROI file, derived from names and a time stamp. We should
% probably store the time stamp in the guidata of the window. Perhaps we
% update the time stamp whenever we change a file name?
[p,str1] = fileparts(get(gd.editFilenameROI1,'String'));
[p,str2] = fileparts(get(gd.editFilenameROI2,'String'));
tString = get(gd.txtTimeStamp,'String');
str = ctrFilename(str1,str2,tString,'roi');
txt = sprintf('%sROI mask: %s\n\n',txt,str);

% Attach the text to the window
set(gd.txtInfo,'String',txt);

return;

% --- Executes on button press in btnTimeStamp.
function btnTimeStamp_Callback(hObject, eventdata, handles)
% Time Stamp button

gd  = guidata(hObject);
str = datestr(now,30);
set(gd.txtTimeStamp,'String',str);
menuEditRefresh_Callback(hObject);

return;

%---------------
function params = getParams(hObject)
% Create the parameters structure that we send to ctInitParams

hndl = guidata(hObject);

% Figure out what the user selected and send it to the script that
% generates the fibers.
params.dt6File      = get(hndl.editFilenameDti6,'String');
params.roi1File     = get(hndl.editFilenameROI1,'String');
params.roi2File     = get(hndl.editFilenameROI2,'String');
params.dSamples     = str2double(get(hndl.editDesiredSamples','string'));
params.maxNodes     = str2double(get(hndl.editMaxNodes','string'));
params.minNodes     = str2double(get(hndl.editMinNodes','string'));
params.stepSize     = str2double(get(hndl.editStepSize','string'));
params.pddpdf       = get(hndl.checkboxPDDPDF,'Value');
params.wm           = get(hndl.checkboxWMFile,'Value');
params.roi1Seed     = get(hndl.checkboxROI1,'Value');
params.roi2Seed     = get(hndl.checkboxROI2,'Value');
params.timeStamp    = get(hndl.txtTimeStamp,'String');

return;

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

function setParams(hObject,ctr)
% Set the interface based on the contrack parameters structure

if notDefined('ctr'), error('contrack structure required.'); end

hndl = guidata(hObject);

%Copy the contrack structure to the interface
%The ctr structure could have more stuff, like this file information

% set(hndl.editFilenameDti6,'String',params.dt6File);
% set(hndl.editFilenameROI1,'String',params.roi1File);
% set(hndl.editFilenameROI2,'String',params.roi2File);
% set(hndl.editWaypointFilename,'String',params.wayPointFile);

set(hndl.editDesiredSamples,'string',num2str(ctrGet(ctr,'desired_samples')));
set(hndl.editMaxNodes,'string',num2str(ctrGet(ctr,'max_nodes')));
set(hndl.editMinNodes,'string',num2str(ctrGet(ctr,'min_nodes')));
set(hndl.editStepSize,'string',num2str(ctrGet(ctr,'step_size')));
set(hndl.checkboxWMFile,'Value',1);
set(hndl.checkboxROI1,'Value',1);
set(hndl.checkboxROI2,'Value',1);
set(hndl.txtTimeStamp,'String',datestr(now,30));

menuEditRefresh_Callback(hObject);

return;

% --- Executes on button press in btnFiberGeneration.
function btnFiberGeneration_Callback(hObject, eventdata, handles)
% hObject    handle to btnFiberGeneration (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

author: vistalab

repo: vistasoft

save_path: github-repos/MATLAB/vistalab-vistasoft


text:
function maximize( varargin )

%MAXIMIZE Specifiies a concave (or affine) objective to be maximized.

if nargin < 1,
    cvx_throw( 'Objective expression missing.' );
elseif iscellstr( varargin ),
    x = evalin( 'caller', sprintf( '%s ', varargin{:} ) );
elseif nargin > 1,
    cvx_throw( 'Too many input arguments.' );
else
    x = varargin{1};
end
evalin( 'caller', 'cvx_verify' );
cvx_pushobj( 'maximize', x );

% Copyright 2005-2014 CVX Research, Inc.
% See the file LICENSE.txt for full copyright information.
% The command 'cvx_where' will show where this file is located.

author: yu-jiang

repo: radpbook

save_path: github-repos/MATLAB/yu-jiang-radpbook


%==============================================================================
% This code is part of the Matlab-based toolbox
% FAIR - Flexible Algorithms for Image Registration.
% For details see
% - https://github.com/C4IR and
% - http://www.siam.org/books/fa06/
%==============================================================================
% function [ML,minLevel,maxLevel,fig] = getMultilevel(IS,omega,m,varargin)
%
% compute a multi-level representation of an image IS or a list of images
% IS = {T,R,MT,MR,...}
%
%  Input:
%      IS        an image or a list of images
%      omega     domain specification
%      m         number of discretization points
%      varargin  optional parameter, see below
%
% Output:
%   ML       struct containing multi-level representation
%                   ML{j} =  {T_j,R_j,omega,m_j}
%   minLevel     coarsest level
%   maxLevel     finest   level
%   fig          handle to the graphical output
%
% A continuous representation of an image serves as a starting point. Using a cell
% centered grid, the data is replace by interpolated values. The representation on
% a coarser level is obtained by averaging over adjacent cells (see the filters for
% options).
%
% ML{level} is a structure containing the image(s), omega, and m,
% where level runs 0:maxLevel, ML{level} = [] for level < minLevel.
% Note that we assume m=2^p, p integer; the original data size is arbitrary.
%==============================================================================

function [ML,minLevel,maxLevel,fig] = getMultilevel(IS,omega,m,varargin)

if nargin == 0,    % show help and provide minimal example
  help(mfilename);
  runMinimalExample;
  ML = 'endOfMinimalExample';
  return;
end;

if nargin == 1,
  % simply return ML, minLevel, maxLevel
  ML = IS; maxLevel = length(ML); fig = 0;
  for minLevel=maxLevel:-1:1,
    if minLevel == 1 || isempty(ML{minLevel-1}), break; end;
  end;
  return;
end;

% start the work
fig      = 2;                  % figure number for output
dopause  = 0;                  % make pause for demonstrations
filter   = 'block';            % discrete smoothing kernel
minLevel = 3;                  % minimal level size
names    = {'T','R','Q','Q1','Q2','Q3','Q4','Q5'};
restrictdim = ones(size(m));   % by default: restrict all dimensions
imgModel = @linearInter;       % use linear interpolation as a default
                               % image model

for k=1:2:length(varargin),    % overwrite defaults
  eval([varargin{k},'=varargin{',int2str(k+1),'};']);
end;

if ~iscell(IS), IS = {IS}; end;% in case of a single image input
lenIS = length(IS);                 % number of images to be handled
dim   = length(omega)/2;            % spacial dimension
maxLevel = ceil(log2(min(m(restrictdim==1))));      % finest level

% messaging
msg = sprintf('%s(%d image(s), filter=%s, %dD data,level=%d:%d)',...
  mfilename,lenIS,filter,dim,minLevel,maxLevel);

fprintf('%s, figure=%d\n',msg,fig);
fprintf('note: model from data of any dimension is sampled on a grid of size m\n');
fprintf('      at present it is assume that m = 2.^p, p an interger\n');
fprintf('      ML is an array of cells, where the i-th entry contains the data\n');
fprintf('      of size 2^p, with some modifications is not all dimensions are equal\n');

% some output to the console
fprintf('%s: [',mfilename);

omegak = @(k)omega(min(1+(k>1),size(omega,1)),:);

% start the loop from maxlevel to minlevel
for level = maxLevel:-1:minLevel
  if level == maxLevel,
    % replace data by the interpolant on x
    for k=1:lenIS,
      sizeData = size(IS{k});
      % note: data is either m1-by-....-by-md
      %                   or m1-by-....-by-md-by-nrVolume
      % but the 1d case is tricky

      if (dim>1) && (length(sizeData) > dim),
        nrVolume = sizeData(end),
        sizeData(end) = [];
        keyboard
      else
        nrVolume = 1;
      end;

      data  = reshape(IS{k},prod(sizeData),nrVolume);
      block = zeros(prod(m),nrVolume);
      xc    = getCellCenteredGrid(omegak(k),m);
      % create dsample of size m from image model (interpolation) obtained from data
      for vol=1:size(data,2),
        block(:,vol) = imgModel(reshape(data(:,vol),sizeData),omegak(k),xc);
      end

      IS{k} = reshape(block,[m,nrVolume]);
    end;
    fprintf('%d',level);
  else
    % restrict the image to the coarser grid
    L = ML{level+1};
    for k=1:lenIS, % run over all images
      for j=1:dim, % run over all dimensions
        if restrictdim(j),
            IS{k} = restrict(IS{k},dim,j,filter);
        end
      end;
    end;
    fprintf(',%d',level);
  end;

  % store the current level data in a struct
  L.m = size(IS{1}); L.m = L.m(1:dim); L.omega = omega;
  for k= 1:lenIS,
    L = setfield(L,names{k},IS{k});
  end;
  ML{level} = L;

  if fig, % do some plots
    if level == maxLevel,
       FAIRfigure(fig,'figname',msg,'position','default');
    end;

    str = @(k) sprintf('%s(level=%d), %s',...
      names{k},level,sprintf('m=[%s]',sprintf(' %d',L.m)));
    p0 = level-minLevel+1;
    dp = (maxLevel-minLevel+1);
    figureh(fig);
    for k=1:lenIS,
      xc = getCellCenteredGrid(omegak(k),L.m);
      subplot(lenIS,dp,p0+(k-1)*dp);
      viewImage(imgModel(IS{k},omegak(k),xc),omegak(k),L.m);
      title(str(k));
    end;
    if dopause, pause; else drawnow; end;
  end;
end;
fprintf('] done\n');

%------------------------------------------------------------------------------
function T = restrict(T,dim,j,filter)
J = [j,setdiff(1:dim,j)]; J = [J dim+1]; % bring dimension j to front
T = permute(T,J);
if rem(size(T,1),2), T(end+1,:,:,:) = T(end,:,:,:); end;
switch filter,
  case 'gaussian',
    Ta = zeros(size(T)+[2,zeros(1,ndims(T)-1)]);
    Ta(2:end-1,:,:,:) = T; Ta([1,end],:,:,:) = Ta([2,end-1],:,:,:);
    T =   0.125 * Ta(1:end-3,:,:,:) + 0.375 * Ta(2:end-2,:,:,:) ...
        + 0.375 * Ta(3:end-1,:,:,:) + 0.125 * Ta(4:end,:,:,:);
    T = T(1:2:end,:,:,:);
  case 'block'
    T = (T(1:2:end,:,:,:)+T(2:2:end,:,:,:))/2;
  case 'harmonic'
    T = 1./(T+eps);
    T = (T(1:2:end,:,:,:)+T(2:2:end,:,:,:))/2;
    T = 1./T;
end;
T = ipermute(T,J); % bring j back home

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

function runMinimalExample
Tdata = double(flipud(imread('hands-T.jpg'))');
Rdata = double(flipud(imread('hands-R.jpg'))');
omega = [0,20,0,25]; % specify physical domain
m     = size(Tdata);

% set view options and interpolation options and initialize viewer and interpolator
viewPara = {'viewImage','viewImage2D','colormap','gray(256)'};
viewImage('reset',viewPara{:});
% create multilevel representation of the data
ML = getMultilevel({Tdata,Rdata},omega,m,'fig',2);

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


text:
function [ p, rank ] = perm_tj_successor ( n, p, rank )

%*****************************************************************************80
%
%% PERM_TJ_SUCCESSOR computes the Trotter-Johnson permutation successor.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    12 January 2011
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    Donald Kreher, Douglas Simpson,
%    Combinatorial Algorithms,
%    CRC Press, 1998,
%    ISBN: 0-8493-3988-X,
%    LC: QA164.K73.
%
%  Parameters:
%
%    Input, integer N, the number of values being permuted.
%    N must be positive.
%
%    Input/output, integer P(N), describes the permutation.
%    P(I) is the item which is permuted into the I-th place
%    by the permutation.
%
%    Input/output, integer RANK, the rank.
%    If RANK = -1 on input, then the routine understands that this is
%    the first call, and that the user wishes the routine to supply
%    the first element in the ordering, which has RANK = 0.
%    In general, the input value of RANK is increased by 1 for output,
%    unless the very last element of the ordering was input, in which
%    case the output value of RANK is 0.
%

%
%  Return the first element.
%
  if ( rank == -1 )
    p = 1 : n;
    rank = 0;
    return
  end
%
%  Check.
%
  missing = perm_check ( n, p );

  if ( missing ~= 0 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'PERM_TJ_SUCCESSOR - Fatal error!\n' );
    fprintf ( 1, '  The input array is illegal.\n' );
    fprintf ( 1, '  Missing element = ', missing );
    error ( 'PERM_TJ_SUCCESSOR - Fatal error!' );
  end

  st = 0;
  q(1:n) = p(1:n);
  done = 0;
  m = n;

  while ( 1 < m && ~done )

    d = 1;
    while ( q(d) ~= m )
      d = d + 1;
    end

    for i = d : m - 1
      q(i) = q(i+1);
    end

    par = perm_parity ( m - 1, q );

    if ( par == 1 )

      if ( d == m )
        m = m - 1;
      else
        t         = p(st+d);
        p(st+d)   = p(st+d+1);
        p(st+d+1) = t;
        done = 1;
      end

    else

      if ( d == 1 )
        m = m - 1;
        st = st + 1;
      else
        t         = p(st+d);
        p(st+d)   = p(st+d-1);
        p(st+d-1) = t;
        done = 1;
      end

    end

  end
%
%  Last element was input.  Return first one.
%
  if ( m == 1 )
    p = 1 : n;
    rank = 0;
    return
  end

  rank = rank + 1;

  return
end

author: johannesgerer

repo: jburkardt-m


text:
classdef TrueAnomalyOptimizationVariable < AbstractOptimizationVariable
    %TrueAnomalyOptimizationVariable Summary of this class goes here
    %   Detailed explanation goes here

    properties
        varObj(1,1) TrueAnomalyTermCondition = TrueAnomalyTermCondition(0);

        lb(1,1) double = 0;
        ub(1,1) double = 0;

        useTf(1,1) = false;
    end

    methods
        function obj = TrueAnomalyOptimizationVariable(varObj)
            obj.varObj = varObj;
            obj.varObj.optVar = obj;

            obj.id = rand();
        end

        function x = getXsForVariable(obj)
            x = [];

            if(obj.useTf)
                x = obj.varObj.tru;
            end
        end

        function [lb, ub] = getBndsForVariable(obj)
            lb = obj.lb(obj.useTf);
            ub = obj.ub(obj.useTf);
        end

        function [lb, ub] = getAllBndsForVariable(obj)
            lb = obj.lb;
            ub = obj.ub;
        end

        function setBndsForVariable(obj, lb, ub)
            obj.lb = lb;
            obj.ub = ub;
        end

        function useTf = getUseTfForVariable(obj)
            useTf = obj.useTf;
        end

        function setUseTfForVariable(obj, useTf)
            obj.useTf = useTf;
        end

        function updateObjWithVarValue(obj, x)
            obj.varObj.tru = x;
        end

        function nameStrs = getStrNamesOfVars(obj, evtNum, varLocType)
            nameStrs = {sprintf('Event %i True Anomaly Termination Condition', evtNum)};
        end

        function varsStoredInRad = getVarsStoredInRad(obj)
            varsStoredInRad = true;
        end
    end
end

author: Arrowstar

repo: ksptot

save_path: github-repos/MATLAB/Arrowstar-ksptot


text:
function output = rand_same_class(imdb,label)
    index = find(imdb.images.label2==label);
    output = index(randi(numel(index)));
end

author: layumi

repo: Image-Text-Embedding


text:
function [ fxx, fxy, fyy ] = f00_f2 ( fi, n, x, y )

%*****************************************************************************80
%
%% F00_F2 returns second derivatives of any function.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    01 August 2012
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer FI, the index of the function.
%
%    Input, integer N, the number of evaluation points.
%
%    Input, real X(N,1), Y(N,1), the evalution points.
%
%    Output, real FXX(N,1), FXY(N,1), FYY(N,1), the second derivatives.
%
  if ( fi == 1 )
    [ fxx, fxy, fyy ] = f01_f2 ( n, x, y );
  elseif ( fi == 2 )
    [ fxx, fxy, fyy ] = f02_f2 ( n, x, y );
  elseif ( fi == 3 )
    [ fxx, fxy, fyy ] = f03_f2 ( n, x, y );
  elseif ( fi == 4 )
    [ fxx, fxy, fyy ] = f04_f2 ( n, x, y );
  elseif ( fi == 5 )
    [ fxx, fxy, fyy ] = f05_f2 ( n, x, y );
  elseif ( fi == 6 )
    [ fxx, fxy, fyy ] = f06_f2 ( n, x, y );
  elseif ( fi == 7 )
    [ fxx, fxy, fyy ] = f07_f2 ( n, x, y );
  elseif ( fi == 8 )
    [ fxx, fxy, fyy ] = f08_f2 ( n, x, y );
  elseif ( fi == 9 )
    [ fxx, fxy, fyy ] = f09_f2 ( n, x, y );
  elseif ( fi == 10 )
    [ fxx, fxy, fyy ] = f10_f2 ( n, x, y );
  elseif ( fi == 11 )
    [ fxx, fxy, fyy ] = f11_f2 ( n, x, y );
  elseif ( fi == 12 )
    [ fxx, fxy, fyy ] = f12_f2 ( n, x, y );
  elseif ( fi == 13 )
    [ fxx, fxy, fyy ] = f13_f2 ( n, x, y );
  else
    fprintf ( 1, '\n' );
    fprintf ( 1, 'F00_F2 - Fatal error!\n' );
    fprintf ( 1, '  Illegal function index FI = %d\n', fi );
    error ( 'F00_F2 - Fatal error!' );
  end

  return
end


text:
% .^  Array power.
%    Z = X.^Y denotes element-by-element powers. X and Y must have
%    compatible sizes. In the simplest cases, they can be the same size or
%    one can be a scalar. Two inputs have compatible sizes if, for every
%    dimension, the dimension sizes of the inputs are either the same or one
%    of them is 1.
%
%    C = POWER(A,B) is called for the syntax 'A .^ B' when A or B is an
%    object.
%
%    See also MPOWER, NTHROOT, REALPOW.
%
%    Reference page in Doc Center
%       doc power
%
%    Other functions named power
%
%       codistributed/power    gpuArray/power    sym/power    ts/power
%       fints/power
%

author: jmaih


text:
function [p] = eeg_peaks(p)

% eeg_peaks - Find peaks in EEG data
%
% Usage: [p] = eeg_peaks(p)
%
% The p structure is explained in eeg_toolbox_defaults.
% For this function, it must contain the field called
% 'volt.data,' which is assumed to be a matrix
% with N rows of EEG sample values from M columns
% of electrodes.
%
% The indices to peak values are returned to
% p.volt.peaks.data (which is size(p.volt.data)).
% The values are:
%                    0   no peak
%                    1   a +ve peak
%                   -1   a -ve peak
%
% This function finds all peaks in a waveform, where
% f(x) changes from an increasing to a decreasing
% function and vice versa). These are local minima
% and maxima that have a first derivate of zero and
% satisfy the first derivative test.
%
% Thus, the +ve/-ve peak refers to peaks and troughs
% in a waveform, regardless of the +ve/-ve value of
% the waveform at that point. The peak values are
% commonly +ve/-ve potentials, but they may not be.
% A +ve peak can have a -ve value where the change
% from increasing f(x) to decreasing f(x) occurs
% entirely in the range of negative values.
%
% All the +ve/-ve peak values are returned in:
%
%   p.volt.peaks.all
%   p.volt.peaks.pos
%   p.volt.peaks.neg
%
% If the p.volt.timeArray is defined, then the timing
% of the peaks is returned in:
%
%   p.volt.peaks.alltimes
%   p.volt.peaks.postimes
%   p.volt.peaks.negtimes
%
% eg,   p.volt.timeArray = linspace(-20,20,200)';
%       p.volt.data = sin(linspace(-20,20,200))';
%       p.volt.data(:,2) = sin(-1 .* linspace(-20,20,200))';
%      [p] = eeg_peaks(p);
%       scatter(p.volt.timeArray(p.volt.peaks.all),p.volt.data(p.volt.peaks.all))
%       hold on
%       plot(p.volt.timeArray,p.volt.data)
%

% $Revision: 1.1 $ $Date: 2009-04-28 22:13:52 $

% Licence:  www.gnu.org GPL, no implied or express warranties
% Created:  02/2002, Darren.Weber_at_radiology.ucsf.edu
%                    - adapted a perl foreach loop to peaks_matrix
%                      to take advantage of matlab speed.
%           05/2002, Darren.Weber_at_radiology.ucsf.edu
%                    - used matlab 'diff' command
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

fprintf('EEG_PEAKS\n...calculating peaks...'); tic;

p.volt.peaks.data = peaks_matrix(p.volt.data);

% Find all the peak values
p.volt.peaks.all = p.volt.data .* (p.volt.peaks.data ~= 0);
p.volt.peaks.pos = p.volt.data .* (p.volt.peaks.data  > 0);
p.volt.peaks.neg = p.volt.data .* (p.volt.peaks.data  < 0);

if ~isfield(p.volt,'timeArray'),
    fprintf('\n...p doesn''t contain ''volt.timeArray''.\n');
    return
elseif isempty(p.volt.timeArray),
    fprintf('\n...p.volt.timeArray is empty.\n');
    return
end

% make sure that volt.timeArray is same size as volt.data
if ~isequal(size(p.volt.timeArray),size(p.volt.data)),
    p.volt.timeArray = repmat(p.volt.timeArray(:,1),1,size(p.volt.data,2));
end
% Find the peak timing
p.volt.peaks.alltimes = p.volt.timeArray .* (p.volt.peaks.data ~= 0);
p.volt.peaks.postimes = p.volt.timeArray .* (p.volt.peaks.data  > 0);
p.volt.peaks.negtimes = p.volt.timeArray .* (p.volt.peaks.data  < 0);

t = toc;
fprintf('done (%5.2f sec)\n',t);
return

function [peaks] = peaks_matrix(data)

    % This function should be quicker than that below.
    % The code below remains as it has been proved to
    % be reliable, but this code effectively replaces it.
    % Both code sets produce the same results.

    [r,c] = size(data);
    % replicate data, but staggered one point
    data2 = zeros(1,c);
    data2(2:r+1,:) = data;
    % add zero values to end of data points
    data(r+1,:) = zeros(1,c);
    % generate boolean difference between two matrices
    dif = data < data2;

    % Now replicate this boolean matrix
    [br,bc] = size(dif);
    dif2 = zeros(1,bc);
    dif2(2:br+1,:) = dif;
    % add zeros to end of boolean dif
    dif(br+1,:) = 0;

    % Calculate the peaks
    peaks = dif(2:r+1,:) - dif2(2:r+1,:);

    % now trim peaks dif back to size of data
    peaks(1,:) = 0; % first point cannot be a peak
    peaks(r,:) = 0; % last point cannot be a peak

    % Could use matlab diff command, like this
    %[r,c] = size(data);
    %decdat = (diff(data)<=0) .* -1;
    %incdat =  diff(data)>=0;

    %zrow = zeros(1,c);
    %dirdat1 = [ zrow; decdat + incdat ];
    %dirdat2 = [ decdat + incdat; zrow ];

    %negpeaks = ((dirdat1 - dirdat2) < 0) .* -1;
    %pospeaks =  (dirdat1 - dirdat2) > 0;
    %peaks = negpeaks + pospeaks;

    %peaks(1,:)   = 0; % first point cannot be a peak
    %peaks(end,:) = 0; % last point cannot be a peak

return

% The function 'peaks_matrix' above replaced the
% following for loop (02/2002, Darren Weber)

% foreach electrode
%for e = 1:size(p.volt.data,2),
%	v = p.volt.data(:,e);
%    n = 1;
%    volt = v(n);
%	while (n <= size(v,1)),
%        if(volt < v(n)),
%            volt = v(n); n = n + 1;
%            if(volt < v(n)), p.volt.peaks.data(n-1,e) = 0;
%            else             p.volt.peaks.data(n-1,e) = 1;
%            end
%        elseif(volt > v(n)),
%            volt = v(n); n = n + 1;
%            if(volt > v(n)), p.volt.peaks.data(n-1,e) = 0;
%            else             p.volt.peaks.data(n-1,e) = -1;
%            end
%        else
%            volt = v(n); n = n + 1;
%           p.volt.peaks.data(n-1,e) = 0;
%        end
%        if isequal(n,size(v,1)),
%            break;
%        end
%    end
%    p.volt.peaks.data(1,e) = 0; % first point cannot be a peak
%    p.volt.peaks.data(n,e) = 0; % last point cannot be a peak
%end
%p.volt.peaks.all = find(p.volt.peaks.data ~= 0);
%p.volt.peaks.pos = find(p.volt.peaks.data  > 0);
%p.volt.peaks.neg = find(p.volt.peaks.data  < 0);
%return

author: PatternRecognition

repo: OpenBMI


text:
function r = tanh(a)
%TANH          Taylor hyperbolic tangent  tanh(a)
%
%Thanks to George Corliss for providing the Taylor expansion
%

% written  06/03/09     S.M. Rump
% modified 08/26/12     S.M. Rump  global variables removed
%

  e = 1e-30;
  if 1+e==1-e                   % fast check for rounding to nearest
    rndold = 0;
  else
    rndold = getround;
    setround(0)
  end

  K = getappdata(0,'INTLAB_TAYLOR_ORDER');

  r = a;
  ct = a.t;
  N = size(a.t,2);
  r.t(1,:) = tanh(a.t(1,:));
  ct(1,:) = 1-r.t(1,:).^2;     % 1-tanh(a)^2
  r.t(2,:) = ct(1,:) .* a.t(2,:);
  for j=2:K
    ct(j,:) = -sum( r.t(1:j,:).*r.t(j:-1:1,:) , 1 );
    r.t(j+1,:) = sum( ct(1:j,:).*a.t(j+1:-1:2,:).*repmat((j:-1:1)',1,N) , 1 )/j;
  end

  if rndold
    setround(rndold)
  end

repo: OpenMAS


text:
classdef TestGraphSegmentation
    %TestGraphSegmentation

    methods (Static)
        function test_1
            img = imread(fullfile(mexopencv.root(),'test','fruits.jpg'));
            gs = cv.GraphSegmentation('K',300);
            assert(isequal(gs.K,300));
            L = gs.processImage(img);
            validateattributes(L, {'int32'}, ...
                {'size',[size(img,1) size(img,2)]});
        end
    end

end

author: kyamagu

repo: mexopencv

save_path: github-repos/MATLAB/kyamagu-mexopencv


text:
function [varargout] = likGumbel(sign, hyp, y, mu, s2, inf, i)

% likGumbel - Gumbel likelihood function for extremal value regression.
% The expression for the likelihood is
%   likGumbel(t) = exp(-z-exp(-z))/be, z = ga+s*(y-t)/be, be = sn*sqrt(6)/pi
% where s={+1,-1} is a sign switching between left and right skewed, ga is the
% Euler-Mascheroni constant, y is the mean, sn^2 is the variance.
% The skewness and kurtosis of likGumbel are 1.14*s and 2.4, respectively.
%
% The hyperparameters are:
%
% hyp = [ log(sn)  ]
%
% Several modes are provided, for computing likelihoods, derivatives and moments
% respectively, see likFunctions.m for the details. In general, care is taken
% to avoid numerical issues when the arguments are extreme.
%
% Copyright (c) by Hannes Nickisch, 2013-11-01.
%
% See also LIKFUNCTIONS.M.

if nargin<4, varargout = {'1'}; return; end   % report number of hyperparameters
if sign=='-', s = -1; else s = 1; end                 % extract sign of skewness

sn2 = exp(2*hyp);                                      % extract hyperparameters
ga = 0.5772156649;                                   % Euler-Mascheroni constant
be = sqrt(6*sn2)/pi;
lZ = -log(be);

if nargin<6                              % prediction mode if inf is not present
  if numel(y)==0,  y = zeros(size(mu)); end
  s2zero = 1; if nargin>4, if norm(s2)>0, s2zero = 0; end, end         % s2==0 ?
  if s2zero                                         % log probability evaluation
    lp = likGumbel(sign, hyp, y, mu, [], 'infLaplace'); s2 = 0;
  else                                                              % prediction
    lp = likGumbel(sign, hyp, y, mu, s2, 'infEP');
  end
  ymu = {}; ys2 = {};
  if nargout>1
    ymu = mu;                                                   % first y moment
    if nargout>2
      ys2 = s2 + sn2;                                          % second y moment
    end
  end
  varargout = {lp,ymu,ys2};
else
  switch inf
  case 'infLaplace'
    z = ga+s*(y-mu)/be; emz = exp(-z);
    if nargin<7                                             % no derivative mode
      dlp = {}; d2lp = {}; d3lp = {};
      lp = lZ -z -emz;
      if nargout>1
        dz = -s/be;                                                     % dz/dmu
        dlp = dz*(emz-1);                    % dlp, derivative of log likelihood
        if nargout>2                    % d2lp, 2nd derivative of log likelihood
          d2lp = -dz^2*emz;
          if nargout>3                  % d3lp, 3rd derivative of log likelihood
            d3lp = dz^3*emz;
          end
        end
      end
      varargout = {lp,dlp,d2lp,d3lp};
    else                                             % derivative w.r.t. log(sn)
      dz = -s/be;                                                       % dz/dmu
      dzs = -s*(y-mu)/be;                                          % dz/dlog(sn)
      lp_dhyp   =  dzs.*(emz-1) -1;
      dlp_dhyp  = dz*(1-emz.*(1+dzs));
      d2lp_dhyp = dz^2*emz.*(2+dzs);
      varargout = {lp_dhyp,dlp_dhyp,d2lp_dhyp};
    end

  case 'infEP'
    if nargout>1
      error('infEP not supported since likT is not log-concave')
    end
    n = max([length(y),length(mu),length(s2)]); on = ones(n,1);
    y = y(:).*on; mu = mu(:).*on; sig = sqrt(s2(:)).*on;          % vectors only
    % since we are not aware of an analytical expression of the integral,
    % we use Gaussian-Hermite quadrature
    N = 20; [t,w] = gauher(N); oN = ones(1,N);
    lZ = likGumbel(sign, hyp, y*oN, sig*t'+mu*oN, []);
    lZ = log_expA_x(lZ,w); % log( exp(lZ)*w )
    varargout = {lZ};

  case 'infVB'
    error('infVB not supported')
  end
end

%  computes y = log( exp(A)*x ) in a numerically safe way by subtracting the
%  maximal value in each row to avoid cancelation after taking the exp
function y = log_expA_x(A,x)
  N = size(A,2);  maxA = max(A,[],2);      % number of columns, max over columns
  y = log(exp(A-maxA*ones(1,N))*x) + maxA;  % exp(A) = exp(A-max(A))*exp(max(A))

author: benfulcher

repo: hctsa

save_path: github-repos/MATLAB/benfulcher-hctsa


text:
function bound = AnglesLocalMaxMin(f,N)

%================================================================
% function bound = AnglesLocalMaxMin(f,N)
%
% This function segments f into a maximum of N supports by taking
% the middle point between the N largest local maxima.
% Note: the detected boundaries are given in term of indices
%
% Inputs:
%   -f: the function to segment
%   -N: maximal number of bands
%
% Outputs:
%   -bound: list of detected boundaries
%
% Author: Jerome Gilles
% Institution: UCLA - Department of Mathematics
% Year: 2013
% Version: 1.0
%===============================================================

locmax=zeros(size(f));
locmin=max(f)*ones(size(f));
% We detect the local maxima
for i=2:length(f)-1
    if ((f(i-1)<f(i)) && (f(i)>f(i+1)))
        locmax(i)=f(i);
    end

    if ((f(i-1)>f(i)) && (f(i)<=f(i+1)))
        locmin(i)=f(i);
    end
end

% We check if the endpoint are local maxima or minima (we work on the torus)
if ((f(end)<f(1)) && (f(1)>f(2)))
    locmax(1)=f(1);
end
if ((f(end)>f(1)) && (f(1)<=f(2)))
    locmin(1)=f(1);
end

if ((f(end-1)<f(end)) && (f(end)>f(1)))
    locmax(end)=f(end);
end
if ((f(end-1)>f(end)) && (f(end)<=f(1)))
    locmin(end)=f(end);
end

% We keep the N-th highest maxima and their index
[lmax,Imax]=sort(locmax,1,'descend');
if length(lmax)>N
    Imax=sort(Imax(1:N));
else
    Imax=sort(Imax);
    N=length(lmax);
end

% We detect the lowest minima between two consecutive maxima
bound=zeros(1,N);
for i=1:N
   if i==N
       [lmin,ind]=sort([locmin(Imax(i):end);locmin(1:Imax(1))]);
       tmp=lmin(1);
       n=1;
       if n<length(lmin)
            n=2;
            while ((n<=length(lmin)) && (tmp==lmin(n)))
                n=n+1;
            end
       end
        bound(i)=Imax(i)+ind(ceil(n/2))-1;
       if bound(i)>length(f)
           bound(i)=bound(i)-length(f);
       end
   else
       [lmin,ind]=sort(locmin(Imax(i):Imax(i+1)));
       tmp=lmin(1);
       n=1;
       if n<length(lmin)
            n=2;
            while ((n<=length(lmin)) && (tmp==lmin(n)))
                n=n+1;
            end
       end
       bound(i)=Imax(i)+ind(ceil(n/2))-1;
   end
end

bound=sort(bound);

author: Sable

repo: mcbench-benchmarks