Sonic Wind Correction Program
This program is called by
prof_file_all2.m
and uses the variable names as assigned there. Note that
two functions are also included "uv.m" and "spddir.m". These
are used to convert the wind speed and direction into
u, v components and vice-versa.
<![CDATA[
%windcorrect.m
%Peter Guest 4/15/95
%read in data
%set all bad value flags to NaN
%calculate relative wind direction
rwd1=wd1-twr_orien;
hi=find(rwd1>360);
lo=find(rwd1<0);
rwd1(hi)=rwd1(hi)-360;
rwd1(lo)=rwd1(lo)+360;
rwd2=wd2-twr_orien;
hi=find(rwd2>360);
lo=find(rwd2<0);
rwd2(hi)=rwd2(hi)-360;
rwd2(lo)=rwd2(lo)+360;
rwd3=wd3-twr_orien;
hi=find(rwd3>360);
lo=find(rwd3<0);
rwd3(hi)=rwd3(hi)-360;
rwd3(lo)=rwd3(lo)+360;
rwd4=wd4-twr_orien;
hi=find(rwd4>360);
lo=find(rwd4<0);
rwd4(hi)=rwd4(hi)-360;
rwd4(lo)=rwd4(lo)+360;
rwd5=wd5-twr_orien;
hi=find(rwd5>360);
lo=find(rwd5<0);
rwd5(hi)=rwd5(hi)-360;
rwd5(lo)=rwd5(lo)+360;
%Defined U,V such that
%wind from relative 0 means U>0 V=0
%wind from relative 90 means V>0 U=0
[u1 v1] = uv(ws1,rwd1);
[u2 v2] = uv(ws2,rwd2);
[u3 v3] = uv(ws3,rwd3);
[u4 v4] = uv(ws4,rwd4);
[u5 v5] = uv(ws5,rwd5);
%correct biases and slopes assume mean bias and slope is "true wind"
bcoru1= .224;%mean bias
scoru1= -0.0365; %mean slope
bcorv1= -.078;%mean bias
scorv1= -0.0372; %mean slope
ib400=find(jd<400);
ia400=find(jd>=400);
%correct u,v components
%cb represents bias correction only
%c represents bias and slope correction
u1cb=u1-bcoru1;
u1c=u1cb+u1cb.*(0-scoru1);
u2cb=u2+0.35-bcoru1;
u2c=u2cb+u2cb.*(0-scoru1);
u3cb(ib400,1)=u3(ib400,1)+0.20-bcoru1;
u3c(ib400,1)=u3cb(ib400,1)+u3cb(ib400,1).*(0-scoru1);
u3cb(ia400,1)=u3(ia400,1)+0.40-bcoru1;
u3c(ia400,1)=u3cb(ia400,1)+u3cb(ia400,1).*(-0.098-scoru1);
u4cb=u4+0.29-bcoru1;
u4c=u4cb+u4cb.*(-0.063-scoru1);
u5cb=u5+0.16-bcoru1;
u5c=u5cb+u5cb.*(-0.050-scoru1);
v1cb=v1+0-bcorv1;
v1c=v1cb+v1cb.*(0-scorv1);
v2cb=v2+0.-bcorv1;
v2c=v2cb+v2cb.*(0-scorv1);
v3cb(ib400,1)=v3(ib400,1)+0.-bcorv1;
v3c(ib400,1)=v3cb(ib400,1)+v3cb(ib400,1).*(-.023-scorv1);
v3cb(ia400,1)=v3(ia400,1)+0.-bcorv1;
v3c(ia400,1)=v3cb(ia400,1)+v3cb(ia400,1).*(-.0665-scorv1);
v4cb=v4-0.11-bcorv1;
v4c=v4cb+v4cb.*(-.050-scorv1);
v5cb=v5-0.16-bcorv1;
v5c=v5cb+v5cb.*(-.034-scorv1);
%convert back to wind speed and direction bias correction only
[wscb1 rwdcb1]=spddir(u1cb,v1cb);
[wscb2 rwdcb2]=spddir(u2cb,v2cb);
[wscb3 rwdcb3]=spddir(u3cb,v3cb);
[wscb4 rwdcb4]=spddir(u4cb,v4cb);
[wscb5 rwdcb5]=spddir(u5cb,v5cb);
%convert back to wind speed and direction bias and slope correction
[wsc1 rwdc1]=spddir(u1c,v1c);
[wsc2 rwdc2]=spddir(u2c,v2c);
[wsc3 rwdc3]=spddir(u3c,v3c);
[wsc4 rwdc4]=spddir(u4c,v4c);
[wsc5 rwdc5]=spddir(u5c,v5c);
%correct true wind directions
wdc1=rwdc1+twr_orien;
hi=find(wdc1>360);
lo=find(wdc1<0);
wdc1(hi)=wdc1(hi)-360;
wdc1(lo)=wdc1(lo)+360;
wdc2=rwdc2+twr_orien;
hi=find(wdc2>360);
lo=find(wdc2<0);
wdc2(hi)=wdc2(hi)-360;
wdc2(lo)=wdc2(lo)+360;
wdc3=rwdc3+twr_orien;
hi=find(wdc3>360);
lo=find(wdc3<0);
wdc3(hi)=wdc3(hi)-360;
wdc3(lo)=wdc3(lo)+360;
wdc4=rwdc4+twr_orien;
hi=find(wdc4>360);
lo=find(wdc4<0);
wdc4(hi)=wdc4(hi)-360;
wdc4(lo)=wdc4(lo)+360;
wdc5=rwdc5+twr_orien;
hi=find(wdc5>360);
lo=find(wdc5<0);
wdc5(hi)=wdc5(hi)-360;
wdc5(lo)=wdc5(lo)+360;
%I suggest correcting ustars by the same proportion since u prime is
%probably in error by the same proportion as ws
%i.e. :
usc1=us1.*wsc1./ws1;
usc2=us2.*wsc2./ws2;
usc3=us3.*wsc3./ws3;
usc4=us4.*wsc4./ws4;
usc5=us5.*wsc5./ws5;
function [u,v] = uv(spd,dir)
% Uv converts Speed, Direction (met convention)
% to U and V components
%changed so that u,v fit met convention 4/12/99
%now u = speed from north; v = speed from east (90 deg)
%works with vectors
dir = dir + 180;
dirr = dir.*pi./180; % convert to radians
s = sin(dirr);
c = cos(dirr);
v = -spd.*s;
u = -spd.*c;
function [s,d] = spddir(u,v)
% Converts U,V components to speed, direction
% using met convention, works with vectors
%changed to correspond to uv function 4/12/99
%a bit messy but works on vectors
uu=v;
vv=-u;
i1 = find(vv>0);
i2 = find(vv<0 & -v<0);
i3 = find(vv<0 & -v>=0);
i4 = find(vv==0 & -v<0);
i5 = find(vv==0 & -v>0);
i6 = find(vv==0 & -v==0);
i7 = find(isnan(u) | isnan(v));
if length(i1) > 0
d(i1) = 180 - (180/pi).*atan(uu(i1)./vv(i1));
end
if length(i2) > 0
d(i2) = -(180/pi).*atan(uu(i2)./vv(i2));
end
if length(i3) > 0
d(i3) = 360 - (180/pi).*atan(uu(i3)./vv(i3));
end
if length(i4) > 0
d(i4) = 90 .* ones(1,length(i4));
end
if length(i5) > 0
d(i5) = 270 .* ones(1,length(i5));
end
if length(i6) > 0
d(i6) = 0 .* ones(1,length(i6));
end
s = sqrt(uu.*uu+vv.*vv);
if length(i7) > 0
s(i7) = NaN .* ones(1,length(i7));
d(i7) = NaN .* ones(1,length(i7));
end
d=d'; % want to end up with a column
]]>
|
|
| Outline Page |
Sonic Wind Page |
|
Last update: 5/4/99
Please send all comments and suggestions to the author, Peter
Guest,