function [dataUnd, matrixUnd] = ApplyUndersampling(kp,data,binning)
% [dataUnd, matrixUnd] = ApplyUndersampling(kp,data,binning)
%
% Inputs: 
% kp        = 192x200, undersampling pattern for the 200 temporal points
% data      = 192x192x200x4 before binning for 200 temporal points, and 4
% coils 
% binning   = 192*200x1, vector that assigns each encoding line to a
% certain frame 
%
% Outputs:  
% dataUnd   = 192x192x8x4 undersampled data after binning into 8 frames
% for 4 coils 
% matrixUnd = 192x192x8, undersampling pattern after binning

N               =  size(data);
aux             = N(1)/2;
NI              = 200;
frames          = 8;

dataUnd         = zeros(N(1),N(2),frames,4);
count           = zeros(N(1),frames);

%Apply undersampling and data binning
for j=1:NI
    for k=1:N(1)
        index= k+(j-1)*N(1);
        a=binning(index);
        if (a ~= 0)
            if (k <= aux)
                pos=2*k-1;
            else
                pos=N(1)-2*(k-1-aux);
            end
            b=kp(pos,j);
            if (b~= 0)
                count(pos,a)=count(pos,a)+1;
                dataUnd(pos,:,a,:)=dataUnd(pos,:,a,:)+data(pos,:,j,:);
            end
        end
    end
end

% Prepare data and undersampling matrix to be reconstructed
matrixUnd=zeros(N(1),N(2),frames);
matrixAux=count>0;

for paux=1:N(1)
    matrixUnd(:,paux,:)=matrixAux(:,:);
end

for j=1:frames
    for k=1:N(1)
        if (count(k,j) ~= 0)
            dataUnd(k,:,j,:)=dataUnd(k,:,j,:)./count(k,j);
        end
    end
end

for p=1:frames
    matrixUnd(:,:,p)=ifftshift(matrixUnd(:,:,p));
    for a=1:4
        dataUnd(:,:,p,a)=ifftshift(dataUnd(:,:,p,a));
    end
end


end