function [dX,dY] = pinvUpdate(Value,Mask,Value_,varargin)
%PINVUPDATE Summary of this function goes here
%   Detailed explanation goes here

p = inputParser; p.KeepUnmatched = true; p.StructExpand = true;
p.addParameter('constrained', false, @islogical);

p.parse(varargin{:});
constrained = p.Results.constrained;

%% Constants
projSize = size(Value);
uref = floor(projSize(1)/2)+1;
vref = floor(projSize(2)/2)+1;

seAng = strel('arbitrary',[0,1,0]|[0,1,0]');
seSpa = strel('arbitrary',reshape([0,1,0]|[0,1,0]',1,1,3,3));

[u,v] = ndgrid(1:projSize(1),1:projSize(2));
u = u-uref;
v = v-vref;

%% Compute difference
Delta = Value_-Value;
DeltaMat = LF.LFToMat(Delta);

%% Extend signal (to avoid cropping the gradient at the boundaries)
ValDil = extend(Value,seSpa,2);
ValDil = extend(ValDil,seAng,2);

%% Compute gradient
[LFx,LFy] = LF.gradient(ValDil);
LFx(~Mask) = nan;
LFy(~Mask) = nan;

%% Compute disparity update
if constrained
    %% Compute jacobian and delta
    Fx = u.*LFx;
    Fy = v.*LFy;
    
    Fx = Fx(Mask);
    Fy = Fy(Mask);
    
    J = [Fx(:),Fy(:)];
    
    Delta = Delta(Mask);
    Delta = Delta(:);

%% Compute disparity update
    dXdY = pinv(J)*Delta;
    
    dX = dXdY(1);
    dY = dXdY(2);
else
    %% Compute jacobian
    LFxMat = LF.LFToMat(LFx);
    LFyMat = LF.LFToMat(LFy);
    
    LFxMat = LF.split(LFxMat,[0,1]);
    LFyMat = LF.split(LFyMat,[0,1]);
    
    LFxMat = cellfun(@(x) x(~isnan(x)),LFxMat,'UniformOutput',false);
    LFyMat = cellfun(@(x) x(~isnan(x)),LFyMat,'UniformOutput',false);
    
    JMat = cellfun(@(x,y) [x,y],LFxMat,LFyMat,'UniformOutput',false);
    
    DeltaMat = LF.split(DeltaMat,[0,1]);
    DeltaMat = cellfun(@(x) x(~isnan(x)),DeltaMat,'UniformOutput',false);
    
    %% Compute disparity update
    dXdY = cellfun(@(j,d) pinv(j)*d,JMat,DeltaMat,'UniformOutput',false);
    dXdY = reshape(dXdY,projSize(1:2));
    
    [dX,dY] = cellfun(@(dxdy) deal(dxdy(1),dxdy(2)),dXdY);
end
end

function Value = extend(Value,se,numIter)
for i=1:numIter
        Mask = ~isnan(Value);
        Temp = Value;
        Temp(~Mask)=-Inf;
        Temp = imdilate(Temp,se);
        Temp(isinf(Temp))=NaN;
        Value(~Mask) = Temp(~Mask);
end
end