function par = vfd_clippar_loop_par1(test_structs, clip_structs, feature_base_dir)
% VFD_CLIPPAR_LOOP_PAR1
%  Computes the vfd_par1 parameter of the VQM_VFD model. The vfd_par1
%  parameter quantifies the perceptual motion distortion of a video stream
%  with variable frame delay (VFD).  This parameter is necessary to capture
%  information that is lost by the VFD feature-based clip parameters, which
%  use an original sequence that is changed to match the processed
%  sequence. The vfd_par1 parameter is described in detail in NTIA
%  TM-11-475. The function takes variable (1) 'test_structs' (of the same
%  format as GTests), which describes each video test and the location of
%  the associated video, and variable (2) 'clip_structs' (of the same
%  format as GClips), which specifies the set of video clips, and (3) the
%  feature directory 'feature_base_dir', which specifies the directory
%  location of the VFD results files (named vfd_test.mat, where test is the
%  name of the video clip test). The returned par will be an ITS parameter
%  structure that contains the parameters for all the processed clips in
%  clip_structs. The sort_clips_by function is used to sort the clips
%  before parameter calculation, and the original clips are skipped if they
%  are encountered.
% SYNTAX
%  par = vfd_clippar_loop_par1(test_structs, clip_structs, feature_base_dir)
% DESCRIPTION
%  See NTIA TM-11-475 for a technical description of the vfd_par1
%  parameter. This function loads the VFD results files (in the same format
%  as generated by the script est_var_frame_delays_script.m in the
%  ITS_VIDEO library), and uses the information in those variable names to
%  extract the vfd_par1 parameter for the user-specified set of video clips
%  in 'clip_structs'.
%

% State variable that gives the last test whose vfd information was loaded.
% This speeds up the reading of the VFD results files as they are only read
% when the test changes.
last_test = '';

% Define the par_name for the returned parameter structure.
par.par_name = {'vfd_par1'};

% Loop to calculate the VFD clip parameter for each clip in this clipset,
% after sorting alphabetically.
offsets = sort_clips_by('none', clip_structs, test_structs);
clip_structs = clip_structs(offsets);
ccnt = 1;  % processed clip counter for returned par array
for loop = 1:length(offsets)
    
    % Skip if an original clip was included.
    if strcmpi(clip_structs(loop).hrc{1},'original'),
        continue;
    end
    
    %  Pick off this clip's test, scene, and HRC name
    clip_test =  clip_structs(loop).test{1};
    clip_scene = clip_structs(loop).scene{1};
    clip_hrc = clip_structs(loop).hrc{1};
    
    %  Load the appropriate VFD delay information for this test.
    %  This code utilizes the variables 'this_clipset' and
    %  'results' from the VFD mat files.  The align_start and
    %  align_stop of the clip will be checked against the VFD info to
    %  make sure the VFD alignment is synched to the clipset's alignment.
    if(~strcmpi(last_test, clip_test))  % Only load the VFD results file if it is different from last time
        last_test = clip_test;
        load([feature_base_dir  'vfd_' clip_test '.mat']);
    end
    vfd_ind = find_clip(this_clipset, clip_test, clip_scene, clip_hrc);
    
    % Check the align_start and align_stop of the video clip as specified
    % by clip_struct and the computed VFD results file (they must match).
    if (clip_structs(loop).align_start ~= this_clipset(vfd_ind).align_start ...
            || clip_structs(loop).align_stop ~= this_clipset(vfd_ind).align_stop)
        fprintf('\tError encountered in vfd_clippar_loop_par1 for processed clip %s:%s(%s)\n', clip_test, clip_scene, clip_hrc);
        error('Invalid VFD information for processed clip: align_start or align_stop')
    end
    
    %  Check the original's alignment with the VFD alignment
    vfd_ind_orig = find_clip(this_clipset, clip_test, clip_scene, 'original');  % original location in VFD variable 'this_clipset'
    ind_orig = find_clip(clip_structs, clip_test, clip_scene, 'original');  % original location in clip_structs
    if (clip_structs(ind_orig).align_start ~= this_clipset(vfd_ind_orig).align_start ...
            || clip_structs(ind_orig).align_stop ~= this_clipset(vfd_ind_orig).align_stop)
        fprintf('\tError encountered in vfd_clippar_loop_par1 for the original clip of %s:%s(%s)\n', clip_test, clip_scene, clip_hrc);
        error('Invalid VFD information for original clip: align_start or align_stop')
    end
    
    % Set the initial alignment point (in frames) for the VFD
    % correction.  first_align is considered frame or field
    % number 1 in the VFD results file.
    first_align = this_clipset(vfd_ind_orig).align_start - this_clipset(vfd_ind_orig).loc_start + 1;
    
    % Set the interlaced flag for the VFD correction.
    if (strcmpi(this_clipset(vfd_ind_orig).video_standard, 'interlace_lower_field_first') == 1)
        is_interlaced = 1;
        field_first = 1;
        first_align = 2*first_align-1;  % convert to fields
    elseif (strcmpi(this_clipset(vfd_ind_orig).video_standard, 'interlace_upper_field_first') == 1)
        is_interlaced = 1;
        field_first = 2;
        first_align = 2*first_align-1;  % convert to fields
    else
        is_interlaced = 0;  % progressive
    end
    
    % Find the VFD original alignment results for this processed
    % clip.  You must eliminate the originals from this_clipset to
    % find the proper results index.
    this_clipset_proc = this_clipset(find_clip(this_clipset, '*', '*', 'original', 'not'));
    results_ind = find_clip(this_clipset_proc, clip_test, clip_scene, clip_hrc);
    this_result = results{results_ind};  % Pick off the alignment results for this clip from the retrieved VFD info
    this_result_fuzzy = results_fuzzy{results_ind};
    
    % Determine the total number of aligned frames that are available
    % (in frames) for the processed clip.  This will be used to
    % generate a gclips alignment result if the VFD algorithm failed.
    nframes = clip_structs(loop).align_stop-clip_structs(loop).align_start+1;
    % Reduce by one if the processed clip is reframed
    if (is_interlaced && mod(clip_structs(loop).spatial.vertical,2))
        nframes = nframes-1;
    end
    
    %  If the VFD alignment algorithm failed, then this_result==0.
    %  In that case, use first_align to generate this_result so there will
    %  be no difference between VFD results and normal results.
    if (this_result == 0)
        
        if (is_interlaced)
            npts = nframes*2;
        else
            npts = nframes;
        end
        
        this_result = first_align:first_align+npts-1;
        
        % Convert the VFD information to a vector that gives abnormal
        % forward field/frame jumps.  See description of algorithm
        % given below.
        this_result = max([0 abs(diff(this_result))-1], 0);
        
    else
        
        npts = length(this_result);
        
        % Convert the VFD information to a vector that gives abnormal
        % forward field/frame jumps.  Subtracting 1 from the abs(diff)
        % and maxing with zero produces a parameter that (1) does not
        % penalize for normal field/frame delivery (where the VFD
        % field/frame indices increase by one from one field/frame to
        % the next), (2) does not penalize for progressive frame
        % repeats (where the VFD frame indices stay fixed from one
        % frame to the next), and (3) does not penalize for interlaced
        % frame repeats (where the VFD field indices jump back in time
        % from one field to the next). A non-impairment value is used
        % for the first field/frame (which must be padded since it's
        % diff is not available).  If one assumes that both the
        % early and late time sides used for the frame jump estimates
        % are absolutely correct (i.e., no fuzzy alignments), then this
        % equation would result:
        %
        % this_result = max([0 abs(diff(this_result))-1], 0);
        %
        % One could also use a modification of the above diff that
        % subtracts the current frame position from the max fuzzy
        % causal alignment of the prior frame.  This is thus a more
        % conservative estimate of the frame jumps (i.e., smaller frame
        % jumps) at any given position since it allows uncertainty on
        % the early time side. However, this allows no uncertainty on
        % the late time side, since that alignment is considered
        % absolutely correct. This equation would result:
        %
        % fuzzy_max = min(max(this_result_fuzzy(:,1:npts-1)),this_result(2:npts));
        % this_result = max([0 abs(this_result(2:npts)-fuzzy_max)-1], 0);
        %
        % A final modification assumes uncertainty on both the early and
        % late time sides so it is the most forgiving of all three.
        % Thus, frame jumps are only penalized when they are absolutely
        % certain to be correct (i.e., no fuzzy overlapping).  Even
        % then, when analyzing this algorithm for the test scene
        % t1y1_kiel12_A1pass, it is clear that the VFD estimation
        % algorithm is overconfident for some alignments (e.g., not
        % enough fuzzy alignments are included for field 511). The
        % final equation is:
        fuzzy_max_early = min(max(this_result_fuzzy(:,1:npts-1)),this_result(2:npts));
        fuzzy_min_late = max(min(this_result_fuzzy(:,2:npts)),fuzzy_max_early);
        this_result = max([0 abs(fuzzy_min_late-fuzzy_max_early)-1], 0);
        
    end
    
    % Calculate vfd_par1 for this clip using log10 of the RMS error.
    this_par = log10(sqrt(mean(this_result.^2))+1);
    
    % Fill the clip_structs information into the returned parameter.
    par.clip_name{ccnt} = sprintf('%s_%s_%s', clip_structs(loop).test{1}, ...
        clip_structs(loop).scene{1}, clip_structs(loop).hrc{1});
    par.inlsa_mos(ccnt) = clip_structs(loop).inlsa_mos;
    par.mos(ccnt) = clip_structs(loop).mos;
    par.data(1, ccnt) = this_par;
    
    % Update clip counter.
    ccnt = ccnt + 1;
    
end







