function [time,swdn,lwdn,swup,lwup,precip,press_pa] = read_rad(path, period)
warning off
z = dir([path 'nsaqcrad1longC1*']);
Archive = [];
for i = 1:length(z)
    x = rd_netcdf([path z(i).name]);
    y = attribute_list([path z(i).name]);
    time = double(x.base_time)/60/60/24 + datenum([1970 1 1 0 0 0]) + double(x.time/60/60/24);
    % QC
    swdn = x.BestEstimate_down_short_hemisp;
    swdn(x.qc_BestEstimate_down_short_hemisp > 0) = NaN;
    lwdn = x.down_long_hemisp;
    lwdn(x.qc_down_long_hemisp > 0) = NaN;
    swup = x.up_short_hemisp;
    swup(x.qc_up_short_hemisp > 0) = NaN;
    lwup = x.up_long_hemisp;
    lwup(x.qc_up_long_hemisp > 0) = NaN;
    precip = x.precip;
    precip(x.qc_precip > 0) = NaN;
    press = x.press;
    press(x.qc_press > 0) = NaN;
    % accumulate
    Archive = [Archive; time double(swdn) double(lwdn) double(swup) double(lwup) double(precip) double(press)];
end

Archive(Archive(:,1) < datenum(period(1,:)) | Archive(:,1) > datenum(period(2,:))+ 1 - 1/24/60,:) = [];
clear mex

% interpret
time = Archive(:,1);
swdn = Archive(:,2);
lwdn = Archive(:,3);
swup = Archive(:,4);
lwup = Archive(:,5);
precip = Archive(:,6);
press = Archive(:,7);

% convert
press_pa = press*1000;
