function [apc_idx, apc_symidx, apc_emesh] = run_apc(apc_in, varargin)

    BIN_APC='~/bin/apcluster';
    CLUS_DIST='correlation';
    force_data = 0;
    verbose = 0;
    run_matlab = 1;
    run_noclean = 0;
    apc_dist = [];
    apc_nondefaults = 0;
    apc_maxits = 500;
    apc_convits = 50;
    apc_dampfact = 0.5;
    sym_idx = [];
    sym_text = [];
    apc_symidx = [];
    emesh = [];
    filename = '_tmp_apc_runapc';
    
    
    if ~isempty(varargin)
        i = 1;
        while i <= length(varargin),
            switch(varargin{i}),
                case 'force_data'
                    force_data = 1;
                case 'matlab'
                    run_matlab = 1;
                case 'unix'
                    run_matlab = 0;
                case 'nocleanup'
                    run_noclean = 1;
                case 'maxits'
                    i = i+1;
                    apc_maxits = varargin{i};
                    apc_nondefaults = 1;
                case 'convits'
                    i = i+1;
                    apc_convits = varargin{i};
                    apc_nondefaults = 1;
                case 'dampfact'
                    i = i+1;
                    apc_dampfact = varargin{i};
                    apc_nondefaults = 1;
                case 'file'
                    i = i+1;
                    filename = varargin{i};
                case 'symidx'
                    i = i+1;
                    sym_idx = varargin{i};
                case 'symtext'
                    i = i+1;
                    sym_text = varargin{i};
                case 'show'
                    i = i+1;
                    emesh = varargin{i};
                case 'distance'
                    i = i+1;
                    CLUS_DIST = varargin{i};
                case 'verbose'
                    i = i+1;
                    verbose = varargin{i};
                 otherwise
                    error('Argument %s not known', varargin{i});
            end
            i = i+1;
        end
    end
    
    if size(apc_in,1) == size(apc_in, 2),
        if force_data == 0,
            if verbose > 0, fprintf('Input is assumed to be distance\n'); end
            apc_dist = apc_in;
        end
    end
    
    if isempty(apc_dist)
        
        if verbose > 0, fprintf('Calculating ''%s'' distance ... ', CLUS_DIST); end
        
        apc_data = apc_in;
    
        if strcmp(CLUS_DIST, 'euclidian'),
            apc_dist = squareform(pdist(normalize(apc_data')));
        else
            apc_dist = squareform(pdist(apc_data', CLUS_DIST));
        end

        if verbose > 0, fprintf('done\n'); end
    end

    N = size(apc_dist,1);

    if run_matlab == 1,

        if verbose > 0, fprintf('Running matlab code\n'); end
        
        s = zeros(N*N-N,3); 
    
        j=1;
        for i=1:N
            for k=[1:i-1,i+1:N]
                s(j,1)=i; s(j,2)=k; s(j,3)=-apc_dist(i,k);
                j=j+1;
            end;
        end

        if verbose > 0, fprintf('Clustering ... '); end

        apc_p = median(s(:,3));
        if apc_nondefaults == 0,
            apc_idx = apcluster(s, apc_p);
        else
            apc_idx = apcluster(s, apc_p, ...
                    'maxits', apc_maxits, 'convits', apc_convits, 'dampfact', apc_dampfact);
        end

        if verbose > 0, fprintf('done\n'); end
    else
        if verbose > 0, fprintf('Running external code'); end
        
        if apc_nondefaults == 0,
            apc_extopt = '';
        else
            apc_extopt = sprintf('%d %d %d', apc_maxits, apc_convits, apc_dampfact);
        end
        
         
        if verbose > 0, fprintf('Exporting to files with names starting with ''%s'' ...', filename); end
        s = zeros(N*N-N,1); 

        fn_s = sprintf('%s_s.txt', filename);

        fid = fopen(fn_s, 'wt');
        j=1;
        for i= 1:N,
            k=[1:i-1,i+1:N];
            fprintf(fid, '%d %d %f\n', [repmat(i, [1, N-1]); k;  -apc_dist(i,k)]);
            s(j) = -apc_dist(i,k);
            j=j+1;
        end
        fclose(fid);
        
        apc_p = median(s);

        fn_p = sprintf('%s_p.txt', filename);
        fid = fopen(fn_p, 'wt');
        for i= 1:N,
            fprintf(fid, '%f\n', repmat(apc_p, [1, N-1]));
        end
        fclose(fid);

 
        apc_bin = sprintf('%s %s %s %s %s', BIN_APC, fn_s, fn_p, fn_out, apc_extopt);

        if verbose > 0, fprintf('done\nRunning ''%s''\n', apc_bin); end

        status = 1;
        % status = unix(apc_bin);
        
        if status == 0,
            fid=fopen(fn_out, 'r');
            apc_idx = fscanf(fid, '%d', [1 inf]);
            fclose(fid);
        else
            apc_idx = [];
            error('External program had status %d, nothing returned', status);
        end

        
    end

    if ~isempty(sym_idx), 
        if ~isempty(sym_text),
            for i=1:length(sym_idx),
                apc_symidx{i} = sym_text{sym_idx(i)};
            end
        else
            apc_symidx = sym_idx(apc_idx);    
        end
    end
    
    
    apc_emesh = [];
    if ~isempty(emesh) && ~isempty(apc_data),
        
        apc_emesh.stain = apc_data;
        apc_emesh.dmesh = [];
        apc_emesh.idx = 1:size(apc_data,2);
        % st4tmesh_ref(st4tpat,1) shows the indices of the genes in st4tmesh_pat corresponding to st4usym
        apc_emesh.symidx = sym_idx;
        for i=1:length(sym_idx), 
            apc_emesh.sym{i} = sym_text{sym_idx(i)}; 
        end
        
        study_clusters(emesh, [], apc_emesh, apc_idx, -1,'display','display_method', 'svd');
    end
    
    
end