function [hed, bad_hed] = read_hed_pitch_HiWinGS_2013(dfl,ddd,hhh)
%{
reads heading/pitch file from PSD Hemisphere system

outputs: [hed, bad_hed]

    3 columns in array hed:
        1   jd_ref      10 Hz timestamp
        2   heding      deg
        3   pitch       deg

    bad_hed = count of bad lines in file

input parameters: dfl = string path to hed0 file
                  ddd = julian date
                  hhh = hour
%}

%% initialize
disp(['Reading hed file for hour ',int2str(hhh)]);
flist = fopen(dfl);

temp = {};    % empty arrays for textscan
stx = [];

%% read file
% 0000009 $PSAT,HPR,025959.80,319.48,,-0.47,N*1F
% 0000118 $PSAT,HPR,025959.90,319.49,,-0.46,N*1E
while feof(flist)==0,	% read entire file into cell array
    try
        % temp will be 1x6 cell array, cells will be ~36000x1 arrays
        temp = textscan(flist,'%2f%2f%3f %*5c %*3c %f %f %*f %f %*[^\n]','delimiter',', ','headerlines',1,'emptyvalue', NaN);
        stx = [stx cell2mat(temp)'];  % 6x~36000 array
    catch
        for ii=1:6      % length of last cell reflects missing values from any/all fields
             if length(temp{1,ii})~=length(temp{1,6})
                 temp{1,ii}(length(temp{1,6})+1)=[];  %truncate length of all cells to length(temp{1,11})
             end;
        end;
        if ~isempty(temp{1,1}); % save good values from last read, if any
            stx = [stx cell2mat(temp)'];
        end;
    end;
end;
fclose(flist);
jd_hed = ddd + (hhh + (stx(1,:)+(stx(2,:)+stx(3,:)/1000)/60)/60)/24;
heading = stx(5,:);
pitch = stx(6,:);

% count bad or missing points in heading file
missing = 36000 - length(stx(5,:));
if missing<0; missing=0; end;
bad_hed = length(find(isnan(heading))) + missing;

% fix NaNs
heading = replace_NaN_nearest_neighbor(heading);
pitch = replace_NaN_nearest_neighbor(pitch);

%% Filter noise from heading & pitch
[bhi,ahi,~,~,~,~] = Filtcoef(10,1/30);
heading = heading - filtfilt(bhi,ahi,heading);
pitch = pitch - filtfilt(bhi,ahi,pitch);

%% format to 3x36000 output array
start = (ddd+hhh/24);
delta = double(1.0/864000);
last = start + 35999*delta;
tref = [start:delta:last];
npts = length(jd_hed);
hed = zeros(36000,3);
hed(:,1) = tref;
if npts>36000
    hed(:,2) = heading(1:36000)';
    hed(:,3) = pitch(1:36000)';
elseif npts<36000
    dn = 36000 - npts;
    hed(:,2) = [heading';ones(dn,1).*heading(npts)]; % dup last few points to fill
    hed(:,3) = [pitch';ones(dn,1).*pitch(npts)];
else
    hed(:,2) = heading';
    hed(:,3) = pitch';
end;

