function hedm = auto_read_hed_day_WHOTS_2019(path_working_ddd,ddd,yyyy)
%{
Reads hourly 2 Hz heading files for one day, averages to 10-min
Loops through 24 hours, calling auto_read_hed_roll_WHOTS_2018 for each hour.

Inputs: path_working_ddd: path to daily folder containing hourly files
        ddd: day-of-year variable
        yyyy: year string

Output: hedm: 1440x3 array of heading and roll angle
        Returns NaN for data gaps

    1   jd_ref      10 min timestamp
    2   hed         heading
    3   roll        roll angle
%}

fclose all;

hedm = NaN(144,5);
delta = double(10/1440);
last = ddd + 144*delta;
jd_bin = ddd:delta:last; 
jd_bin = jd_bin'; % 10 min bin edges
hedm(:,1) = jd_bin(1:end-1);
jd = sprintf('%03i',ddd);

for hhh = 0:23              % cycle thru 24 hourly gprm files
    hr = sprintf('%02i',hhh);
    dfl = fullfile(path_working_ddd,['hed0' yyyy(3:4),jd,hr,'_raw.txt']);
    if exist(dfl,'file')==2
        hed = auto_read_hed_pitch_WHOTS_2019(dfl,ddd,hhh);
    else
        continue % no data this hour
    end
    
    [hed_n, hed_e] = pol2cart(hed(:,2)*pi/180,1);  % N/E heading components
    hed = [hed,hed_n,hed_e]; % add N/E components as cols 4,5
    
    % average this hour's data into 10min hedm array
    start = ddd + hhh/24.0; % jd start time this hour
    diff = jd_bin - start;  % look for closest time stamp to start
    [~,ii] = min(abs(diff));
    temp = interval_avg(hed(:,1), hed(:,2:5), jd_bin(ii:ii+6));
    hedm(ii:ii+5,2:5) = temp(:,2:5);
end

% recompute cog/sog from averaged E/N heading components
hed_10 = atan2(hedm(:,5), hedm(:,4))*180/pi; % from avg components
hed_10 = mod(hed_10 + 360*3, 360);
hedm(:,2) = hed_10;
hedm = hedm(:,1:3);

end



