Skip to content

pvl_getISDdata.m and pvl_readtmy3.m update #1

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 3 commits into from
Jun 23, 2016
Merged

pvl_getISDdata.m and pvl_readtmy3.m update #1

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
c0c46f4
Update to pvl_getISDdata to use Haversine formula and include input p…
DanRiley Jun 14, 2016
4bd7bfe
Update to the default delta T calculation within the SPA code.
DanRiley Jun 23, 2016
3abd593
Update pvl_spa.m
DanRiley Jun 23, 2016
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
61 changes: 43 additions & 18 deletions pvl_getISDdata.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,19 +1,21 @@
function outfilen = pvl_getISDdata(lat,lon,yr,archive)
function outfilen = pvl_getISDdata(latitude,longitude,year,archive)
%
% getISDdata: fetch data for one year of integrated surface data for location
%
% Syntax
% data = getISDdata(latitude,longitude,year,archive)
%
% Inputs:
% latitude - latitude of location of interest
% longitude - longitude of location of interest
% latitude - latitude of location of interest in decimal degrees
% longitude - longitude of location of interest in decimal degrees
% year - desired year
% archive - string specifying path to local ISD archive. Current
% directory is used if archive is not specified.
%
% Output:
% outfilen - name for the returned data file
% outfilen - name for the returned data file in the form "USAF-WBAN-yyyy"
% where USAF is the USAF ID of the station and WBAN is the WBAN ID of
% the system, as provided by the ISD.
%
% Notes:
% This function copies up to two files from NOAA's ftp site to the
Expand All @@ -40,18 +42,30 @@
mod = vers(pos+6);
if versyr<2013 || ((versyr==2013) && mod=='a')
display(['Matlab ' vers ' detected'])
display('Function pvl_getISDdata requires R2013b or later')
return
error('Function pvl_getISDdata requires R2013b or later')
end

%% Get index data
% Change pwd to a local directory if you want to use an archive
if nargin < 4
%% Check the input data
p = inputParser;
p.addRequired('latitude',@(x) all(isscalar(x) & isnumeric(x) & x<=90 & x>=-90));
p.addRequired('longitude', @(x) all(isscalar(x) & isnumeric(x) & x<=180 & x>=-180));
p.addRequired('year',@(x) all(isscalar(x) & isnumeric(x) & (x >= 1901)));
p.addOptional('archive', @(x) ischar(x));
p.parse(latitude, longitude, year, varargin{:});

defaultchecker = {'archive'};
if ~any(strcmp(defaultchecker,p.UsingDefaults))
archive = pwd;
end

lat = p.Results.latitude;
lon = p.Results.longitude;
yr = p.Results.year;

%% Get index data
% Change pwd to a local directory if you want to use an archive

isdfile = 'isd-history.csv';
fname = [archive '\' isdfile];
fname = [archive filesep isdfile];
if ~exist(fname,'file')
% Get index file
ftphandle = ftp('ftp.ncdc.noaa.gov');
Expand All @@ -69,17 +83,23 @@
index.Properties.VariableNames{'ELEV_M_'} = 'ELEV';
index.USAF = char(index.USAF);
index.WBAN = char(index.WBAN);
index.LAT(strcmp('',index.LAT)) = {'99999'};
index.LON(strcmp('',index.LON)) = {'99999'};
index.ELEV(strcmp('',index.ELEV)) = {'99999'};
index.LAT(strcmp('',index.LAT)) = {'NaN'};
index.LON(strcmp('',index.LON)) = {'NaN'};
index.ELEV(strcmp('',index.ELEV)) = {'NaN'};
index.LAT = str2num(char(index.LAT)); %#ok<*ST2NM>
index.LON = str2num(char(index.LON));
index.ELEV = str2num(char(index.ELEV));
index.BEGIN = datenum(char(index.BEGIN),'yyyymmdd');
index.END = datenum(char(index.END),'yyyymmdd');

%% Find row for closest location (might iterate through several locations)
dist = hypot(index.LAT-lat,index.LON-lon);
%% Find row for closest location (might iterate through several locations) using haversine formula
% calculate the argument for the haversine formula, make sure the argument
% is between -1 and 1
arg = sqrt((sind((index.LAT - lat) ./ 2)).^2 + cosd(index.LAT).*cosd(lat).*(sind((index.LON-lon)./2)).^2);
arg(arg < -1) = -1;
arg(arg > 1) = 1;
% compute the distance in meters between the provided lat/lon and the each lat/lon from the isd.
dist = 2 .* 6371008.8 .* asin(arg); % mean earth radius is 6371.0088 km
[~,idx] = sort(dist);

% Find first row that contains year (might iterate through several years)
Expand All @@ -99,10 +119,15 @@

%
% get file from ftp server if not in archive
fname = [archive '\' outfilen];
fname = [archive filesep outfilen];
if ~exist(fname,'file')
addr = sprintf('ftp://ftp.ncdc.noaa.gov/pub/data/noaa/%4d/',yr);
gunzip([addr '/' outfilen '.gz'],archive);
try
gunzip([addr '/' outfilen '.gz'],archive);
catch
error(['Unable to reach the ftp site for the integrated surface database.' ...
' The data file cannot be accessed.']);
end
%% Alt approach if gunzip doesn't seem to be working
% % cd(ftphandle,num2str(yr));
% % gunzip(mget(ftphandle,[fname '.gz']));
Expand Down
125 changes: 123 additions & 2 deletions pvl_readtmy3.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,9 @@
% If a FileName is not provided, the user will be prompted to browse to
% an appropriate TMY3 file.
%
% TMY3 format was changed slightly on January 19, 2015 ([2]) to include
% Present weather variables along with changes to other variables.
%
% Input
% FileName - an optional argument which allows the user to select which
% TMY3 format file should be read. A file path may also be necessary if
Expand Down Expand Up @@ -106,12 +109,23 @@
% TMYData.Lprecipquantity - The period of accumulation for the liquid precipitation depth field, hour
% TMYData.LprecipSource - See [1], Table 1-5, 8760x1 cell array of strings
% TMYData.LprecipUncertainty - See [1], Table 1-6
% TMYData.PresWth - Present weather code, see [2]. Only available in
% TMY3 files created after the January 2015 TMY3 update.
% TMYData.PresWthSource - Present weather code source, see [2]. Only
% available in TMY3 files created after the Jan. 2015 TMY3 update.
% TMYData.PresWthUncertainty - Present weather code uncertainty, see
% [2]. Only available in TMY3 files after the Jan. 2015 update.
%
% Reference
% [1] Wilcox, S and Marion, W., 2008. Users Manual for TMY3 Data Sets,
% NREL/TP-581-43156, National Renewable Energy Laboratory. Available at
% <http://www.nrel.gov/docs/fy08osti/43156.pdf>.
%
% [2] National Solar Radiation Data Base, 1991-2005 Update: Typical
% Meteorological Year 3.
% http://rredc.nrel.gov/solar/old_data/nsrdb/1991-2005/tmy3/ retrieved June
% 6, 2016.
%
%
% See also
% DATEVEC PVL_MAKELOCATIONSTRUCT PVL_MAKETIMESTRUCT PVL_READTMY2
Expand All @@ -136,15 +150,118 @@
Header2 = textscan(FileID, '%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s',1,'Delimiter',',');
DataLines = textscan(FileID, '%s%s%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%f%s%f',8760,'Delimiter',',');
ST = fclose(FileID);
%% Read in the data from the struct created from the textscan

% Parse out the header information since it is common to both formats
TMYData.SiteID = Header1{1};
TMYData.StationName = Header1{2};
TMYData.StationState = Header1{3};
TMYData.SiteTimeZone = Header1{4};
TMYData.SiteLatitude = Header1{5};
TMYData.SiteLongitude = Header1{6};
TMYData.SiteElevation = Header1{7};
%% Determine if the file is from pre-January 19, 2015 or post-20150119
if isempty(strfind(char(DataLines{1,1}{1,1}), 'PresWth'))
% The string 'PresWth' was not found in DataLines{1,1}{1,1}, so we have a
% pre-2015 TMY3 file
TMYData = ParsePre2015TMY(DataLines);
else
% The string 'PresWth' was found in DataLines{1,1}{1,1}, so we have a
% post-2015 TMY3 file. We'll need to re-import the file to accomodate
% the additional columns.
FileID = fopen(FilePathAndNameAndExt);
Header1 = textscan(FileID, '%f%q%q%f%f%f%f', 1,'Delimiter',',');
Header2 = textscan(FileID, '%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s%s',1,'Delimiter',',');
DataLines = textscan(FileID, '%s%s%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%s%f%f%f%s%f%f%s%f',8760,'Delimiter',',');
ST = fclose(FileID);

TMYData = ParsePost2015TMY(DataLines);
end
end

function TMYData = ParsePre2015TMY(DataLines)

%% Read in the data from the struct created from the textscan for TMY3
% This fuction will read pre-2015 TMY3 files and sort them into column
% vectors. A different format is used for TMY3 files after the January 2015
% update, and thus requires a different function.

TMYData.DateString = DataLines{1};
TMYData.TimeString = DataLines{2};
TMYData.ETR = DataLines{3};
TMYData.ETRN = DataLines{4};
TMYData.GHI = DataLines{5};
TMYData.GHISource = DataLines{6};
TMYData.GHIUncertainty = DataLines{7};
TMYData.DNI = DataLines{8};
TMYData.DNISource = DataLines{9};
TMYData.DNIUncertainty = DataLines{10};
TMYData.DHI = DataLines{11};
TMYData.DHISource = DataLines{12};
TMYData.DHIUncertainty = DataLines{13};
TMYData.GHillum = DataLines{14};
TMYData.GHillumSource = DataLines{15};
TMYData.GHillumUncertainty = DataLines{16};
TMYData.DNillum = DataLines{17};
TMYData.DNillumSource = DataLines{18};
TMYData.DNillumUncertainty = DataLines{19};
TMYData.DHillum = DataLines{20};
TMYData.DHillumSource = DataLines{21};
TMYData.DHillumUncertainty = DataLines{22};
TMYData.Zenithlum = DataLines{23};
TMYData.ZenithlumSource = DataLines{24};
TMYData.ZenithlumUncertainty = DataLines{25};
TMYData.TotCld = DataLines{26};
TMYData.TotCldSource = DataLines{27};
TMYData.TotCldUnertainty = DataLines{28};
TMYData.OpqCld = DataLines{29};
TMYData.OpqCldSource = DataLines{30};
TMYData.OpqCldUncertainty = DataLines{31};
TMYData.DryBulb = DataLines{32};
TMYData.DryBulbSource = DataLines{33};
TMYData.DryBulbUncertainty = DataLines{34};
TMYData.DewPoint = DataLines{35};
TMYData.DewPointSource = DataLines{36};
TMYData.DewPointUncertainty = DataLines{37};
TMYData.RHum = DataLines{38};
TMYData.RHumSource = DataLines{39};
TMYData.RHumUncertainty = DataLines{40};
TMYData.Pressure = DataLines{41};
TMYData.PressureSource = DataLines{42};
TMYData.PressureUncertainty = DataLines{43};
TMYData.Wdir = DataLines{44};
TMYData.WdirSource = DataLines{45};
TMYData.WdirUncertainty = DataLines{46};
TMYData.Wspd = DataLines{47};
TMYData.WspdSource = DataLines{48};
TMYData.WspdUncertainty = DataLines{49};
TMYData.Hvis = DataLines{50};
TMYData.HvisSource = DataLines{51};
TMYData.HvisUncertainty = DataLines{52};
TMYData.CeilHgt = DataLines{53};
TMYData.CeilHgtSource = DataLines{54};
TMYData.CeilHgtUncertainty = DataLines{55};
TMYData.Pwat = DataLines{56};
TMYData.PwatSource = DataLines{57};
TMYData.PwatUncertainty = DataLines{58};
TMYData.AOD = DataLines{59};
TMYData.AODSource = DataLines{60};
TMYData.AODUncertainty = DataLines{61};
TMYData.Alb = DataLines{62};
TMYData.AlbSource = DataLines{63};
TMYData.AlbUncertainty = DataLines{64};
TMYData.Lprecipdepth = DataLines{65};
TMYData.Lprecipquantity = DataLines{66};
TMYData.LprecipSource = DataLines{67};
TMYData.LprecipUncertainty = DataLines{68};
%% Create a MATLAB datenum from the string date and time
TMYData.DateNumber = datenum(strcat(TMYData.DateString,'-',TMYData.TimeString),'mm/dd/yyyy-HH:MM'); %MATLAB datenum format
end



function TMYData = ParsePost2015TMY(DataLines)

%% Read in the data from the struct created from the textscan
TMYData.DateString = DataLines{1};
TMYData.TimeString = DataLines{2};
TMYData.ETR = DataLines{3};
Expand Down Expand Up @@ -213,5 +330,9 @@
TMYData.Lprecipquantity = DataLines{66};
TMYData.LprecipSource = DataLines{67};
TMYData.LprecipUncertainty = DataLines{68};
TMYData.PresWth = DataLines{69};
TMYData.PresWthSource = DataLines{70};
TMYData.PresWthUncertainty = DataLines{71};
%% Create a MATLAB datenum from the string date and time
TMYData.DateNumber = datenum(strcat(TMYData.DateString,'-',TMYData.TimeString),'mm/dd/yyyy-HH:MM'); %MATLAB datenum format
TMYData.DateNumber = datenum(strcat(TMYData.DateString,'-',TMYData.TimeString),'mm/dd/yyyy-HH:MM'); %MATLAB datenum format
end
Loading