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

function [tcArrayCell, responseIndArrayCell, rSquareArrayCell, fitParamsCell, ...
          flyIndsCell, tcExtCellArray, tcArgExtCellArray, varargout] = getStimArraysForSameRois(StimTableCell, trimArray)

if nargin < 2
    trimArray = 1;
end    
[tcArrayCell, responseIndArrayCell, rSquareArrayCell, fitParamsCell, ...
 tcExtCellArray, tcArgExtCellArray] = cellfun(@getTuningCurveArrayFromTable, ...
                                              StimTableCell, 'uni', 0);

nStims = numel(StimTableCell);

for iStim = 1: nStims
    flyIndsCell{iStim} = [];                       
    for jFly = 1: size(StimTableCell{iStim})
        nRoisTemp = numel(StimTableCell{iStim}(jFly, 'responseIndex').responseIndex{1});
        flyIndsCell{iStim} = [flyIndsCell{iStim} repelem(StimTableCell{iStim}{jFly, 1}, nRoisTemp)];
    end
end
% Big ugly chunk to get rois that have seen all stimuli (bars in this case)
if trimArray
    % find mismatches of number of ROIs accross stimuli
    nRoisCell = cellfun(@(y) cellfun(@numel, y.responseIndex), ...
                        StimTableCell, 'uni', 0);
    nRoisArray = [nRoisCell{:}];
%     % use the fitting parameters, because using only response quality neglects
%     % the fitting failures.
%     nRoisCell = cellfun(@(y) cellfun(@numel, y.paramTuning), ...
%                         StimTableCell, 'uni', 0);
%     nRoisArray = cat(1, nRoisCell{:});
    for iFly = 1: size(nRoisArray, 1)
        equalNoCells(iFly) = all(nRoisArray(iFly, :) == nRoisArray(iFly, 1));
    end
    % Check if there are any stimuli with distinct number of rois.
    if any(~equalNoCells)
        for iStim = 1: nStims
            roisMeta = StimTableCell{iStim}.roiMeta;
            bgRois = cellfun(@(x) x.backgroundRois, roisMeta);
            nValidRois = cellfun(@numel, StimTableCell{iStim}.responseIndex);
            % Last roi picked is the background.
            roiInds = arrayfun(@(x) 1: x, bgRois, 'uni', 0);    
            bgOrInvalidRois = cellfun(@(x) [x.backgroundRois x.invalidRois], ...
                                      roisMeta, 'uni', 0);
            validRoiInds{iStim} = cellfun(@(x, y) setdiff(x, y), ...
                                          roiInds, bgOrInvalidRois, 'uni', 0);
            nFlies = size(StimTableCell{iStim}, 1);
        end
        validRoiIndsCellArray = [validRoiInds{:}];
        % Check what rois to remove only for flies with different roi
        % number.
        for iFly = find(~equalNoCells)
            if nStims == 4
                offInds = intersect(validRoiIndsCellArray{iFly, 1: 2});
                onInds = intersect(validRoiIndsCellArray{iFly, 3: 4});
                allInds = intersect(offInds, onInds);
            else
                allInds = validRoiIndsCellArray{iFly, 1};
                for iStim = 2: nStims
                    allInds = intersect(allInds, validRoiIndsCellArray{iFly, iStim});
                end
            end
            for jStim = 1: nStims
                % Find roi index of the ones not in the intersection.
                validInds = validRoiIndsCellArray{iFly, jStim};
                invalidRoisInds = setdiff(validInds, allInds); 
                tempInd = arrayfun(@(x) find(validInds == x), ...
                                   invalidRoisInds, 'uni', 0);
                roisToRemoveInd{iFly, jStim} = tempInd;
            end
        end
        for jStim = find(~all(cellfun(@isempty, roisToRemoveInd), 1))
            for iFly = find(~equalNoCells) 
            % Quick fixup for only one invalid roi per fly
                if ~isempty(roisToRemoveInd{iFly, jStim})
                    deleteInds = arrayfun(@(x) find(validRoiIndsCellArray{iFly, jStim} == x), ...
                                          roisToRemoveInd{iFly, jStim}{1}, 'uni', 1);
                    validRoiIndsCellArray{iFly, jStim}(deleteInds) = 0;
                end
            end
        end
        transposedValidRois = cellfun(@(x) x', validRoiIndsCellArray, 'uni', 0);
        vectorValidRois = arrayfun(@(x) cat(1, transposedValidRois{:, x}), 1: size(validRoiInds, 2), 'uni', 0);
        for jStim = find(~all(cellfun(@isempty, roisToRemoveInd), 1))
            tempInd = vectorValidRois{jStim} == 0;
            responseIndArrayCell{jStim}(tempInd, :) = [];
            rSquareArrayCell{jStim}(tempInd) = [];
            fitParamsCell{jStim}.a1(tempInd) = [];
            fitParamsCell{jStim}.b1(tempInd) = [];
            fitParamsCell{jStim}.c1(tempInd) = [];
            flyIndsCell{jStim}(tempInd) = [];
            if ~isempty(tcArrayCell{jStim})
                tcArrayCell{jStim}(tempInd, :) = [];
    %           tcExtCellArray{jStim}(tempInd, :) = [];
    %           tcArgExtCellArray{jStim}(tempInd, :) = [];
                tcExtCellArray{jStim}(:, tempInd) = [];
                tcArgExtCellArray{jStim}(:, tempInd) = [];
            end
        end
                
    else
        roisToRemoveInd = {};
    end
    varargout{1} = roisToRemoveInd;
    % Check that all stimuli have same number of rois, it could be by
    % chance, but here we assume the cells all saw same stimulus. Actually
    % it should have an unique roisID but so far its not implemented.
    sizeCell = cellfun(@(x) size(x, 1), tcArrayCell, 'uni', 0);
    assert(isequal(sizeCell{:}));
end

end