function [out_val]=interp_wrf(in_val,in_hgt,out_hgt,nl)

zk=in_hgt;
p=in_val;
hgtnow=out_hgt;

if hgtnow> zk(nl); % EXTRAPOLATE ALOFT! if the model surface is greater than the top of sounding
    %disp('extrap aloft');
    w2 = (zk(nl)-hgtnow)./(zk(nl)-zk(nl-1)); % how many zk levels fit into the excess height
    w1 = 1- w2;
    p_surf=(w1.*(p(nl)))+w2.*(p(nl-1));
elseif hgtnow < zk(1) % EXTRAPOLATE BELOW
    %disp('extrap below'); 
    w2 = (zk(2)-hgtnow)./(zk(2)-zk(1)); % how many zk levels fit into the excess height
    w1 = 1- w2;
    p_surf=(w1.*(p(2)))+w2.*(p(1));
else % we need to search through the column to find where the ground resides
    intp='false';
    kp=nl;
    while (strcmp(intp,'false') && kp >=2)
        if (zk(kp)>= hgtnow && zk(kp-1) <= (hgtnow)) % the ground is between these layers
            w2=(hgtnow-zk(kp))./(zk(kp-1)-zk(kp));
            w1=1-w2;
            p_surf=w1.*(p(kp))+w2.*(p(kp-1));
            intp='true';
            %disp('Interp'); 
        end
        kp=kp-1; % step down through the column until we find the ground
    end
end


out_val=p_surf;

end