<div dir="ltr"><div>If your forecast and verification data are on the same grid, you can use the maskout and const functions to construct fields of ones and zeros and then use the asum function to add them up to give you a contingency table. I believe it can be done in a relatively small number of lines. For example, to compute the correct yes forecasts (i.e., hits), you might use<br>
<br></div>*For each threshold<br><div>'define binforecast = const(const(maskout(rain,rain-threshold),1),0,-u)'<br></div><div>'define binobs = const(const(maskout(prec,prec-threshold,1),0,-u)'<br></div><div>
'define a = const(binforecast+binobs,1)'<br></div><div>'define hits = asum(a,[lon/lat range])'<br><br></div><div>The first two lines do the same thing, one for the forecast rain, the other for the verification (which I've renamed observed) rain. Since that first line is a bear, I'll go through it piece by piece:<br>
</div><div>***maskout(rain,rain-threshold)<br></div><div>The maskout function returns the rain field, except values that fall below a threshold (that you pick) are converted to missing data. We want missing data because those will be ignored in final calculations.<br>
</div><div>***const(maskout(rain-rain-threshold),1)<br></div><div>The inner nested const function takes the output from the maskout function and converts all values that are NOT missing to 1. It does not touch the undefined data.<br>
***const(const(maskout(rain,rain-threshold),1),0,-u)<br></div><div>The outer const function isn't technically necessary, but it completes the construction of a binary field (just ones and zeros) of the precip field where the value is 1 at all grid points where the forecast rain exceeds the threshold and zero otherwise. This const function takes the output from the previous const function and converts all missing values (the flag -u) to zero.<br>
<br></div><div>We do the same thing for the observed precip field.<br><br></div><div>Now for the third line. I want to create a field that is basically 1 for hits and zero (or undefined) otherwise.<br>Since I have two binary fields, I can just add them. Where there is a hit, I will get 2. When adding fields, adding a real value to an undefined value at a grid point will turn the whole grid point undefined. So this field really computes an intersection of areas. I use the const function to convert the remaining valid values to 1. Finally, I can compute the total number of hits in the domain by using the asum function. Make sure you don't use atot() because the latter function incorporates latitude-weighting, which will result in a non-integer sum that isn't representative of your forecast result.<br>
<br></div><div>I wanted to create a true binary field (strictly ones and zeros). If I changed the valid values to 0.5 instead of to 1 in the first two lines, I wouldn't really need the third line of code. But this was just to illustrate the process.<br>
<br></div><div>Now use the same general algorithm to code misses, false alarams, and correct "no"s.<br><br></div><div>Good luck.<br><br>Jeff Duda<br></div></div><div class="gmail_extra"><br><br><div class="gmail_quote">
On Thu, May 1, 2014 at 5:39 AM, Amulya Chevuturi <span dir="ltr"><<a href="mailto:amulya.chevuturi@gmail.com" target="_blank">amulya.chevuturi@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div dir="ltr">Dear Grads Users,<div><br></div><div>Recently I asked a question about correction in a script building contingency tables (Asnwered by Dr. Eric Altshuler). I am attaching the script at the end of the mail.</div>
<div><br></div><div>But this script is for when the model and observation temporal and spatial extent are same. Even the model and observations resolution is same for this script to run. Can anyone tell me how to build 2*2 contingency tables for when the model and observations datasets are at different resolutions and extent. </div>
<div><br></div><div>With regards,</div><div><div style="font-family:arial,sans-serif;font-size:13px">Amulya Chevuturi</div><div style="font-family:arial,sans-serif;font-size:13px">School of Environmental Sciences</div><div style="font-family:arial,sans-serif;font-size:13px">
Jawaharlal Nehru University</div><div style="font-family:arial,sans-serif;font-size:13px">New Delhi, India</div></div><div><br></div><div><br></div><div>******Script for contingency table*************</div><div><div><div style="font-family:arial,sans-serif;font-size:13px">
'reinit'</div><div style="font-family:arial,sans-serif;font-size:13px"><br></div><div style="font-family:arial,sans-serif;font-size:13px">'open model.ctl'</div><div style="font-family:arial,sans-serif;font-size:13px">
'set t 3'</div><div style="font-family:arial,sans-serif;font-size:13px">'set y 57 208'</div><div style="font-family:arial,sans-serif;font-size:13px">'set x 72 228'</div><div style="font-family:arial,sans-serif;font-size:13px">
'define a = rain'</div><div style="font-family:arial,sans-serif;font-size:13px">'close 1'</div><div style="font-family:arial,sans-serif;font-size:13px"><br></div><div style="font-family:arial,sans-serif;font-size:13px">
'open observation.ctl'</div><div style="font-family:arial,sans-serif;font-size:13px">'set t 3'</div><div style="font-family:arial,sans-serif;font-size:13px">'set y 57 208'</div><div style="font-family:arial,sans-serif;font-size:13px">
'set x 72 228'</div><div style="font-family:arial,sans-serif;font-size:13px">'define b = prec'</div><div style="font-family:arial,sans-serif;font-size:13px"><br></div><div style="font-family:arial,sans-serif;font-size:13px">
n1=0</div><div style="font-family:arial,sans-serif;font-size:13px">n2=0</div><div style="font-family:arial,sans-serif;font-size:13px">n3=0</div><div style="font-family:arial,sans-serif;font-size:13px">n4=0</div><div style="font-family:arial,sans-serif;font-size:13px">
j=72<br></div><div style="font-family:arial,sans-serif;font-size:13px">while (j<=228)</div><div style="font-family:arial,sans-serif;font-size:13px">'set y 'j</div><div style="font-family:arial,sans-serif;font-size:13px">
i=57</div><div style="font-family:arial,sans-serif;font-size:13px">while (i<=208)</div><div style="font-family:arial,sans-serif;font-size:13px">'set x 'i</div><div style="font-family:arial,sans-serif;font-size:13px">
'd a'</div>c=subwrd(result,4)<div style="font-family:arial,sans-serif;font-size:13px">'d b'</div><div style="font-family:arial,sans-serif;font-size:13px"><span style="font-family:arial;font-size:small">d=subwrd(result,4)</span><br>
</div><div style="font-family:arial,sans-serif;font-size:13px">if(c>0&d>0)</div><div style="font-family:arial,sans-serif;font-size:13px">n1=n1+1</div><div style="font-family:arial,sans-serif;font-size:13px">else</div>
<div style="font-family:arial,sans-serif;font-size:13px">if(c<=0&d<=0)</div><div style="font-family:arial,sans-serif;font-size:13px">n2=n2+1</div><div style="font-family:arial,sans-serif;font-size:13px">else</div>
<div style="font-family:arial,sans-serif;font-size:13px">if(c<=0&d>0)</div><div style="font-family:arial,sans-serif;font-size:13px">n3=n3+1</div><div style="font-family:arial,sans-serif;font-size:13px">else</div>
<div style="font-family:arial,sans-serif;font-size:13px">if(c>0&d<=0)</div><div style="font-family:arial,sans-serif;font-size:13px">n4=n4+1</div><div style="font-family:arial,sans-serif;font-size:13px">endif</div>
<div style="font-family:arial,sans-serif;font-size:13px">endif</div><div style="font-family:arial,sans-serif;font-size:13px">endif</div><div style="font-family:arial,sans-serif;font-size:13px">endif</div><div style="font-family:arial,sans-serif;font-size:13px">
i=i+1</div><div style="font-family:arial,sans-serif;font-size:13px">endwhile</div><div style="font-family:arial,sans-serif;font-size:13px">j=j+1</div><div style="font-family:arial,sans-serif;font-size:13px">endwhile</div>
<div style="font-family:arial,sans-serif;font-size:13px">say n1</div><div style="font-family:arial,sans-serif;font-size:13px">say n2</div><div style="font-family:arial,sans-serif;font-size:13px">say n3</div><div style="font-family:arial,sans-serif;font-size:13px">
say n4</div></div><div style="font-family:arial,sans-serif;font-size:13px">'close 1'</div></div></div>
<br>_______________________________________________<br>
gradsusr mailing list<br>
<a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>
<a href="http://gradsusr.org/mailman/listinfo/gradsusr" target="_blank">http://gradsusr.org/mailman/listinfo/gradsusr</a><br>
<br></blockquote></div><br><br clear="all"><br>-- <br>Jeff Duda<br>Graduate research assistant<br>University of Oklahoma School of Meteorology<br>Center for Analysis and Prediction of Storms<br>
</div>