function mvClusteringOnly(seq_in, mv_path, seg_path, pC, pR, nFrm, REJRATE, MAXITER)

%% input file
[fid_seq_in message]= fopen(seq_in,'rb');

%% parameters of a synthetic motion field
%
blkSiz = 8; % block size
bC=pC/blkSiz; % number of column in blocks
bR=pR/blkSiz; % number of row in blocks

% get the coordinates
[coorBlkY,coorBlkX]=ndgrid(1:bR,1:bC);
coorX = coorBlkX.*blkSiz-blkSiz/2;
coorY = coorBlkY.*blkSiz-blkSiz/2;

HALFPIX=2; % Half pixel

for ii = 1:nFrm
    % get one frame from YUV sequence (YUV: 411 format) 
    [yCurr, uCurr, vCurr] = readOneFrame(fid_seq_in, pR, pC);
    if ii == 1
        iniH=[];
        fgMap = zeros(bR, bC);
    else
        %% get input noisy motion field
        [npx, npy] = readMVs(mv_path,ii, bC, bR, HALFPIX);

        %% combined GME and segmentation
        [fgMap, MM, iniH] = mrfInGme(npx, npy, coorX, coorY, ...
            blkSiz, iniH, fgMap, REJRATE, MAXITER, ii, seg_path, yCurr, uCurr, vCurr);
    
    end
end

fclose(fid_seq_in);



function [fgMap, MM, iniH] = mrfInGme(npx, npy, coorX, coorY, ...
    blkSiz, iniH, fgMap, REJRATE, nITER, IDFRM, SEG_PATH, img_y, img_u, img_v)

REM_MODE = 3;

GM_TRAN = 1;  % translational model
GM_ISOT = 2;  % isotripic model
GM_AFFI = 3;  % affine model
GM_PERS = 4;  % perspective model
REM_MODE = 3;
MRF_ENA = 1;

PO1 = 0.05;      % type-1 outlier --> 10% of MVF
PO2 = 0.40;      % type-2 outleir --> 20% of MVF
PO = PO1 + PO2;  % the overall outleirs --> 30% of MVF

TYPE1PO = 0.15;
type1Map = OutlierRem(npx, npy, TYPE1PO, PO2, fgMap, blkSiz, 1);

for ii=1:nITER
   
    % iteratively refine GME
    if ii == 1
        % GME based on all inliers using Newton-Raphson method  
        if IDFRM == 2  % for the first frame, we use more iterations
            % the MV outliers removal cascade
            iMap = OutlierRem(npx, npy, PO1, PO2, fgMap, blkSiz, 3);
            [MM, iniH] = gmeLs(coorX, coorY, npx, npy, iniH, iMap, REJRATE, 0.001);
        else
            % the MV outliers removal cascade
            iMap = 1-fgMap;
            [MM, iniH] = gmeLs(coorX, coorY, npx, npy, iniH, iMap, REJRATE, 0.001);
        end
        
        % object segmentation
        fgMap = resSeg(npx, npy, MM, coorX, coorY, blkSiz, PO2, MRF_ENA, 2);
    elseif ii == nITER
        
        % the MV outliers removal cascade
        iMap = 1-fgMap;
        [MM, iniH] = gmeLs(coorX, coorY, npx, npy, iniH, iMap, REJRATE, 0.001);
        
        % object segmentation
        fgMap = resSeg(npx, npy, MM, coorX, coorY, blkSiz, PO2, MRF_ENA, 1);
        
        % segmentation test
        motionSegTest(npx, npy, MM, coorX, coorY, blkSiz, MRF_ENA, SEG_PATH, img_y, img_u, img_v, IDFRM);
        
    else
        % the MV outliers removal cascade
        iMap = 1-fgMap;
        [MM, iniH] = gmeLs(coorX, coorY, npx, npy, iniH, iMap, REJRATE, 0.001);
        
        % object segmentation
        fgMap = resSeg(npx, npy, MM, coorX, coorY, blkSiz, PO2, MRF_ENA, 2);
    end
end

% segmentation performed on residual MV field
function fgMap = resSeg(npx, npy, MM, coorX, coorY, blkSiz, PO2, MRF_ENA, mode)

% object segmentation
[mmpx, mmpy] = mvGen_f(MM, coorX, coorY);% SNR check-up
resMVFx = npx - mmpx;
resMVFy = npy - mmpy;

% moving region threshold
t = 1.2;
thX = t*(mean(resMVFx(:).^2)^0.5);
thY = t*(mean(resMVFy(:).^2)^0.5);

% MV quantization
[segMap] = initCluster(resMVFx, resMVFy, thX, thY);

% segmentation using Markov Random Field classificaiton
if MRF_ENA
    [segMap] = segMRF(resMVFx, resMVFy, segMap);
end

% check whether there is a foreground present in MV
% field
maxReg = max(segMap(:));
if maxReg > 1
    % if so, then, predict the foreground map
    fgMap = predFgMap(npx, npy, segMap, blkSiz, PO2, mode);
else
    fgMap = zeros(size(segMap));
end

function motionSegTest(npx, npy, MM, coorX, coorY, blkSiz, MRF_ENA, SEG_PATH, img_y, img_u, img_v, IDFRM)

% object segmentation
[mmpx, mmpy] = mvGen_f(MM, coorX, coorY);% SNR check-up
resMVFx = npx - mmpx;
resMVFy = npy - mmpy;

% moving region threshold
t = 1.2;
thX = t*(mean(resMVFx(:).^2)^0.5);
thY = t*(mean(resMVFy(:).^2)^0.5);

% MV quantization
[segMap] = initCluster(resMVFx, resMVFy, thX, thY);

DispInitMap(img_y, img_u, img_v, '', segMap, blkSiz);        
 