function [xnew, S] = BP_LOT(x,F,G,b,s,band)
%%% BLOT



[m,n] = size(F); % = size(G)

S = [];
Region = 1:n;


q = 2;

for iter=1:s
    % Step1
    xk = x(1:n);
    xkp = x(n+1:2*n);
    abs_x_band = abs(xk(Region)) + abs(xkp(Region))/(2*pi* q);    
    [temp,lambda] = max(abs_x_band);
    NewIndex = Region(lambda);
    Region = setdiff(Region,NewIndex-band:NewIndex+band);
    
    % Step2
    S = [S; NewIndex];
   
end
%S
%S = LO([F,G],b,S,band);
%S

%A = [F,G];
A_S = [F(:,S), G(:,S)];

xnew = zeros(2*n,1);
temp = lscov(A_S,b);

xnew(S) = temp(1:s);
xnew(S+n) = temp(s+1:2*s);
%x = (A_S'*A_S)\(A_S'*b);
    
%%

function S = LO(A,b,S0,eta)
%%% Local Optimization
display('ya')
S = sort(S0);

min(diff(S))

for k=1:length(S0)    
    ik = S(k);
    S = setdiff(S,ik);

    res = 1e20;
    for j=max([ik-eta,1]):min([ik+eta,size(A,2)])

        A_S = A(:,[S;j]);
        x = lscov(A_S,b);
        r = b - A_S*x;
        
        if(norm(r,2) < res)
            res = norm(r,2);
            j0 = j;
        end
    end
    
    S = sort([S;j0]);   
end



%%





