Unsure on how to achieve goal - Function output isn't right
1 view (last 30 days)
Show older comments
Hello,
I'm trying to work out what I'm doing wrong with my approach to a question. I am trying to create a function to calculate solar radiation absorption as a function of time of day&year, cloud cover, and latitude.
My code is as follows:
function [abs_sol_rad]=solarrad(nc,lat,t,I)
% Function to calculate absorbed solar radiation given cloud cover,
% latitude, seconds from midnight and day in the year.
% Inputs:
% I= day number (365 days)
%lat = latitude
%t = seconds from mightnight
%nc = cloud cover from 0 (cloudless) to 100 (full cloud cover)
% Outputs are abs_sol_rad=absorption solar radiation with units W/m^2
d0=23.45; %angle between equatorial plane and earth rotational axis
wd=(2*pi)/86400;wy=(2*pi)/365; %wd= day rotational angle frequency, wy=year rotation angle frequency
I=365;
hr=24;
abs_sol_rad=nan(hr,I);
for I=1:365 %
for hr=1:24
t=hr*3600;
end
delta=d0*(pi/180)*(cos(wy*(I+9))); %declination of the sun from winter solstice (day 356)
S=1357+45*(cos(wy*I)); %solar constant
sin_h=max(0,(sind(lat)*sind(delta))+(cosd(lat)*cos(delta)*cos(wd*(t+43200)))); %sun height (angle) equation
m=asin(sin_h); %optical length
Cext=0.128-(0.0235.*log(m)); %integrated extinction coefficient
for h=m
i=((pi/2)-h); %angle of incidence
j=asin(0.75*(sin(i))); %angle of reflection
r=0.5*(abs((sin(i-j).^2)/(sin(i+j).^2)+(tan(i-j).^2/tan(i+j).^2))); %reflectance
end
if h>0
Insd=S*(exp(-Cext*m))*(sin_h*(1-0.71*nc));
end
if h<=0
Insd=0;
end
Insg=0.52*nc*Insd; %Global radiation
Qsun=(Insd*(1-r))+(0.97*Insg); %absorbed radiation in water
abs_sol_rad(hr,I)=Qsun;
end
end
The output I'm getting is just for the last time step I'm wanting (want a time step of 1 hour each day, 365 days). I need the output to be a matrix 24x365. Which I believe shows absorption when given latitude and cloud cover I guess.
Do I need to create a secondary script to use this function to create an array?
Should I just create a meshgrid to produce an 24x365 array and then just apply the function to each cell?
2 Comments
Answers (1)
Miriam
on 17 Nov 2018
Have the second for loop end after assigning the value of abs_sol_rad, ie:
function [abs_sol_rad]=solarrad(nc,lat,t,I)
%% ...
for I=1:365 %
for hr=1:24
t=hr*3600;
delta=d0*(pi/180)*(cos(wy*(I+9)));
%% ...
abs_sol_rad(hr,I)=Qsun;
end
end
end
2 Comments
Miriam
on 17 Nov 2018
Are you still using hr=1:24? I'm assuming you've adjusted your code slightly since your indexing has reversed. If you also, say, adjusted hr=0:23 that would explain why you need to add one, as MATLAB indexing begins from 1.
See Also
Categories
Find more on Environment and Clutter in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!