NCEP fnl

zlpyl

由于fnl资料不提供各层的比湿，最近用NCEP 每日4次的fnl资料计算了各层比湿，公式看了很久没发现问题，但计算出的比湿，就是比实况小很多，不存在数量级问题，比如700hpa基本计算的比实况少一半。有没有计算过的，知道这是什么问题吗，附上计算比湿的公式
do i=0,20
    if(tmp1(i).lt.258.16)then
     es1(i)=6.1078*exp(21.8745*(tmp1(i)-273.16)/(tmp1(i)-7.66))
    else
     es1(i)=6.1078*exp(17.27*(tmp1(i)-273.16)/(tmp1(i)-35.86))
    end if
  
   qs1(i)=0.62197*es1(i)/(lev(i)/100-0.378*es1(i))
   rs1(i)=0.622*es1(i)/(lev(i)/100-es1(i))
   

  q1(i)=rh1(i)*qs1(i)*1000/100
  end do