Re: confused on the issue of calculation result using NCL

From: Mary Haley <haley_at_nyahnyahspammersnyahnyah>
Date: Mon Aug 05 2013 - 10:12:20 MDT

Thomas,

[Sorry that this question slipped through the cracks.]

I'm having trouble debugging your code, because there are a *lot* of calculations and I'm not sure if the units are all being handled correctly. It would take me awhile to go through every line of your code to determine if each function is being called correctly.

Some things of note:

1. I'm getting a warning that says your data is already of type "float" , so I would recommend not using "short2float".

2. You are setting "_FillValue" in some places that seems unnecessary. For example:

U@_FillValue = -999
V@_FillValue = -999
T@_FillValue = -999

Fp_int@_FillValue = 0
Fphi_int@_FillValue = 0
EPdiv_int@_FillValue = 0

These six variables already have missing values, so there's no need to reset them.

Also, in general, I would highly caution against setting any _FillValue to 0. It's better to use a value that you know is far outside the range of the expected values. One possibility is to use the "default_fillvalue" function. So, for example, instead of:

x@_FillValue = 0

use:

x@_FillValue = default_fillvalue(typeof(x))

3. Why are you doing the following:

U = where(ismissing(U), 0.0, U)
V = where(ismissing(V), 0.0, V)
T = where(ismissing(T), 0.0, T)

? It doesn't seem to make a difference in the calculations that follow, but I'm not sure you want to do this.

4. You have some code like this:

THETA = T*(conform(T,level,1)/1000)^(-0.286)

To the casual observer, I'm not sure what you are doing here, and what units are expected.
It would help if you commented your code, and/or immediately add a "units" attribute so
you can tell what this calculation is supposed to be doing.

5. Since you are creating a lot of variables in this script, my suggestion is to do a "printMinMax" on all of these the minute they are created or calculated, so you can be sure everything looks correct.

I like using "printMinMax", but this function is more useful if you have a "long_name" attribute, so I created a "printMinMax2", which uses the name of the variable, instead of the "long_name" attribute. Also, if the variable contains a "units" attribute, then this is included. I included this modified script in your code so you can use it.

See the attached "script_mod.ncl" and the output that I got:

(0) U min/max: -154.313/151.585
(0) V min/max: -151.698/142.518
(0) T min/max: 183.654/308.712
(0) THETA min/max: 229.309/2124.05
(0) THETAzm min/max: 235.924/2110.06
(0) loglevel min/max: 0/6.90776
(0) THETAptemp min/max: -5.9165/0.00181306
(0) THETAp min/max: -5.9165/0.00181306
(0) THETAptemp min/max: -5.9165/0.00181306
(0) Uza min/max: -130.228/135.084
(0) Vza min/max: -144.9/148.049
(0) THETAza min/max: -248.343/252.288
(0) UV min/max: -9593.2/11725.2
(0) UVzmtm min/max: -180.951/1145.19
(0) VTHETA min/max: -7391.56/32627.4
(0) VTHETAzmtm min/max: -133.385/7841.23
(0) THETAptm min/max: -5.9165/0.00181306
(0) phi min/max: -1.5708/1.5708
(0) asphi min/max: -6.37122e+06/6.37122e+06
(0) acphi min/max: -0.278495/6.37122e+06
(0) omega min/max: 7.2921e-05/7.2921e-05
(0) f min/max: -0.000145842/0.000145842
(0) latfac min/max: 1.21734e-08/6.37122e+06
(0) Fphi min/max: -2.32857e+09/3.85735e+08
(0) Fp min/max: -1.16072e+09/9.11285e+09
(0) EPdiv1 min/max: -1581.09/2822.4
(0) EPdiv2 min/max: -1.21428e+06/795020
(0) EPdiv min/max: -1.21429e+06/795019
(0) lev (millibar) min/max: 1/1000
(0) lat_int (degrees_north) min/max: -87.5/87.5
(0) Fphi_int2 min/max: -2.32786e+09/3.37465e+08
(0) Fp_int2 min/max: -1.16072e+09/3.60732e+08
(0) EPdiv_int2 min/max: -149035/54246.5
(0) EPdiv_int min/max: -149035/54246.5
(0) dudt_int min/max: -2991.55/1088.88

This is not all your variables, but it may be a good start. Note that many of these don't have any units, so I'm not sure if these are all correct.

Some possible suspicions might be your asphi/acphi variables. The scales seem off.

--Mary

> On Wed, Jul 31, 2013 at 10:59 AM, Thomas Tobian <thomastobian@gmail.com> wrote:
> Dear NCL Community,
>
> I recently created my own codes using NCL to calculate flux.
> Everything was OK, the codes just run well without any errors.
> However the results were not realistic, it was twice bigger than I was expecting, for instance
> the result of (dudt) is with magnitude of ten to the power 2, in reality is normally less than 80 m/s2 .
> I did many time debugging procedures and could not find what exactly the problem.
> I have already made sure that all the type of input are float and there are no missing values (has been interpolated).
>
> For this, I expect a help from the NCL experts to verify what exactly the problem that
> happened in my input dataset. I have already put the "script.ncl, input.nc and output.nc" in ftp
> incoming by name "comp.tar.gz", extracted this *tar, and run : "ncl<script.ncl" you will directly
> find the output in a second. I just wonder whether this issue related to the type of coordinate, or file structures etc.
> I really thank you for your suggestion.
>
> Best regards,
> Thomas Tobian
>
>
>

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk

Received on Mon Aug 5 10:12:31 2013

This archive was generated by hypermail 2.1.8 : Wed Aug 07 2013 - 13:59:00 MDT