[ncl-talk] Discrepancies between Python and NCL results

[Bill Ladwig]

Hi Guido,

Just at a quick glance, I see this in your NCL script:

w=data->w(:,1:150,ind_dry)
>

In Python, you have this:

 w=np.array(nc2.variables["w"][:][:,0:150,ind_dry])


You've got difference slices going (both NCL and Python use 0-based
indexes), was that intentional?

Bill


On Wed, Feb 21, 2018 at 2:58 AM, Guido Cioni <guidocioni at gmail.com> wrote:

> Hey all,
> I have a weird question addressed to experts of both NCL and Python.
>
> In my effort to speed up a NCL script which computes the Moist Static
> Energy of updraft points and of the environment in a LES model simulation
> of diurnal convection, I've tried to translate the entire script into
> Python. I was interested to see whether Python was faster than NCL: long
> story short, it's not :) The problem is that my data is really big (72 x
> 150 x 1280000 points) so that a single assignment a = file->var takes
> approx. 40 minutes. A single plot takes about 3 hours to be produced. And
> I'm using no single loop! :)
>
> But let's come to the problem. I have to select updrafts, which I defined
> as grid points having a vertical velocity of at least 1 m/s and a cloud
> water content of at least 1E-4 kg/kg. For these points I want to construct
> a vertical profile of moist static energy. Thus, it seems to me that using
> the where function was going to be my best shot. This is what I tried to do
> to select the grid points:
>
>    w_cloud=where(qc.gt.th_qc .and. w.ge.1, mse, mse at _FillValue)
>    w_envir=where(qc.le.th_qc .or. w.lt.1, mse, mse at _FillValue)
>
> And then averaged over time and space to get a single vertical profile (I
> know that NCL will take care of the missing values using these operations).
> The result is not bad..
>
>
> The red line shows the difference w_cloud-w_envir which is, as expected,
> positive meaning that updrafts are moister and warmer than the surrounding.
>
> Coming to Python I decided to use a similar approach with numpy masked
> arrays
>
> mse_cloud_3D=np.ma.masked_where( ((w<1) | (qc<1E-04)), mse)
> mse_cloud=np.ma.mean(mse_cloud_3D, axis=(0,2))
>
> mse_envir_3D=np.ma.masked_where( ((w>=1) & (qc>=1E-04)), mse)
> mse_envir=np.ma.mean(mse_envir_3D, axis=(0,2))
>
> For the average MSE the profile looks almost exactly the same (good!) but
> the difference instead is completely different, also considering the order
> of magnitude.
>
> So my question is...what I'm doing wrong? How could the average MSE be
> fine and the difference between the updraft and the environment not? I
> would tend to trust more NCL as the masked arrays of numpy are behaving
> pretty weird sometimes.
>
> I'm attaching the scripts cleaned as much as possible, although it is
> impossible (for obvious reasons) to send the data along :)
>
> If anyone has opinion/comments...
>
> Cheers
>
>
>
> Guido Cioni
> http://guidocioni.altervista.org
>
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20180221/9d4fef41/attachment.html>