Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

sum hourly precipiation data into individual storm events

Asked by Oli on 26 Apr 2012

Hi all. I have hourly precip data that I am trying to sum into separate storm events. 24 hours of no precipitation (or 24 zero values) signify the end of the storm event. I'm trying to create a code that sums the precip values until the 24 hours of no rainfall and then starts summing the following storm event. The output will be an array with individual storm event depths of the time series. Does anyone have any suggestions the best method to do this?

0 Comments

Oli

Products

2 Answers

Answer by Image Analyst on 26 Apr 2012

Well if you have the Image Processing Toolbox you can get it in 4 lines: a line to identify rain-free hours. Another line to get rid of small stretches of no rain and combine storms on either side into a single storm. The third line to identify stretches of hourly measurements as individually numbered storms, and the fourth line to actually do the measurements of each numbered storm. Here's the code:

rain = rand(1, 50)
zeroIndices = rain<0.5;
rain(zeroIndices) = 0
% Now we have some sample data
% Let's start the analysis:
%========================================================
% Key part right here:
% Find out where it's zero (no rain):
binaryData = rain == 0
% Get rid of small stretches 2 or less in length:
% Change 3 to 24 if you want to combine storms 23 hours or closer together.
rainFreeHours = bwareaopen(binaryData, 3)
% Now rain-free = 1 and raining = 0.  
% Invert it to find rainy stretches, then label it to find individual storms.
[labeledStorms numStorms] = bwlabel(~rainFreeHours)
% Now call regionprops to get the amount of rain over 
% all hours of each storm:
measurements = regionprops(labeledStorms, rain, 'PixelValues');
%========================================================
% We're done!
%Now report results
for storm = 1 : numStorms
	rainInThisStorm(storm) = sum(measurements(storm).PixelValues);
	fprintf('Rainfall total in storm #%d = %.4f\n', ...
		storm, rainInThisStorm(storm));
end
% For fun, report the rain over ALL the hours.
totalRainOverAllHours1 = sum(rain)
totalRainOverAllHours2 = sum(rainInThisStorm) % Will be the same.

2 Comments

Image Analyst on 26 Apr 2012

Roger just posted on to your duplicate post in the newsgroup: http://groups.google.com/group/comp.soft-sys.matlab/browse_frm/thread/f34210f9aaf3da6e/da9dfdef7e004370?hl=en#da9dfdef7e004370

Oli on 27 Apr 2012

I don't have the Image processing toolbox so I need a code for basic matlab. But thank you so much for helping and for pointing out the newsgroup reply! It is my first time using mathworks so I'm not sure how I duplicated this post. The newsgroup suggestion worked for me.

Image Analyst
Answer by Geoff on 26 Apr 2012

You can use basic MatLab stuff for this too.

Convolution will detect stretches of no rainfall.

dry = conv(rainfall, ones(1,24), 'valid') == 0;

You can then detect rain->dry events with diff (dry changes from 0 to 1)

stormEnd = find(diff(dry) == 1) + 1;

That will give you the index of the first dry hour after each storm.

Now you just treat your data as a series of events. Chuck the first and last index of rainfall in there.

evt = [1, stormEnd, numel(rainfall)];

And sum the rainfall between each event:

totals = arrayfun( @(n) sum(rainfall(evt(n-1):evt(n))), 2:numel(evt) );

2 Comments

Oli on 27 Apr 2012

I attempted this method but got an error at evt = [1, stormEnd, numel(rainfall)]; Error using ==> horzcat
CAT arguments dimensions are not consistent. I am not familiar with some of the operations above (somewhat new to matlab too), so I'm not sure what the problem it. I ended up using the method from the newsgroup reply. Thank you for helping though!

Geoff on 29 Apr 2012

Oh, your data was probably in columns and my code assumed rows. You should transpose 'stormEnd'. Never mind =)

Geoff

Contact us