Creating a matrix after giving initial conditions

22 ビュー (過去 30 日間)
Damith
Damith 2015 年 3 月 12 日
コメント済み: dpb 2015 年 3 月 31 日
Hi,
I need to create a matrix when given the begin year, month, day, time and end year, month, day and time.
for example; (I do not need to see the header in the matrix)
I would like to enter:
begin year = 1998
begin month = 1
begin day = 1
begin time = 0
end year = 1998
end month = 2
end day = 15
end time = 12
(as shown in the table below)
But, if I enter leap years (e.g. 2000,2004, etc) it should include Feb 29 and do the same time steps (0 to 21) as shown in the table.
Then, the program should generate the output for the time periods given above. See the sample below.
Year Month Day Time
1998 1 1 0
1998 1 1 3
1998 1 1 6
1998 1 1 9
1998 1 1 12
1998 1 1 15
1998 1 1 18
1998 1 1 21
1998 1 2 0
.
.
.
1998 2 15 12
Any help is appreciated.
Thanks in advance.

採用された回答

dpb
dpb 2015 年 3 月 12 日
編集済み: dpb 2015 年 3 月 12 日
yr=2000;
dn=datevec(datenum(yr,1,1,[0:3:24*(365+isleapyr(yr))-1].',0,0))
>> [dn(1:5,:);dn(470:485,:);dn(length(dn)-10,:)]
ans =
2000 1 1 0 0 0
2000 1 1 3 0 0
2000 1 1 6 0 0
2000 1 1 9 0 0
2000 1 1 12 0 0
...
2000 2 28 15 0 0
2000 2 28 18 0 0
2000 2 28 21 0 0
2000 2 29 0 0 0
2000 2 29 3 0 0
2000 2 29 6 0 0
2000 2 29 9 0 0
2000 2 29 12 0 0
2000 2 29 15 0 0
2000 2 29 18 0 0
2000 2 29 21 0 0
2000 3 1 0 0 0
2000 3 1 3 0 0
2000 3 1 6 0 0
2000 3 1 9 0 0
2000 3 1 12 0 0
...
2000 12 30 15 0 0
2000 12 30 18 0 0
2000 12 30 21 0 0
2000 12 31 0 0 0
2000 12 31 3 0 0
2000 12 31 6 0 0
2000 12 31 9 0 0
2000 12 31 12 0 0
2000 12 31 15 0 0
2000 12 31 18 0 0
2000 12 31 21 0 0
>>
I introduced the ellipses in the output to show the breaks clearly...
Salt to suit...
Oh, isleapyr is one of my utilities...
function is=isleapyr(yr)
% returns T for input year being a leapyear
is=(datenum(yr+1,1,1)-datenum(yr,1,1))==366;
there may be a TMW-supplied equivalent since R2012b which lacks it.
  11 件のコメント
Damith
Damith 2015 年 3 月 18 日
yes, ignore the days that are incomplete.
Will see whether your plan works or not :)

サインインしてコメントする。

その他の回答 (2 件)

dpb
dpb 2015 年 3 月 18 日
編集済み: dpb 2015 年 3 月 20 日
Since dir returns a sorted list, I think I'd just process it similarly as you are but simply test the days within the loop by parsing the date/time string within a double loop. One could, as mentioned, get fancy and do some serious attempts at matching datenums but just doesn't seem to be needed here...
I'd start with something like the following which builds an array of the summed third column value in matching days by year. Each day is in a given column, only a small additional logic is needed to keep track of which day is present, a month and day row filled in by column while in the loop...
d=dir('TRMM_*_newntcl.csv');
ymdh=cell2mat(textscan([d.name], ...
'TRMM_%4d_%2d_%2d%2d_newntcl.csv', ...
'collectoutput',true));
idx=find(ymdh(:,end)==0); % get the beginning day 0 hr locations
N=length(idx);
idx=[idx;length(ymdh)+1];
precip=nan(51681,fix(N)); % preallocate on num whole days
k=0;
for i=1:length(idx)-1
isOK=~any(diff(ymdh(idx(i):idx(i+1)-1,3))) & idx(i+1)-idx(i)==8;
if ~isOK, continue, end % either aren't eight or it's a wrong day mixed in
data=zeros(51681,8); % preallocate based on num lat/lon points
for j=0:7 % ok, process each of these, they're one day's worth
data(:,j+1)=csvread(d(idx(i)+j).name,0,2); % add hourly data
end
k=k+1;
precip(:,k)=sum(data,2); % the daily sums are in an array
fn=d(idx(i)).name; fn(16:17)=[]; % fn(1:1)='O';
disp(fn)
csvwrite(fn,precip(:,k))
end
precip(:,k+1:end)=[]; % remove columns to leave complete days only
This snippet reads only the third column on the assumption you can easily add the constant position columns as you see fit.
I've never used a netCDF file so on casting it into that form I'm not of much use but the above should provide the basic building block to get the data consolidated by day. Note I specifically avoided cell arrays here and did the summing in place while accumulating the data to keep total memory down to minimal amount...as you've already done, use this inside the higher level loop to traverse the directories.
Unless there aren't multiple partial days, the final value of k will be the same as the floor(L/8) calculation. If there are there will be that many zero columns that can simply remove.
Hope that helps...
  30 件のコメント
dpb
dpb 2015 年 3 月 24 日
Go back an reread the comments starting on Mar 20@17:45 where I described how to build the 3D array of daily totals that you can subsequently simply sum by column for each plane and by plane for monthly and yearly totals for each year. What's the issue in doing that, it's simply modifying the size of the precip array initially and then using the ymdh array info on day and month to set the indices instead of just incrementing a counter k as doing now.

サインインしてコメントする。


dpb
dpb 2015 年 3 月 24 日
編集済み: dpb 2015 年 3 月 29 日
d=dir('TRMM_*_newntcl.csv');
ymdh=cell2mat(textscan([d.name], ...
'TRMM_%4d_%2d_%2d%2d_newntcl.csv', ...
'collectoutput',true));
idx=find(ymdh(:,end)==0); % get the beginning day 0 hr locations
N=length(idx);
idx=[idx;length(ymdh)+1];
precip=nan(51681,31,12); % preallocate on num whole days
%fprintf('Building daily sums for %d: \n',ymdh(1,1))
dispstat('','init')
for i=1:length(idx)-1
isOK=~any(diff(ymdh(idx(i):idx(i+1)-1,3))) & idx(i+1)-idx(i)==8;
if ~isOK, continue, end % either aren't eight or it's a wrong day mixed in
%msg=sprintf(' Month: %2d Day: %2d',ymdh(idx(i),2:3));
%dispstat(msg,'keepthis')
data=zeros(51681,8); % preallocate based on num lat/lon points
for j=0:7 % ok, process each of these, they're one day's worth
%msg=sprintf('Reading hourly file: %s: ',d(idx(i)+j).name);
%dispstat(msg)
data(:,j+1)=csvread(d(idx(i)+j).name,0,2); % add hourly data
end
mo=ymdh(idx(i),2); % this month
da=ymdh(idx(i),3); % this day
precip(:,da,mo)=sum(data,2);
fn=d(idx(i)).name; fn(16:17)=[]; fn(1:1)='O';
%disp(fn)
csvwrite(fn,precip(:,k))
end
mTot=squeeze(sum(precip,2)); % monthly totals in each column 1-12
yTot=sum(mTot,2); % yearly total
Now just write mTot and yTot to files as you want them. As noted before, you'll probably want to do some checking also on whether there's data in every day of the month before accumulating that total and then for the months before the yearly total or you may have totals that aren't actually representative of the entire month/year.
Using nan to initialize instead of zero is one way, then any sums with NaN as a value are indicative altho there you've got to fix up the detail of how many days are in each month but that's pretty simple bookkeeping...
  30 件のコメント
dpb
dpb 2015 年 3 月 31 日
BTW, if you haven't yet figured this out, you could save tons of time if you were to run the above script once for each year subdirectory and then
save precip
in that working directory. Then you can simply
load precip
for that year and experiment with how to create and save the monthly and yearly sums to your heart's content and build whatever files from there without having to reread the same values over and over and over which is quite slow.

サインインしてコメントする。

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by