<html><head></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; ">Interesting question! From the man page on 'pow':&nbsp;pow(x, y) returns a NaN and raises the "invalid" floating-point exception for finite x &lt; 0 and finite non-integer&nbsp;y.<div><div><br></div><div>Here are some GrADS calls that illustrate this functionality of 'pow':</div><div><div>ga-&gt; d pow(2,0.5)&nbsp;</div><div>Result value = 1.41421&nbsp;</div><div>ga-&gt; d pow(-2,2) &nbsp;</div><div>Result value = 4&nbsp;</div><div>ga-&gt; d pow(-2,-1)</div><div>Result value = -0.5&nbsp;</div><div>ga-&gt; d pow(-2,0.9999999999)</div><div>Result value = nan&nbsp;</div><div>ga-&gt; d pow(-2,1.1) &nbsp; &nbsp; &nbsp; &nbsp;&nbsp;</div><div>Result value = nan&nbsp;</div><div><br></div></div><div>The problem is that GrADS is returning the result of the call to 'pow' without checking whether it is a NaN and turning it into a proper missing data value. I will add this to my list of bugs to fix. Meanwhile, wrap a maskout around the expression 1000/(lev-e) and that should any negative values for the first argument to 'pow'.&nbsp;</div><div>--Jennifer</div><div><br></div><div><br></div><div><br><div><div>On Jan 18, 2013, at 1:21 PM, Jeff Duda wrote:</div><br class="Apple-interchange-newline"><blockquote type="cite">Hello,<br>I am trying to compute some thermodynamic parameters, shown below.<br clear="all"><br>&nbsp;'define e = rhprs*6.11*exp((2.5e6/461.5)*((1/273)-(1/tmpprs)))'<br>&nbsp;'define r = 622*e/(lev-e)'<br>&nbsp;'define condtemp = 55 + 1 / (1/(tmpprs-55) - log(rhprs/100)/2840)'<br>
&nbsp;'define thetaDL = tmpprs*pow(1000/(lev-e),0.2854)*pow(tmpprs/condtemp,0.0028*r)'<br>&nbsp;'define thetae = thetaDL*exp(((3.036/condtemp)-0.00178)*r*(1+0.448+r/1000))'<br>&nbsp;'define thetaw = 45.114 - 51.489*pow(thetae/273.15,-3.504)'<br>
&nbsp;'define LSI = max(thetaw,lev=1000,lev=500) - max(thetaw,lev=1000,lev=900)'<br><br>I get a lot of NaNs.&nbsp; They are coming from this line:<br>&nbsp;'define thetaDL = tmpprs*pow(1000/(lev-e),0.2854)*pow(tmpprs/condtemp,0.0028*r)'<br>
<br>And in fact, from this piece:<br>pow(1000/(lev-e),0.2854)<br><br>I know this because I tested this for a model grid column (from 1000 - 500 mb). <br>ga-&gt; d 1000/(lev-e)&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; <br>Printing Grid -- 13 Values -- Undef = -9.99e+08<br>
-3.55595 -4.04332 -4.69364 -9.98891 39.3747 -5.3341 -11.5849 4.19346 <br>2.07427 2.42358 2.11214 2.1287 2.24423 <br>ga-&gt; d pow(1000/(lev-e),0.2854)<br>Printing Grid -- 13 Values -- Undef = -9.99e+08<br>-nan -nan -nan -nan 2.85282 -nan -nan 1.5055 <br>
1.2315 1.28743 1.23787 1.24063 1.25949 <br><br>For whatever reason, the pow() function doesn't like negative bases, because for example, I can use a calculator to compute -3.55595^0.2854 as -1.43629.&nbsp; Why, then, does the pow() function give a -nan as a result?<br>
<br>Jeff Duda<br><br>-- <br>Jeff Duda<br>Graduate research assistant<br>University of Oklahoma School of Meteorology<br>Center for Analysis and Prediction of Storms<br>
_______________________________________________<br>gradsusr mailing list<br><a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>http://gradsusr.org/mailman/listinfo/gradsusr<br></blockquote></div><br><div apple-content-edited="true">
<span class="Apple-style-span" style="font-size: 12px; "><div>--</div><div>Jennifer M. Adams</div><div>IGES/COLA</div><div>4041 Powder Mill Road, Suite 302</div><div>Calverton, MD 20705</div><div><a href="mailto:jma@cola.iges.org">jma@cola.iges.org</a></div><div><br class="khtml-block-placeholder"></div><br class="Apple-interchange-newline"></span>
</div>
<br></div></div></body></html>