Thank you Dennis,<br>I will try to use that, or maybe I will find an other quicker solution.<br><br>Giorgio<br>
<br>
----Messaggio originale----<br>Da: shea@ucar.edu<br>Data: 30-set-2015 17.52<br>A: "Giorgio Graffino"&lt;giorgio.graffino@alice.it&gt;<br>Ogg: Re: [ncl-talk] R: NCL uv2vr_cfd function on Python<br><br>Attached..<br><br>Not sure why you want to do python. Mimicking the NCL function could<br>be cumbersome.<br>Also, 'for/do' loops in any interpreted language can be slow.<br><br>Good Luck<br><br>On Wed, Sep 30, 2015 at 9:49 AM, Dennis Shea &lt;shea@ucar.edu&gt; wrote:<br>&gt; The underlying code used by NCL is a f77 subroutine that handles the<br>&gt; boundaries (half widths). Also, it allows both the lat and lon arrays<br>&gt; to be unequally spaced. ie<br>&gt;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; dlat_a = lat(2)-lat(0)&nbsp;&nbsp; , dlat_b = lat(10)-lat(8)<br>&gt;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; dlat_a does not need to be the same as dlat_b<br>&gt;<br>&gt; Allowing the above, makes the code a bit more involved than your simple script.<br>&gt;<br>&gt; The f77 code will be set to you<br>&gt;<br>&gt;<br>&gt;<br>&gt;<br>&gt;<br>&gt;<br>&gt;<br>&gt; On Tue, Sep 29, 2015 at 8:29 AM, Giorgio Graffino<br>&gt; &lt;giorgio.graffino@alice.it&gt; wrote:<br>&gt;&gt; Ok, I think I've found the problem: dx2 and dy2 are distance in meters, not<br>&gt;&gt; degrees.<br>&gt;&gt; Since u and v are given in latitude/longitude, there must be an implicit<br>&gt;&gt; conversion inside the uv2vr_cfd function.<br>&gt;&gt;<br>&gt;&gt; Anyone knows how this conversion is performed?<br>&gt;&gt;<br>&gt;&gt; Giorgio<br>&gt;&gt;<br>&gt;&gt;<br>&gt;&gt; ----Messaggio originale----<br>&gt;&gt; Da: giorgio.graffino@alice.it<br>&gt;&gt; Data: 29-set-2015 11.31<br>&gt;&gt; A: &lt;ncl-talk@ucar.edu&gt;<br>&gt;&gt; Ogg: [ncl-talk] NCL uv2vr_cfd function on Python<br>&gt;&gt;<br>&gt;&gt;<br>&gt;&gt; Dear NCL Users,<br>&gt;&gt; I'm trying to adapt an NCL script I wrote to compute water vertical speed on<br>&gt;&gt; Python.<br>&gt;&gt;<br>&gt;&gt; In particular, I'm finding some difficulties on "translating" the NCL<br>&gt;&gt; built-in function uv2vr_cfd into Python; I used that function to compute the<br>&gt;&gt; vertical component of wind stress. On the web page<br>&gt;&gt; (https://www.ncl.ucar.edu/Document/Functions/Built-in/uv2vr_cfd.shtml) it's<br>&gt;&gt; stated that<br>&gt;&gt;<br>&gt;&gt; According to H.B. Bluestein [Synoptic-Dynamic Meteorology in Midlatitudes,<br>&gt;&gt; 1992, Oxford Univ. Press p113-114], let D represent the partial derivative,<br>&gt;&gt; a the radius of the earth, phi the latitude and dx2/dy2 the appropriate<br>&gt;&gt; longitudinal and latitudinal spacing, respectively. Then, letting j be the<br>&gt;&gt; latitude y-subscript, and i be the longitude x-subscript:<br>&gt;&gt;<br>&gt;&gt;&nbsp;&nbsp;&nbsp;&nbsp; rv = Dv/Dx - Du/Dy + (u/a)*tan(phi)<br>&gt;&gt;<br>&gt;&gt;&nbsp;&nbsp;&nbsp;&nbsp; rv(j,i) = (v(j,i+1)-v(j,i-1))/dx2(j)<br>&gt;&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; - (u(j+1,i)-u(j-1,i))/dy2(j)<br>&gt;&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; + (u(j,i)/a)*tan(phi(j))<br>&gt;&gt;<br>&gt;&gt; Now, I tried to write that function on Python in this way<br>&gt;&gt;<br>&gt;&gt; def Curl(u,v,nlat,nlon):<br>&gt;&gt;<br>&gt;&gt; rv = np.empty(shape=(nlat,nlon))<br>&gt;&gt;<br>&gt;&gt; rv.fill(np.nan)<br>&gt;&gt;<br>&gt;&gt; M0 = np.arange(1,nlat-1) # M0 = [1,...,nlat-1]<br>&gt;&gt;<br>&gt;&gt; N0 = np.arange(1,nlon-1) # N0 = [1,...,nlon-1]<br>&gt;&gt;<br>&gt;&gt; for m in M0:<br>&gt;&gt;<br>&gt;&gt; for n in N0:<br>&gt;&gt;<br>&gt;&gt; rv[m,n] = (v[m,n+1]-v[m,n-1])/(lons[n+1]-lons[n-1]) - \<br>&gt;&gt;<br>&gt;&gt; (u[m+1,n]-u[m-1,n])/(lats[m+1]-lats[m-1]) + \<br>&gt;&gt;<br>&gt;&gt; (u[m,n]/radius)*tans[m]<br>&gt;&gt;<br>&gt;&gt; return rv<br>&gt;&gt;<br>&gt;&gt;<br>&gt;&gt; but, even if the pattern is similar, there are four orders of magnitude<br>&gt;&gt; between the values. I attach results using NCL (MedDailyINDEX.png file) and<br>&gt;&gt; Python (WindStressCurl.png) for the same day.<br>&gt;&gt;<br>&gt;&gt; Could it be possible to know how exactly rv is computed inside uv2vr_cfd<br>&gt;&gt; function? Because I think there is something missing in my implementation...<br>&gt;&gt;<br>&gt;&gt; Thanks for your attention.<br>&gt;&gt;<br>&gt;&gt; Regards,<br>&gt;&gt; Giorgio<br>&gt;&gt;<br>&gt;&gt;<br>&gt;&gt;<br>&gt;&gt;<br>&gt;&gt;<br>&gt;&gt; _______________________________________________<br>&gt;&gt; ncl-talk mailing list<br>&gt;&gt; ncl-talk@ucar.edu<br>&gt;&gt; List instructions, subscriber options, unsubscribe:<br>&gt;&gt; http://mailman.ucar.edu/mailman/listinfo/ncl-talk<br>&gt;&gt;<br><br><br>