function [h_mean,ws_mean,wd_mean,u_mean,v_mean]=mean_profile(hurricane)
hurricane='frances'  %'swp.104090504*.dat';
rep='modified\vert\';
radar_files=dir(['C:\Documents and Settings\lorsolo\My Documents\research\',hurricane,'\data\',rep,'swp.104090500*.dat']);
ws=[];
wd=[];
h=[];
aaa=15;
bbb=400;
L=length(radar_files)
for i=1:length(radar_files)
    data=load(['C:\Documents and Settings\lorsolo\My Documents\research\',hurricane,'\data\total_wind\',radar_files(i).name(1:length(radar_files(i).name)-4),'_wind.dat']);
    el=str2num(radar_files(i).name(30:32));
    if el~=0.5 &el~=0.8
        ws=[ws;data(aaa:bbb,1)];
        wd=[wd;data(aaa:bbb,2)];
        rge=33.35:66.7:399*66.7+33.35;
        rge=rge(aaa:bbb);
        
        h=[h;rge'.*sin(deg2rad(el))];
    end
end

[h,index]=sort(h);
ws=ws(index);
wd=wd(index);

u_sl=find(h<=260);
ws_sl=nanmean(ws(u_sl))
[u_sl,v_sl]=wind_decomposition(wd(u_sl));
u_sl=nanmean(u_sl);
v_sl=nanmean(v_sl);
[wd_sl]=wind_direction_deg(u_sl,v_sl)

ws_mean=[];
wd_mean=[];
u_mean=[];
v_mean=[];
h_mean=[];

segment=0:20:800;
for j=1:length(segment)-1%0:299%1:length(segment)-1
    u_segment=find(h > segment(j) & h <= segment(j+1));
    %u_segment=j*L+1:(j+1)*L;
    temp=ws(u_segment);
    u_NaN=find(~isnan(temp));
    temp_mean=mean(temp(u_NaN));
    sigmaWS=std(temp(u_NaN));
    outliers = abs(temp - temp_mean) > 1.5*sigmaWS;
    u_outl=find(outliers==1);
    temp(u_outl)=NaN;
    u_NaN=find(~isnan(temp));
    ws_mean=[ws_mean;mean(temp(u_NaN))];
    h_mean=[h_mean;(segment(j)+segment(j+1))/2];
    %h_mean=[h_mean;h(j*L+1)];
    temp=wd(u_segment);
    temp(u_outl)=NaN;
    u=[];
    v=[];
    ang=[];
    alpha=[];
    for i=1:length(temp)
        ang(i)=(temp(i)*pi)/180;
        if ang(i)<=pi/2
            u(i)=-sin(ang(i));
            v(i)=-cos(ang(i));
            
            
        elseif ang(i) > pi/2 & ang(i) <= pi
            ang(i)=ang(i)-pi/2;
            u(i)=-cos(ang(i));
            v(i)=sin(ang(i));
        elseif ang(i) > pi & ang(i) <= 3*pi/2
            ang(i)=ang(i)-pi;
            u(i)=sin(ang(i));
            v(i)=cos(ang(i));
        elseif ang(i) > 3*pi/2 & ang(i) <=2* pi
            ang(i)=ang(i)-3*pi/2;
            u(i)=cos(ang(i));
            v(i)=-sin(ang(i));
        else
            u(i)=NaN;
            v(i)=NaN;
        end
        
    end
    u_NaN=find(~isnan(u));
    if length(u_NaN)==0
        wd_mean=[wd_mean;NaN];
        u_mean=[u_mean;NaN];
        v_mean=[v_mean;NaN];
    else
        U=mean(u(u_NaN));
        V=mean(v(u_NaN));
        
        sigmaU=std(u(u_NaN));
        sigmaV=std(v(u_NaN));
        outliers = abs(u - U) > 1.5*sigmaU;
        u_outl=find(outliers==1);
        u(u_outl)=NaN;
        v(u_outl)=NaN;
        outliers = abs(v - V) > 1.5*sigmaV;
        u_outl=find(outliers==1);
        u(u_outl)=NaN;
        v(u_outl)=NaN;
        u_NaN=find(~isnan(u));
        U=mean(u(u_NaN));
        V=mean(v(u_NaN));
        u_mean=[u_mean;U];
        v_mean=[v_mean;V];
        
        if U >0 & V >0
            alpha=atan(U/V);
            alpha=pi+alpha;
            
        elseif U >0 & V <0
            alpha=-atan((V/U));
            alpha=3*pi/2+alpha;
            
        elseif U <0 & V >0
            alpha=-atan((V/U));
            alpha=pi/2+alpha;
            
        elseif U <0 & V <0
            alpha=atan(U/V);
        else
            alpha=NaN;
        end
        alpha=rad2deg(alpha);
        wd_mean=[wd_mean;alpha];
    end
end

