Fixing history:

<ol>
<li>[2004/06/21] Make fcurl actually do the curl (vrtec) not inverse. Make rws fit this.
<li>[ditto] Swap u,v in various routines. The spherepack stuff uses v,u which I find odd.
<li>[2004/06/22] Make fdiv call div1 not laplacian :-))))
</ol>

<h1>Now, just the RWS tested</h1>

OK, lets try to see if the great <a href=rws.pro>rws</a> code works.

<p>
We'll be guided by James, p269. That shows multi-y ECMWF for DJF. We'll just use Jan, 850 hPa winds.

<pre>
; We want the 200 hPa field (well, 150 ideally...)
;uv=gm('era-40/0.01',m=1,y=1979+indgen(11),sc=15201+[0,1],lev=200)                          
;uv2=reform(uv,11,2)
;uv3=pp_avg1(uv2,1)                                                                         
;uv4=pp_regrid(red2reg(uv3),get_orog(),/pole)                                               
;pp_write,uv4,'era-40-uv-200-hpa-79_89'                                                     

uv1=get_field('era-40-uv-200-hpa-79_89',nf=2)

; pp_plot,/glob,uv1,/vec,skip=[4,3],col=3,lev=indgen(8)*2+1,th=2
pp_plot,/glob,uv1,/vec,skip=[4,3],col=3,lev=(indgen(10)+1)*5,th=2,model=0   

gettwogifs,out='s1'
</pre>

<a href=s1.png><img src=s1_s.png></a>

<p>
So, now lets generate the RWS etc terms:
<pre>
z=rws(uv1,div=div,sfvp=sfvp,vort=vort,abs_vort=abs_vort,pvort=pvort,rot=rot)

pp_plot,abs_vort,/glob,model=0,lev=(indgen(15)-7)*2e-5
pp_arrow,div,skip=[4,3],col=1,th=2
gettwogifs,out='james-fig-8.6-a'
</pre>

First, lets look at the abs vorticity and the vectors of the divergent wind:
<p>
<a href=james-fig-8.6-a.png><img src=james-fig-8.6-a_s.png></a>
<p>
This looks correct at last! Hurrah. The Abs vort has the right magnitude and structure, though a bit more structure than James
gives - presumably his was contoured by hand.

<pre>
pp_plot,/glob,wpp_fc(z,!met_values.radius_earth,2),model=0,lev=(indgen(15)-7)*5e-11,yr=[-80,80]
gettwogifs,out='james-fig-8.6-b'
</pre>

<a href=james-fig-8.6-b.png><img src=james-fig-8.6-b_s.png></a>

<p>
This at last looks right: it has the correct magnitudes, and it has the correct pattern, noce
you allow for James's human contouring. Hurrah!


<h4>Case 2</h4>

(see also <a href=sfvp.html>sfvp</a>)
<p>
Vera at al plot the RWS term (the diff of WC and WW cases) in the S Pacific in their fig 8 b:

<P>
<pre>
rwsc=rws(uvc)
rwsw=rws(uvw)
pp_plot,/glob,wpp_fc(pp_diff(rwsc,rwsw),!met_values.radius_earth,2),model=0,lev=(indgen(15)-7)*5e-11,yr=[-80,80],l0=180,yr=[-60,20],xr=[90,300]
gettwogifs,out='vera-fig-8-b'
</pre>
<a href=vera-fig-8-b.png><img src=vera-fig-8-b_s.png></a>
<p>
This substantially replicates their figure. I have used mean winds; they appear to say that they use time_avg(rws(winds)) whereas
I am using rws(time_avg(winds)).

<p>
We can try using the other averaging (though its slow...)
<pre>
u1=pp_regrid(/pole,red2reg(u=gm('era-40/0.01',m=mm,sc=15201,y=y,lev=200)),get_orog())   
v1=pp_regrid(/pole,red2reg(u=gm('era-40/0.01',m=mm,sc=15202,y=y,lev=200)),get_orog())   
rwss=u1
for i=0,n_elements(u1)-1 do rwss(i)=rws([u1(i),v1(i)])   

s1=pp_avg(rwss(is_in(wc,rwss.lbyr)))                                                                                           
s2=pp_avg(rwss(is_in(ww,rwss.lbyr)))
cw=pp_diff(s1,s2)                                                                                                              
pp_plot,/glob,wpp_fc(cw,!met_values.radius_earth,2),model=0,lev=(indgen(15)-7)*5e-11,yr=[-80,80],l0=180,yr=[-60,20],xr=[90,300]
gettwogifs,out='vera-fig-8-b1'

</pre>
<a href=vera-fig-8-b1.png><img src=vera-fig-8-b1_s.png></a>
<p>
But this looks much the same as above... almost disturbingly similar... fortunately its not identical. A fairly linear process, then.