<div dir="ltr">Hello,<div><br></div><div>I believe that the NCL function wgt_volrmse has a bug -- namely, that the formula is wrong. I am using NCL V6.3.0, but downloaded the source code of V6.4.0 and confirmed that the relevant source code is unchanged. I also could not find any reports on a bug in wgt_volrmse on the NCL website.</div><div><br></div><div>I first became suspicious of a bug when my code using wgt_volrmse were an order of magnitude different than results that a colleague of mine obtained using raw Fortran code.</div><div>So I downloaded the source code for versions 6.3.0 and 6.4.0 to look at wgt_volrmse. The relevant fortran file is under ./ncl_ncarg-6.4.0/ni/src/lib/nfpfort/volRmse.f</div><div><br></div><div>I looked at this code and it appears that a square is being taken twice:</div><div>First, SUMTQ = (T - Q)**2</div><div>Then, VARTQ = (SUMTQ/SUMWZ)**2</div><div>A square root is taken at the end. So that "negates" only one of the **2 operations.</div><div><br></div><div><br></div><div>I then looked at the code for wgt_arearmse. In theory, these two functions should provide identical answers if the third-from-rightmost dimension is of size 1 and no weighting is applied to this dimension in wgt_volrmse.</div><div>The relevant source code for wgt_arearmse is under ./ncl_ncarg-6.4.0/ni/src/lib/nfpfort/areaRmse.f</div><div><br></div><div>This code takes a square once (SUMD) then takes a square root at the end. So I have no reason to believe that wgt_arearmse has a bug.</div><div><br></div><div><br></div><div>As I said earlier, in theory wgt_volrmse and wgt_arearmse should provide identical answers if the third-from-rightmost dimension is of size 1 and no weighting is applied to this dimension in wgt_volrmse.</div><div>I tested this out with a script, and with data samples. I upload a .tar file (demo_RMSEbug.tar, <b>filesize 530K</b>) via ftp. The contents of this .tar file are:</div><div><a href="http://fcst.nc">fcst.nc</a> (a file of size (1, 181 lat, 360 lon) of forecast 2-m temperatures): <b>filesize 260K</b></div><div><a href="http://obs.nc">obs.nc</a> (a file of size (1, 181 lat, 360 lon) of observed 2-m temperatures): <b>filesize 261K</b></div><div>demo_bug.ncl (a small sample script): <b>filesize 1.5K</b></div><div><br></div><div>In my test script, I calculate RMSE using both wgt_arearmse and wgt_volrmse, and print the output. As might be expected with the extra **2 operation in wgt_volrmse, the RMSE I get with this function is much higher than the RMSE I get when using wgt_arearmse. Values from both are printed out within the code.</div><div><br></div><div>I believe that this conclusively demonstrates a bug in wgt_volrmse (or at least an inconsistency with wgt_arearmse). If not, I apologize.</div><div><br></div><div>Other relevant information:</div><div>ncl -V returns 6.3.0</div><div>uname -a returns <span style="font-variant-ligatures:no-common-ligatures;font-family:Menlo;font-size:11px">Linux tfe10 3.10.0-693.11.1.el7.x86_64 #1 SMP Fri Oct 27 05:39:05 EDT 2017 x86_64 x86_64 x86_64 GNU/Linux</span></div>
<div><br></div><div><br></div><div>Please let me know if you have any further questions.</div><div><br></div><div>Thank you very much for all of the work you do for NCL -- it is a great tool.</div><div><br></div><div>Sincerely,</div><div><br></div><div>Benjamin Green</div><div><br></div><div>-- <br><div class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr">Benjamin Green<div>CIRES Research Associate at NOAA/ESRL/GSD</div><div>David Skaggs Research Center, Room 2B-501</div><div>Phone: 303-497-5132</div></div></div></div></div></div></div>
</div></div>