off-motion-receptive-fields/receptiveFields/getAlignedTuningCurveArrayFromTable.m

function [alignedRFsCell, ...
          alignedRFsCellNaN, varargout] = getAlignedTuningCurveArrayFromTable(StimTableCell, ...
                                                                   thresholdResponseInd, ...
                                                                   thresholdRSquare, ...
                                                                   varargin)
%% Input should be a cell array of tables for XOFF, YOFF, XON, YON.
switch nargin
    case 3
        alignToOff = false;
        alignToAll = false;
    case 4
        alignToOff = varargin{1};
        aignToAll = false;
    case 5
        alignToOff = varargin{1};
        alignToAll = varargin{2};        
    otherwise
        warning('Too many arguments, cannot process them meaningfully!');
end

%%
nStims = numel(StimTableCell);
% Retrieve table data in arrays.
 [tcArrayCell, responseIndArrayCell, rSquareArrayCell, fitParamsCell, ...
  flyIndsCell, tcExtCellArray, tcArgExtCellArray] = getStimArraysForSameRois(StimTableCell, alignToOff);
%% Set thresholds.
isAboveThresholdCell = cellfun(@(x, y) ...
                               x(:) > thresholdResponseInd & ...
                               y(:) > thresholdRSquare, ...
                               responseIndArrayCell, rSquareArrayCell, 'uni', 0); 

if alignToOff
    offThreshold = isAboveThresholdCell{1} & isAboveThresholdCell{2};
    onThreshold = isAboveThresholdCell{3} & isAboveThresholdCell{4};
    if alignToAll
        [isAboveThresholdCell{:}] = deal(offThreshold & onThreshold);
        [current2OffInds{1:4}] = deal(1: sum(isAboveThresholdCell{1}));
    else
        [isAboveThresholdCell{1:2}] = deal(offThreshold);
        isAboveThresholdCell{3} = offThreshold & isAboveThresholdCell{3};
        isAboveThresholdCell{4} = offThreshold & isAboveThresholdCell{4};    
        offInds = find(offThreshold);
        [current2OffInds{1:2}] = deal(offInds);
        current2OffInds{3} = arrayfun(@(x) find(x == offInds), ...
                                      find(isAboveThresholdCell{3}));
        current2OffInds{4} = arrayfun(@(x) find(x == offInds), ...
                                      find(isAboveThresholdCell{4}));
    end
end
% Fill in in case no cell passed the threshold.
emptyStimulusInd = cellfun(@(x) ~any(x), isAboveThresholdCell);   
if all(emptyStimulusInd)
    fprintf('No cell fulfill requested conditions!!! Change thresholds \n')
    return
end
firstNonEmptyStim = isAboveThresholdCell{find(~emptyStimulusInd, 1)};
for iStim = find(emptyStimulusInd)
    isAboveThresholdCell{iStim} = ones(sum(firstNonEmptyStim), 1);
end

% Now get the thresholded arrays. We still need the case for all stimuli
% passing the threshold.
tcArrayThres = cellfun(@(x, y) x(y, :), ...
                       tcArrayCell, isAboveThresholdCell, 'uni', 0);
respIndArrayThres = cellfun(@(x, y) x(y), ...
                            responseIndArrayCell, isAboveThresholdCell, 'uni', 0);
rSquareArrayThres = cellfun(@(x, y) x(y), ...
                            rSquareArrayCell, isAboveThresholdCell, 'uni', 0);
a1Thres = cellfun(@(x, y) x.a1(y), ...
                  fitParamsCell, isAboveThresholdCell, 'uni', 0);
b1Thres = cellfun(@(x, y) x.b1(y), ...
                  fitParamsCell, isAboveThresholdCell, 'uni', 0);
c1Thres = cellfun(@(x, y) x.c1(y), ...
                  fitParamsCell, isAboveThresholdCell, 'uni', 0);
flyIndsThres = cellfun(@(x, y) x(y), ...
                  flyIndsCell, isAboveThresholdCell, 'uni', 0);
tcExtArrayThres = cellfun(@(x, y) x(y), ...
                          tcExtCellArray, isAboveThresholdCell, 'uni', 0);      
tcArgExtArrayThres = cellfun(@(x, y) x(y), ...
                          tcArgExtCellArray, isAboveThresholdCell, 'uni', 0);  
padWithNaNs = 01;
polarityCell = cellfun(@getStimPolarity, StimTableCell, 'uni', 0);
%% Align the thresholded tuning curves.

for iStim = 1: nStims
    if alignToOff && iStim > 2
        lagsToUse = lagDiffsCell{iStim - 2}(current2OffInds{iStim});
    else
        lagsToUse = [];
    end
    [alignedRFsCell{iStim}, lagDiffsCell{iStim}] = alignRFsToAbsMax(...
                                                        tcArrayThres{iStim}, ...
                                                        padWithNaNs * 0, ...
                                                        lagsToUse, ...
                                                        polarityCell{iStim});
    % Set to zero all array for stimuli without any rois aboive threshold.
    if ~any(responseIndArrayCell{iStim} > thresholdResponseInd & ...
            rSquareArrayCell{iStim} > thresholdRSquare)
        alignedRFsCell{iStim} = zeros(size(tcArrayCell{iStim}));
    end
end
% This is just to get the nans.
for iStim = 1: nStims
    if alignToOff && iStim > 2
        lagsToUse = lagDiffsCell{iStim - 2}(current2OffInds{iStim});
    else
        lagsToUse = [];
    end
    [alignedRFsCellNaN{iStim}, lagDiffsCell{iStim}] = alignRFsToAbsMax(...
                                                        tcArrayThres{iStim}, ...
                                                        padWithNaNs, ...
                                                        lagsToUse, ...
                                                        polarityCell{iStim});
    % Set to zero all array for stimuli without any rois aboive threshold.
    if ~any(responseIndArrayCell{iStim} > thresholdResponseInd & ...
            rSquareArrayCell{iStim} > thresholdRSquare)
        alignedRFsCellNaN{iStim} = zeros(size(tcArrayCell{iStim}));
    end
end

fitParamsThres.a1 = a1Thres;
fitParamsThres.b1 = b1Thres;
fitParamsThres.c1 = c1Thres;
varargout{1} = fitParamsThres;
varargout{2} = flyIndsThres;
varargout{3} = isAboveThresholdCell;
varargout{4} = tcExtArrayThres;
varargout{5} = tcArgExtArrayThres;
varargout{6} = respIndArrayThres;
varargout{7} = rSquareArrayThres;

end


function polarity = getStimPolarity(StimTable)

% Get stim polarity.
tokens = regexp(StimTable(1, :).stimParamFileName{1}, '.*_([mp]).*Con_.*', 'tokens');
switch tokens{1}{1}
    case 'm'
        polarity = 'off';
    case 'p'
        polarity = 'on';
    otherwise
        polarity = 'other';
end

end