tag:blogger.com,1999:blog-63487100026390485272018-09-17T21:40:28.821+12:00GUM StewPosts about the GUM Tree Calculator (GTC) - a programmable calculator that handles measurement uncertainty.Il Cuoconoreply@blogger.comBlogger7125tag:blogger.com,1999:blog-6348710002639048527.post-53713122669511645872013-03-31T11:14:00.001+13:002013-03-31T11:34:03.265+13:00RF power measurement<h2>An example of how uncertainties combine</h2>This post discusses one aspect of an uncertainty calculation for a radio frequency (RF) power measurement. It builds on the earlier post <a href="http://gumstew.blogspot.co.nz/2013/03/formulating-problems-3.html" target="_blank">Formulating uncertainty (3)</a> in which a simple voltmeter model was described.<br /><br />We consider a measurement setup that is often used in calibration laboratories. A signal is applied to the input to a power splitter (see <a href="http://gumstew.blogspot.co.nz/2013/03/lazy-about-differentiation.html" target="_blank">lazy about differentiation</a>) and a pair of power sensors are connected to the two output ports.<br /><br /><div class="separator" style="clear: both; text-align: center;"><a href="http://2.bp.blogspot.com/-Pl5JLe4X5C8/UVdMcxvBMbI/AAAAAAAAAFU/PQ6ah_PHICI/s1600/Setup.png" imageanchor="1" style="margin-left: 1em; margin-right: 1em;"><img border="0" height="152" src="http://2.bp.blogspot.com/-Pl5JLe4X5C8/UVdMcxvBMbI/AAAAAAAAAFU/PQ6ah_PHICI/s320/Setup.png" width="320" /></a></div><br />The sensors convert the RF signal into a voltage that is measured with a voltmeter (not shown).<br /><br />To measure the power at one port, two voltage measurements are required: one with the power switched on, the other with the power off. The port power is then<br />\[<br />P_\mathrm{sens} = \frac{(V_\mathrm{off}-V_\mathrm{on})(V_\mathrm{off}+V_\mathrm{on})} {R_\mathrm{std}}<br />\]<br />where \(R_\mathrm{std}\) is the resistance of a precision standard resistor inside the power sensor circuit. \( P_\mathrm{sens} \) is the power registered by the sensor.<sup><a href="http://www.blogger.com/blogger.g?blogID=6348710002639048527#fn1" id="ref1">1</a></sup><br /><br />However, the ratio of port powers has to be calculated to benefit from using a splitter. So we are interested in<br />\[<br />R_\mathrm{sens} = \frac{P_\mathrm{sens \cdot 1}}{P_\mathrm{sens \cdot 2}} \;.<br />\]<br />Assuming different sensor electronics (and hence different resistors) at each port <br />\[<br /> R_\mathrm{sens} = \frac{R_\mathrm{std\cdot 2}}{R_\mathrm{std\cdot 1}}<br />\frac{(V_\mathrm{off\cdot 1}-V_\mathrm{on\cdot 1})(V_\mathrm{off\cdot 1}+V_\mathrm{on\cdot 1})}{(V_\mathrm{off\cdot 2}-V_\mathrm{on\cdot 2})(V_\mathrm{off\cdot 2}+V_\mathrm{on\cdot 2})}<br />\;.<br />\]<br />So, what is the uncertainty in the power ratio?<br /><br /><h3>Reusing a model</h3>The voltmeter model developed earlier will serve us here too. But a class for power sensors will by keeping an uncertain number for the precision resistor and defining the power calculation.<br /><pre class="brush: python;">class Sensor(object):<br /> def __init__(self,R,u_R,tag=None):<br /> """Initialise a new sensor object"""<br /> label = "R"<br /> if tag is not None: label += "_%s" % tag<br /> self.R = ureal(R,u_R,label=label)<br /><br /> def power(self,v_off,v_on):<br /> """Return an estimate of power"""<br /> return (v_off-v_on)*(v_off+v_on)/self.R<br /></pre>Using this class and the <span style="font-family: "Courier New",Courier,monospace;">Meter</span> class defined <a href="http://gumstew.blogspot.co.nz/2013/03/formulating-problems-3.html#MeterClass">here</a>, we can proceed as follows.<br /><br />First define objects for a voltmeter and two sensors <a href="http://www.blogger.com/blogger.g?blogID=6348710002639048527" name="anchor"></a><br /><pre class="brush: python;">m = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7,tag='m1')<br />s1 = Sensor(R=200,u_R=2E-4,tag='1')<br />s2 = Sensor(R=200,u_R=2E-4,tag='2')<br /></pre>Then convert raw voltage readings into uncertain numbers, obtain power estimates at each port and finally calculate the power ratio<br /><pre class="brush: python;">V1, V2, = m.voltage(2.6,tag='1'), m.voltage(2.4,tag='2')<br />V3, V4, = m.voltage(2.55,tag='3'), m.voltage(2.41,tag='4')<br /><br />P1 = s1.power(V1,V2)<br />P2 = s2.power(V3,V4)<br /><br />X = P1/P2<br /><br /></pre><h3>The results</h3>First, the uncertainty budget of the power ratio<br /><pre class="brush: python;">print summary(X)<br />for cpt in rp.budget(X,trim=0):<br /> print "%s: %G" % cpt<br /></pre>The output is <br /><pre class="brush: plain text;">1.4400922, u=4.9E-06, df=inf<br />e_ran_3: 2.69706E-06<br />e_ran_4: 2.40904E-06<br />e_ran_1: 1.947E-06<br />e_ran_2: 1.65899E-06<br />R_1: 1.44009E-06<br />R_2: 1.44009E-06<br />e_zero_m1: 4.64546E-09<br />e_gain_m1: 5.0822E-21<br /></pre>Random errors associated with the voltage measurements dominate this estimate.<br /><br />A different picture emerges if we look at just one of the power estimates<br /><pre class="brush: python;">print summary(P1)<br />for cpt in rp.budget(P1,trim=0):<br /> print "%s: %G" % cpt<br /></pre>The output is dominated by uncertainty in the residual meter gain error<br /><pre class="brush: plain text;">0.005000000, u=3.2E-08, df=inf<br />e_gain_m1: 3E-08<br />e_ran_1: 6.76E-09<br />e_ran_2: 5.76E-09<br />R_1: 5E-09<br />e_zero_m1: 2E-09<br /></pre>If instead we look at a single voltage reading<br /><pre class="brush: python;">print summary(V1)<br />for cpt in rp.budget(V1,trim=0):<br /> print "%s: %G" % cpt<br /></pre>We see that the random errors were dominated by both the residual meter gain and offset uncertainties.<br /><pre class="brush: plain text;">2.6000000, u=7.9E-06, df=inf<br />e_gain_m1: 7.8E-06<br />e_zero_m1: 1E-06<br />e_ran_1: 2.6E-07<br /></pre><br /><h3>In conclusion </h3>One of the advantages of object-oriented programming is that code written in one context can often be reused in another. GTC can capture the benefits from this. In this post we have reused the <span style="font-family: "Courier New",Courier,monospace;">Meter</span> class definition. <br /><br />This post has also illustrated how to combine classes of objects in data processing: a <span style="font-family: "Courier New",Courier,monospace;">Meter</span> object was used <a href="http://www.blogger.com/blogger.g?blogID=6348710002639048527#anchor">here</a> to create uncertain numbers for the four voltage estimates and these uncertain numbers were then directly used to obtain estimates of power by the two <span style="font-family: inherit;"><span style="font-family: "Courier New",Courier,monospace;">Sensor</span> </span>objects. The use of classes makes data processing code easier to read, understand and maintain.<br /><br />Finally, the way in which different sources of uncertainty propagate into the power ratio is not simple. A full analysis using calculus would be quite a lot more work than the GTC calculation here. With GTC, if we want to, we can track the relative importance of influence quantities through the calculation, from one intermediate step to the next. However, there is no need to audit all these details. GTC ensures that uncertainties are correctly propagated.<br /><br /><br /><sup id="fn1">1. There are other factors that need to be considered to estimate the signal generator power. In this post, we are only interested in the power level registered by the sensor<a href="http://www.blogger.com/blogger.g?blogID=6348710002639048527#ref1" title="Jump back to footnote 1 in the text.">↩</a></sup>Il Cuoconoreply@blogger.com0tag:blogger.com,1999:blog-6348710002639048527.post-60373428442590593452013-03-16T19:20:00.003+13:002013-03-22T11:12:06.402+13:00Lazy about differentiation<h2>Calculating derivatives</h2>While preparing a workshop about measurement uncertainty the other day, I decided to look for discussion topics in an Application Note about an instrument designed for high frequency measurements (at radio and microwave frequencies) . <br /><br />I found a section where the authors made a reasonable sounding decision, choosing between two alternative circuit configurations, one called a power splitter the other called a power divider.<br /><br />I decided to check their assumptions, which eventually led me to calculate some partial derivatives with GTC.<br /><h4> </h4><h4>The splitter and divider networks</h4>The power splitter has two resistors, each of impedance \( Z_0 \).<br /><br /><div class="separator" style="clear: both; text-align: center;"><a href="http://4.bp.blogspot.com/-lPJVU0TkHRw/UUP_ele6mwI/AAAAAAAAAE8/OjaStHCwkiE/s1600/Power+splitter.png" imageanchor="1" style="margin-left: 1em; margin-right: 1em;"><img border="0" src="http://4.bp.blogspot.com/-lPJVU0TkHRw/UUP_ele6mwI/AAAAAAAAAE8/OjaStHCwkiE/s1600/Power+splitter.png" /></a></div><br />As this diagram shows, there are three "ports" (pairs of lines). The input port is on the left and the two output ports are on the right.<br /><br />The power divider has three resistors, each of impedance \( \frac{Z_0}{3} \).<br /><br /><div class="separator" style="clear: both; text-align: center;"><a href="http://1.bp.blogspot.com/-6DGKfNIYriA/UUQARc-g9-I/AAAAAAAAAFE/o4gNn2oXTjU/s1600/Power+divider.png" imageanchor="1" style="margin-left: 1em; margin-right: 1em;"><img border="0" src="http://1.bp.blogspot.com/-6DGKfNIYriA/UUQARc-g9-I/AAAAAAAAAFE/o4gNn2oXTjU/s1600/Power+divider.png" /></a></div><br /><br />There are also three ports. The input is shown here on the left and the two outputs are on the right. <br /><br /><h4>Ideal behaviour and sensitivity to errors</h4>Now, when these networks are connected to an "ideal" circuit, the two lines of each output port are effectively connected together by an impedance of \( Z_0 \). In that case, the impedance at the input port can be calculated and is equal \( Z_0 \) too.<br /><br />That is an ideal configuration, and both the power splitter and power divider behave the same way.<br /><br />However, the Application Note claimed that if the terminating impedance at an output port was slightly different from \( Z_0 \), the resulting change of impedance seen at the input port would be smaller for a splitter than for a divider.<br /><br />In other words, the power splitter would be less sensitive to variation of the impedance at the output ports. <br /><br />This seems reasonable, because the impedances in the splitter are three times bigger of those in the divider. However, I decided to check.<br /><br /><h4>Checking the assumption</h4>We need equations for the input impedance \( Z_\mathrm{in} \), which will depend on \( Z_0 \) and the non-ideal impedance at one port \( Z \).<br /><br />For the power splitter, the equation is<sup><a href="https://www.blogger.com/blogger.g?blogID=6348710002639048527#fn1" id="ref1">1</a></sup> <br />\[<br />Z_\mathrm{in} = \frac{2 Z_0 (Z_0 + Z)}{3Z_0 + Z}<br />\]<br />and for the power divider, it is<br />\[<br />Z_\mathrm{in} = \frac{Z_0}{3} + \frac{\frac{4 Z_0}{3}(\frac{Z_0}{3} + Z)}{\frac{4 Z_0}{3} + \frac{Z_0}{3} + Z}<br />\] <br />To find the sensitivity of \( Z_\mathrm{in} \) to changes, these equations must be differentiated with respect to \( Z \).<br /><br />GTC can do these calculations quite easily (I'm a bit out of practice with differentiation).<br /><br />The following code differentiates the equations above. It first prints out the value of \( Z_\mathrm{in} \), which we expect to be equal to \( Z_0 = 50 \) in this case. Then the function <span style="font-family: "Courier New",Courier,monospace;">rp.sensitivity()</span> obtains the first partial derivative of \( Z_\mathrm{in} \) with respect to \( Z \).<br /><br />It is worth noting, that when <span style="font-family: "Courier New",Courier,monospace;">Z</span> is defined as an uncertain number, in line 14, the uncertainty argument is set to unity, although any number could have been used except for zero. A value of zero would have fooled GTC into thinking that <span style="font-family: "Courier New",Courier,monospace;">Z</span> was a constant, so <span style="font-family: "Courier New",Courier,monospace;">rp.sensitivity()</span> would have returned a value of zero.<br /><br /><pre class="brush: python; gutter: True">"""<br />Compare the sensitivity of the input impedance to departures from Z0.<br /><br />The Application Note, claims that the 2-resistor splitter is less <br />sensitive at the input port to changes in the output port impedance. <br />It isn't!!!<br /><br />"""<br />Z0 = 50 <br />Z0_3 = Z0/3.0<br /><br /># This is the impedance of the output port. <br /># The second parameter can be anything other than zero. <br />Z = ureal(Z0,1)<br /><br /># Splitter expression for input impedance<br />Z_in = 2*Z0*(Z + Z0)/(Z + 3*Z0)<br /><br /># <br />print value(Z_in)<br />print "sensitivity of Z_in to Z = %G" % rp.sensitivity(Z_in,Z)<br />print<br /><br /># Divider expression for input impedance<br />num = 4*Z0_3 * (Z0_3 + Z)<br />den = 4*Z0_3 + (Z0_3 + Z)<br />Z_in = Z0_3 + num/den<br /><br />print value(Z_in)<br />print "sensitivity of Z_in to Z = %G" % rp.sensitivity(Z_in,Z)</pre><h4 class="brush: python; gutter: True">And the results are ... </h4><pre class="brush: plain text">Input impedance Z0 = 50<br />sensitivity of Z_in to Z = 0.25<br /><br />Input impedance Z0 = 50<br />sensitivity of Z_in to Z = 0.25<br /></pre>In both cases, the value of derivative of \( Z_\mathrm{in} \) with respect to \( Z \) is \( \frac{1}{4} \). <br /><br />So the Application Note was wrong!<br /><br />It does not matter whether you use a divider or a splitter, the effect of an imperfect termination at the output ports is the same. <br /><br />I felt obliged to check this result by hand before posting, but it was certainly a lot harder to remember the calculus rules than to write a little GTC script like the one above.<br /><br /><hr /><sup id="fn1">[1] These equations can be obtained fairly easily if you know the rules for combining series and parallel impedances. Just remember to connect another \( Z_0 \) impedance across the output port terminals.<a href="https://www.blogger.com/blogger.g?blogID=6348710002639048527#ref1" title="Jump back to footnote 1 in the text.">↩</a></sup>Il Cuoconoreply@blogger.com0tag:blogger.com,1999:blog-6348710002639048527.post-52087048484623704962013-03-04T15:14:00.005+13:002013-03-31T10:12:08.615+13:00Formulating problems (3)<h2>Object-oriented programming and systematic errors </h2>In the previous post <a href="http://gumstew.blogspot.com/2013/02/formulating-problems-2.html" target="_blank">(March 3, 2013)</a>, I showed how to use GTC when both random and systematic errors contribute to the overall measurement uncertainty.<br /><br />The key is that an estimate of each measurement error in the model must be associated with one, and only one, uncertain number. <br /><br />In the last post, we saw that the idea of matching one uncertain number to one error estimate offers a straightforward way of using GTC. However, it becomes hard to keep track of all the different uncertain numbers as problems become more complicated.<br /><br />In this post, I will show how to define a class of objects to represent the voltmeter used earlier. This class of objects will encapsulate the data and data processing needed by the voltmeter error model, making data processing much more straightforward.<br /><br /><h3>A simple voltmeter</h3>An equation for the meter reading \(x_\mathrm{dis} \) in terms of the applied voltage \( V \) and instrument errors is<br />\[<br />x_\mathrm{dis} = (1 + e_\mathrm{gain} + e_\mathrm{ran}) V + e_\mathrm{zero}<br />\]<br />where, \( e_\mathrm{gain} \) is the residual error in the gain setting, \( e_\mathrm{zero} \) is the residual error in the zero offset and \( e_\mathrm{ran}\) is an error due to system noise.<br /><br />To carry out data processing, we invert this equation and replace the errors by their estimates (indicated by the hat symbol). In this way we can estimate the applied voltage given a reading<br />\[<br />\widehat{V} = \frac{x_\mathrm{dis} - \widehat{e}_\mathrm{zero}}{1 + \widehat{e}_\mathrm{gain} + \widehat{e}_\mathrm{ran}}<br />\]<br />Remember, \( e_\mathrm{gain} \) and \( e_\mathrm{zero} \) are systematic errors: they endure from one measurement to the next. Whereas \( e_\mathrm{ran} \) changes at every meter reading. So we usually estimate \( e_\mathrm{gain} \) and \( e_\mathrm{zero} \) just once, but need a fresh estimate of \( e_\mathrm{ran} \) for every reading.<br /><br />I'll now define a Python class to implement the required uncertain-number data processing.<br /><br /><h3>A voltmeter class</h3>In object-oriented languages like Python, the structure and behaviour of new types of objects can be defined by a 'class'. An object of a given class can have data (that other objects of the same class do not share) and functions that implement certain behaviour (shared among all objects of the same class).<br /><br />Here is a class definition for the voltmeter<br /><a name="MeterClass"></a><br /><pre class="brush: python;">class Meter(object):<br /> def __init__(self,u_e_gain,u_e_zero,u_e_ran=0,tag=None): <br /> """Initialise a new voltmeter object"""<br /><br /> self.u_e_ran = u_e_ran<br /> <br /> tag = "_%s" % tag if tag is not None else ""<br /> label = 'e_gain%s' % tag<br /> self.e_gain = ureal(0,u_e_gain,label=label)<br /> <br /> label = 'e_zero%s' % tag<br /> self.e_zero = ureal(0,u_e_zero,label=label)<br /><br /> self.u_e_ran = u_e_ran<br /><br /> def voltage(self,x,tag=None):<br /> """Return an uncertain number for the voltage"""<br /><br /> tag = "_%s" % tag if tag is not None else ""<br /> if self.u_e_ran != 0:<br /> label = 'e_ran%s' % tag<br /> e_ran = ureal(0,self.u_e_ran,label=label)<br /> else:<br /> e_ran = 0<br /><br /> return (x - self.e_zero)/(1 + self.e_gain + e_ran)<br /></pre>When objects of the Meter class are created, their data is initialised by the <span style="font-family: "Courier New",Courier,monospace;">__init__</span> function. We see here that uncertain numbers for the residual gain and zero setting errors are created during initialisation. The standard uncertainty associated with random noise is also saved.<br /><br />The <span style="font-family: "Courier New",Courier,monospace;">voltage</span> function uses the object-data to return an uncertain number for an applied voltage estimate. When <span style="font-family: "Courier New",Courier,monospace;">voltage</span> is called, an uncertain number associated with the random error is created first. This is combined with the uncertain numbers already defined to estimate the applied voltage. <br /><br /><h3>Resistance measurement with one meter</h3>As in the previous post, I will calculate an estimate of the unknown resistance \( R_\mathrm{x} \) from measurements of the voltage ratio \( V_\mathrm{x} / V_\mathrm{s} \) and an estimate of \( R_\mathrm{s} \) <br /><br /><div class="separator" style="clear: both; text-align: center;"><a href="http://3.bp.blogspot.com/-ZfOGQJUro54/USl4O_OG8JI/AAAAAAAAAEs/kdUJiZonhjU/s1600/resistance_ratio.png" imageanchor="1" style="margin-left: 1em; margin-right: 1em;"><img border="0" src="http://3.bp.blogspot.com/-ZfOGQJUro54/USl4O_OG8JI/AAAAAAAAAEs/kdUJiZonhjU/s1600/resistance_ratio.png" /></a></div><br /><br />To do this, I create an <span style="font-family: inherit;">object </span>of the Meter class by supplying numbers for the standard uncertainty parameters: <span style="font-family: "Courier New",Courier,monospace;">u_e_gain, u_e_zero <span style="font-family: inherit;">and </span>u_e_ran:</span> <br /><pre class="brush: python;">m = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7)</pre>Then use the meter object to obtain a pair of uncertain numbers for the voltage readings<br /><pre class="brush: python;">Vs = m.voltage(5.0100,tag='Vs')<br />Vx = m.voltage(4.9885,tag='Vx')<br /></pre><br />The code to calculate the unknown resistance and display the results (apart from the Meter class definition) now looks like this <br /><pre class="brush: python; gutter: True">m = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7)<br /><br />Rs = 1000 + ureal(0,1E-3,label='Rs')<br /><br />Vs = m.voltage(5.0100,tag='Vs')<br />Vx = m.voltage(4.9885,tag='Vx')<br /><br />Rx = Rs * Vx/Vs<br /><br />print "Rx = %G, u(Rx) = %G" % (value(Rx), uncertainty(Rx))<br />for cpt in reporting.budget(Rx,trim=0):<br /> print "%s: %G" % cpt<br /></pre>we obtain (as before)<br /><pre class="brush: plain text;">Rx = 995.709, u(Rx) = 0.00100562<br />Rs: 0.000995709<br />e_ran_Vs: 9.95709E-05<br />e_ran_Vx: 9.95709E-05<br />e_zero: 8.5657E-07<br />e_gain: 0<br /></pre>This is easier to read, write and understand than the code I posted earlier. <br /><br /><h3>Resistance measurement with two meters </h3>As might be expected, if two physically different meters are used to measure the different voltages, then two different Meter objects should be created to do the data processing.<br /><br />So with just a few simple changes to the code we have<br /><pre class="brush: python; gutter: True">m1 = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7,tag='m1')<br />m2 = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7,tag='m2')<br /><br />Vs = m1.voltage(5.0100,tag='Vs')<br />Vx = m2.voltage(4.9885,tag='Vx')<br /><br />Rs = 1000 + ureal(0,1E-3,label='Rs')<br /><br />Rx = Rs * Vx/Vs<br />print "Rx = %G, u(Rx) = %G" % (value(Rx), uncertainty(Rx))<br />for cpt in reporting.budget(Rx,trim=0):<br /> print "%s: %G" % cpt<br /></pre>which gives the same results as we obtained last time<br /><pre class="brush: plain text;">Rx = 995.709, u(Rx) = 0.0043516<br />e_gain_m1: 0.00298713<br />e_gain_m2: 0.00298713<br />Rs: 0.000995709<br />e_zero_m2: 0.000199601<br />e_zero_m1: 0.000198744<br />e_ran_Vs: 9.95709E-05<br />e_ran_Vx: 9.95709E-05<br /></pre><br />This simple and intuitive behaviour is due to the way that uncertain numbers are created during object initialisation. When one Meter object is initialsed, a pair of uncertain numbers are created; but when two objects are initialsed, each has its own pair of uncertain numbers representing independent estimates of the different gain and zero setting errors in each meter.<br /><br />Hence the definition of a Meter class is able to capture the creation of uncertain numbers we need in a natural and intuitive way that matches the model of the measurement system.Il Cuoconoreply@blogger.com0tag:blogger.com,1999:blog-6348710002639048527.post-51263142325482024512013-03-03T15:34:00.002+13:002013-03-17T10:56:29.492+13:00Formulating problems (2)<h2>Systematic errors </h2>In the previous post <a href="http://gumstew.blogspot.com/2013/02/formulating-problems-i.html" target="_blank">(Feb 10, 2013)</a>, I showed how GTC lets you build up an uncertainty calculation step by step. In this post, I will explain how to deal with systematic errors.<br /><br /><h3>A simple voltmeter</h3>Last time, I used an equation for the reading displayed by the meter \(x_\mathrm{dis} \) in terms of the actual applied voltage \( V \) and instrumental errors<br />\[<br />x_\mathrm{dis} = (1 + e_\mathrm{gain}) V + e_\mathrm{zero}<br />\]<br />where, \( e_\mathrm{gain} \) is the residual error in the gain setting and \( e_\mathrm{zero} \) is the residual error in the zero offset. I did not include a term for random noise, because that was estimated from a sample of readings. In this post we will add a term for random noise to the model instead. So the model becomes<br />\[<br />x_\mathrm{dis} = (1 + e_\mathrm{gain} + e_\mathrm{ran}) V + e_\mathrm{zero}<br />\]<br />Now, \( e_\mathrm{gain} \) and \( e_\mathrm{zero} \) are systematic errors: they endure from one measurement to the next. So, for example, the contribution to uncertainty from random errors can be reduced by averaging readings from the same meter, but not the uncertainty contributions associated with \( e_\mathrm{gain} \) and \( e_\mathrm{zero} \).<br /><br />How can this be handled with GTC? <br /><br />Simple. <br /><br />We associate one, and only one, uncertain number with the estimate of any quantity in the meter model. So, there will be one uncertain number for the estimate of \( e_\mathrm{gain} \) and one for the estimate of \( e_\mathrm{zero} \), regardless of how many readings are made with the voltmeter. <br /><br />There will be a different uncertain number for each reading representing the random noise \( e_{\mathrm{ran}\cdot i}\), because the error changes from one reading to the next.<br /><br /><br /><h3>A voltage ratio measurement</h3>We will consider a simple measurement in which a resistance value \( R_\mathrm{x} \) is estimated by measuring a voltage ratio. <br /><br />The circuit is shown in this figure,<br /><br /><div class="separator" style="clear: both; text-align: center;"><a href="http://3.bp.blogspot.com/-ZfOGQJUro54/USl4O_OG8JI/AAAAAAAAAEs/kdUJiZonhjU/s1600/resistance_ratio.png" imageanchor="1" style="margin-left: 1em; margin-right: 1em;"><img border="0" src="http://3.bp.blogspot.com/-ZfOGQJUro54/USl4O_OG8JI/AAAAAAAAAEs/kdUJiZonhjU/s1600/resistance_ratio.png" /></a></div><br />where a single meter is used to measure the potential difference across each resistor. We will just assume that the current through the circuit does not vary.<br /><br />With \( R_\mathrm{s} \) known, \( R_\mathrm{x} \) can be calculated <br />\[<br />R_\mathrm{x} = R_\mathrm{s}\frac{V_\mathrm{x}}{V_\mathrm{s}}<br />\]<br />Suppose that we measure \( V_\mathrm{x} = 4.9885 \, \mathrm{V}\) and \( V_\mathrm{s} =5.0100 \, \mathrm{V}\). Suppose also that \( R_\mathrm{s} = 1000.000 \; u(R_\mathrm{s}) = 1 \times 10^{-3} \, \Omega\) and that the relative standard uncertainty in the meter gain is \( u(e_\mathrm{gain}) = 3 \times 10^{-6} \; \mathrm{V}/ \mathrm{V} \), that the relative standard uncertainty due to meter noise is \( u(e_\mathrm{ran}) = 1 \times 10^{-7} \; \mathrm{V}/ \mathrm{V} \) and the standard uncertainty in the zero setting is \( u(e_\mathrm{zero}) = 1 \times 10^{-6} \; \mathrm{V}\).<br /><br />Common sense tells us that the uncertainty due to gain error cancels when the voltage ration is calculated. So, one might simply ignore that term. The zero setting uncertainty will not have very much effect on the combined uncertainty \( R_\mathrm{s}\) either.<br /><br />However, I am not going to make those assumptions. Instead, I will show that GTC arrives at the correct result anyway, with no need to think about simplifying the mathematical model. <br /><br />I will define uncertain numbers for the estimates of \( R_\mathrm{s} \), \( V_\mathrm{s} \) and \( V_\mathrm{x} \) and then work out<br />\[<br />R_\mathrm{x} = R_\mathrm{s} \frac{V_\mathrm{x}}{V_\mathrm{s}}<br />\]<br />Here is the standard resistance<br /><pre class="brush: python;">R_s = 1000 + ureal(0,1E-3,label='R_s')<br /></pre>and here are the voltmeter errors<br /><pre class="brush: python;">e_zero = ureal(0,1E-6,label='e_zero')<br />e_gain = ureal(1,3E-6,label='e_gain')<br /></pre>Uncertain numbers representing the individual voltage measurements are then (note an uncertain number for the meter noise is created independently for each reading)<br /><pre class="brush: python;">e_ran_1 = ureal(1,1E-7,label='e_ran_1')<br />V_s = (5.0100 - e_zero) / (e_gain * e_ran_1) <br /><br />e_ran_2 = ureal(1,1E-7,label='e_ran_2')<br />V_x = (4.9885 - e_zero) / (e_gain * e_ran_2) <br /></pre>And the unknown resistance is<br /><pre class="brush: python;">R_x = R_s * V_x/V_s<br /></pre>The value, the uncertainty and the uncertainty budget of the result are displayed by <br /><pre class="brush: python;">print "R_x = %G, u(R_x) = %G" % (value(R_x), uncertainty(R_x))<br />for cpt in reporting.budget(R_x,trim=0):<br /> print "%s: %G" % cpt<br /></pre>We obtain the following results:<br /><pre class="brush: plain text;">R_x = 995.709, u(R_x) = 0.00100562<br />R_s: 0.000995709<br />e_ran_1: 9.95709E-05<br />e_ran_2: 9.95709E-05<br />e_zero: 8.5657E-07<br />e_gain: 0<br /></pre>As expected, the uncertainty due to an error in the voltmeter gain setting is zero and the uncertainty due to an error in the zero offset is also small. The uncertainty of this measurement is dominated by the standard resistor uncertainty, followed in importance by the meter noise.<br /><br />Now, suppose that two meters of the same type are used to simultaneously measure \( V_\mathrm{s} \) and \( V_\mathrm{x} \). This would ensure that the current flowing in the circuit elements was indeed the same, rather than just assuming it to be constant. <br /><br />Different meters will have different residual gain and zero setting errors. So the calculation must now define uncertain numbers for each meter. <br /><br />The code now looks like this<br /><pre class="brush: python; gutter: true;">R_s = 1000 + ureal(0,1E-3,label='R_s')<br /><br /># Meter 1<br />e_zero_1 = ureal(0,1E-6,label='e_zero_1')<br />e_gain_1 = ureal(1,3E-6,label='e_gain_1')<br /><br /># Meter 2<br />e_zero_2 = ureal(0,1E-6,label='e_zero_2')<br />e_gain_2 = ureal(1,3E-6,label='e_gain_2')<br /><br /># Readings<br />e_ran_1 = ureal(1,1E-7,label='e_ran_1')<br />V_s = (5.0100 - e_zero_1) / (e_gain_1 * e_ran_1) <br /><br />e_ran_2 = ureal(1,1E-7,label='e_ran_2')<br />V_x = (4.9885 - e_zero_2) / (e_gain_2 * e_ran_2) <br /><br /># Result<br />R_x = R_s * V_x/V_s<br /><br /># Display results<br />print "R_x = %G, u(R_x) = %G" % (value(R_x), uncertainty(R_x))<br />for cpt in reporting.budget(R_x,trim=0):<br /> print "%s: %G" % cpt<br /></pre>We obtain the following results:<br /><pre class="brush: plain text;">R_x = 995.709, u(R_x) = 0.0043516<br />e_gain_1: 0.00298713<br />e_gain_2: 0.00298713<br />R_s: 0.000995709<br />e_zero_2: 0.000199601<br />e_zero_1: 0.000198744<br />e_ran_1: 9.95709E-05<br />e_ran_2: 9.95709E-05<br /></pre>This is quite a different result. The combined standard uncertainty \( u(R_\mathrm{x}) \) is a lot bigger and the gain error uncertainties now dominate the uncertainty budget. <br /><br /><h3>In summary</h3>Every error in a measurement model has to be estimated in the GUM process of calculating uncertainty. Often residual errors are approximately zero, or unity, but these estimates still carry uncertainty.<br /><br />Using GTC, every estimate of a quantity in the model is associated with one, and only one, uncertain number. When this rule is understood, GTC can be used to handle quite complicated data processing problems.<br /><br />The next post will revisit this example. I will show then how object-oriented programming techniques can be used to implement the voltmeter model. <br /><br />Il Cuoconoreply@blogger.com0tag:blogger.com,1999:blog-6348710002639048527.post-24403150396617912182013-02-10T17:16:00.002+13:002013-03-17T10:58:09.881+13:00Formulating problems (1) <h3>Developing an uncertainty calculation with GTC</h3>The traditional approach to uncertainty calculation described in the GUM requires a mathematical model of the measurement<br />\[<br />Y = f(X_1, X_2,\cdots) <br />\]<br />where \( Y \) is the measurand and \( X_1, X_2,\cdots \) are all the quantities that can influence the outcome of the measurement. Usually, it is difficult to write down this function directly.<br /><br />GTC also works with a model, but instead of writing down a single function that depends on everything, we can build models in smaller, convenient steps.<br /><br />An example should make this clearer.<br /><br />Consider a measurement of current through a standard resistor (also called current shunt) made using a voltmeter, as shown<br /><br /><div class="separator" style="clear: both; text-align: center;"><a href="http://3.bp.blogspot.com/-lZJ8F--M7LI/URHpYhTMKFI/AAAAAAAAAEc/gCPSdRUJN28/s1600/current_shunt.png" imageanchor="1" style="margin-left: 1em; margin-right: 1em;"><img border="0" src="http://3.bp.blogspot.com/-lZJ8F--M7LI/URHpYhTMKFI/AAAAAAAAAEc/gCPSdRUJN28/s1600/current_shunt.png" /></a></div><br />If the resistance \( R \) of the current shunt is known, \( I \) can be calculated using Ohm's law<br />\[<br />I = \frac{V}{R} <br />\]<br />Four voltage readings have been recorded (all in mV)<br />\[ <br />100.6512 \,, \; 100.6401 \,, \; 100.6420 \,, \; 100.6101 \,,<br />\]<br />and the resistor has been calibrated <br />\[<br />R = 0.009011 \, {\Omega} \quad \frac{u(R)}{R} = 0.0003<br />\]<br />So the current can be calculated using uncertain numbers for the voltage and resistance estimates<br /><pre class="brush: python;">x_obs = [100.6512E-3, 100.6401E-3, 100.6420E-3, 100.6101E-3]<br />V = type_a.estimate(x_obs,label='V_type_a')<br /><br />e_R_cal = ureal(1,0.0003,label='e_R_cal')<br />R = 0.009001 * e_R_cal<br /><br />I = V/R<br />I.label = 'I'<br />print summary(I)<br /></pre>The result (in amperes) is <br /><pre class="brush: plain text;">I: 11.1805, u=0.0035, df=465.6</pre><h3>Enhancing the model: temperature dependence</h3>Suppose now that, during measurements, the temperate was controlled only to within \( \pm 3^\circ \mathrm{C}\) of the reported calibration temperature. This means that the value of \( R \) may be slightly different from the one reported. To allow for this, an additional component of uncertainty in the resistance due to temperature should be added to the model.<br /><br />Given the current-shunt temperature coefficient \( \alpha = 4 \times 10^{-6}\; \Omega/ {^\circ \mathrm{C}}\), and the temperature range, we may modify the model of the measurement by simply writing<br /><pre class="brush: python;">e_R_T = ureal(1,4E-6*type_b.uniform(3),label='e_R_T')<br />R = 0.009001 * e_R_cal * e_R_T<br /></pre>Repeating the calculation <br /><pre class="brush: python;">I = V/R<br />I.label = 'I'<br />print summary(I)<br /></pre>we obtain the same result<br /><pre class="brush: plain text;">I: 11.1805, u=0.0035, df=465.6</pre>Apparently, the temperature does not have a big influence. To investigate further, we can look at the uncertainty budget for \( I \)<br /><pre class="brush: python;">for cpt in reporting.budget(I,trim=0):<br /> print "%s: %G" % cpt<br /></pre>This displays the magnitude of the component of uncertainty in \( I \) due to each of the three influence quantities<br /><pre class="brush: plain text;">e_R_cal: 0.00335416<br />V_type_a: 0.000990883<br />e_R_T: 7.74609E-05</pre><h3 class="brush: python;">Enhancing the model: voltmeter errors</h3><div class="brush: python;">Another source of measurement error will be the voltmeter accuracy. In other words, the voltage indicated will not be exactly the same as the applied voltage. Voltmeter errors can be represented as</div>\[<br />x_\mathrm{dis} = (1 + e_\mathrm{gain}) V + e_\mathrm{zero}<br />\] <br />where \( x_\mathrm{dis}\) is the reading shown when a voltage \( V \) is applied: \( e_\mathrm{zero} \) is an offset error (i.e., the reading when \(V=0\) ) and \( e_\mathrm{gain} \) is a gain (i.e., the reading is not exactly \( 1 \times V\) ).<br /><br />The voltmeter specifications give an voltage offset error uncertainty of \( 500 \,\mu V\) and a gain (scale) error uncertainty of \( 25 \times 10^{-6} \,V/V \). Each error contributes to the uncertainty in \( x_\mathrm{dis} \) as an estimate of the voltage \( V \).<br /><br />The measurement model can be extended by adding two more uncertain numbers <br /><pre class="brush: python;">e_gain = ureal(0,25E-6,label='e_V_gain')<br />e_off = ureal(0,500E-6,label='e_V_off')<br />V = (type_a.estimate(v_obs,label='V_type_a') - e_off ) /(1 + e_gain)<br /></pre>Note that these new uncertain numbers have been combined with the type-A estimate of voltage, obtained from the observations. This is OK, because the errors are systematic (they do not change during the four repeat voltage readings), so the effect they have on the average voltage is the same as on each individual reading (however, a future post will show another way of introducing handling this type of problem). <br /><br />Recalculating the current now<br /><pre class="brush: python;">I = V/R<br />I.label = 'I'<br />print summary(I)<br />for cpt in reporting.budget(I,trim=0):<br /> print "%s: %G" % cpt<br /></pre>we obtain<br /><pre class="brush: plain text;">I: 11.181, u=0.056, df=inf<br /><br />e_V_off: 0.0555494<br />e_R_cal: 0.00335416<br />V_type_a: 0.000990883<br />e_V_gain: 0.000279513<br />e_R_T: 7.74609E-05<br /></pre>In this case, the voltmeter offset uncertainty makes a significant contribution to the current measurement.<br /><br /><h3>Putting it all together</h3>The GTC code now looks like this:<br /><pre class="brush: python; gutter: true;">e_R_cal = ureal(1,0.0003,label='e_R_cal')<br />e_R_T = ureal(1,4E-6*type_b.uniform(3),label='e_R_T')<br /><br />R = 0.009001 * e_R_cal * e_R_T<br /><br />x_obs = [100.6512E-3, 100.6401E-3, 100.6420E-3, 100.6101E-3]<br /><br />e_gain = ureal(0,25E-6,label='e_V_gain')<br />e_off = ureal(0,500E-6,label='e_V_off')<br /><br />V = (type_a.estimate(x_obs,label='V_type_a') - e_off ) /(1 + e_gain)<br /><br />I = V/R<br />I.label = 'I'<br /><br />print summary(I)<br /><br />for cpt in reporting.budget(I,trim=0):<br /> print "%s: %G" % cpt<br /></pre>Try to imagine writing down the measurement model for this problem as a simple function \( f() \). It is not straight forward. <br /><br />Using GTC, on the other hand, builds the model implicitly, as the different contributions to measurement uncertainty are introduced. Also, as shown in this post, the model can be easily expanded, or simplified, by adding or removing uncertain numbers from the calculation. <br /><br />That makes GTC a very convenient tool for doing simple uncertainty calculations.Il Cuoconoreply@blogger.com0tag:blogger.com,1999:blog-6348710002639048527.post-57223368820845566672013-01-20T16:16:00.000+13:002013-03-17T10:59:49.228+13:00Reporting results<h2>Obtaining numerical results from uncertain numbers</h2>Every uncertain real number has the attributes:<br /><ul><li>value</li><li>uncertainty</li><li>degrees of freedom</li></ul><i>Value</i> is the estimate of a quantity of interest. <i>Uncertainty</i> is the standard uncertainty of this estimate and the <i>degrees of freedom</i> is, well, the degrees of freedom associated with the uncertainty.<br /><br />The numerical values of these attributes are used to calculate an <i>expanded uncertainty</i> for a quantity of interest. <br /><br />GTC functions can obtain attribute values from an uncertain number, or there is a set of Python object-attributes that do the same thing. For example, if<br /><pre class="brush: python;">y = ureal(1.1,2.5,6)<br /></pre>then<br /><pre class="brush: python;">print value(y)<br />print uncertainty(y)<br />print dof(y)<br /></pre>yields<br /><pre class="brush: plain text;">1.1<br />2.5<br />6<br /></pre>as does the more succinct form using object attributes<br /><pre class="brush: python;">print y.x<br />print y.u<br />print y.df<br /></pre><h3>Expanded uncertainty </h3>To calculate an expanded uncertainty, the standard uncertainty of a result is multiplied by a coverage factor \( k_p \) that depends on the number of degrees of freedom (\( p \) is the level of confidence, or coverage probability).<br /><br />For the uncertain number above (and for a 95% level of confidence)<br /><pre class="brush: python;">k_95 = reporting.k_factor(y.df,p=95)<br />U_95 = k_95 * y.u<br /></pre>Alternatively, the lower and upper bounds of an uncertainty interval can be calculated directly <br /><pre class="brush: python;">Y_lb, Y_ub = reporting.uncertainty_interval(y,p=95)<br /></pre>The numbers obtained here are -5.017 and 7.217<br /><br /><h3>Uncertainty budgets</h3>An uncertainty budget breaks down the various contributions to the combined uncertainty of a result.<br /><br />For example, the following three uncertain numbers are associated with measurements of voltage, current and phase in an electrical circuit <br /><pre class="brush: python;">V = ureal(4.999,3.2E-3,label='V') # volt<br />I = ureal(19.661E-3,9.5E-6,label='I') # amp<br />phi = ureal(1.04446,7.5E-4,label='phi') # radian<br /></pre>Based on these measurements, an uncertain number for the resistance \( R \) is<br /><pre class="brush: python;">R = V * cos(phi) / I<br /></pre>We obtain a budget for this by writing<br /><pre class="brush: python;">for cpt in reporting.budget(R):<br /> print "%s: %G" % cpt<br /></pre>which generates a list of names and components of uncertainty in the same units as the measurement result and in order of decreasing magnitude<br /><pre class="brush: plain text;">phi: 0.164885<br />V: 0.0817649<br />I: 0.0617189<br /></pre>This shows that phase is the most important source of measurement error. The contribution to uncertainty from a phase error is likely to be about twice as big as voltage and current errors.<br /><br /><h3>Correlation coefficients</h3>Measurement results may become correlated if they share influence quantities, or are calculated from common data sets. <br /><br />Continuing with the example above, an estimate of reactance is<br /><pre class="brush: python;">X = V * sin(phi) / I<br /></pre>Now, the estimates associated with <span style="font-family: "Courier New",Courier,monospace;">R</span> and <span style="font-family: "Courier New",Courier,monospace;">X</span> will not be independent, because they were calculated using the same set of data. Their correlation coefficient can be calculated by<br /><pre class="brush: python;">get_correlation(R,X)<br /></pre>(the value obtained is 0.058)<br /> <br /><br />Il Cuoconoreply@blogger.com0tag:blogger.com,1999:blog-6348710002639048527.post-74714355089354004462013-01-13T11:42:00.001+13:002013-02-06T17:10:59.728+13:00Defining uncertain numbers<h2>Representing estimates</h2>Measurements estimate quantities of interest, but there is always some uncertainty. So GTC has built-in types, called uncertain numbers, that represent estimates and keep track of their uncertainty at the same time.<br /><br />The first task then, when using GTC, is to define uncertain numbers.<br /><br /><h3>Absolute uncertainty </h3>When an estimate and standard uncertainty are given in the same units, the definition of an uncertain number is straight forward. <br /><br />Suppose, a length is reported <br />\[<br />\widehat{l} = 10.1\,\mathrm{cm},\; u(\widehat{l}) = 0.3\,\mathrm{cm}\;,<br />\]<br />an uncertain number for \( \widehat{l} \) may be defined by writing<br /><pre class="brush: python;">l = ureal(10.1,0.3) </pre>although it would be better to write<br /><pre class="brush: python;">l = 10.1 - ureal(0,0.3)<br /></pre>The slightly longer form is preferred because the uncertainty is associated with the unknown measurement error (estimated as zero), rather than the number 10.1 (which is known). This corresponds to the mathematical relation between a result \( \widehat{l} \) and the quantity that was intended to be measured <br />\[<br />l = \widehat{l} - e \;,<br />\]<br />with \( e \) being the (unknown) measurement error.<br /><br /><h3>Relative uncertainty</h3>Sometimes an uncertainty is expressed relative to the magnitude of an estimate (i.e., as a standard relative uncertainty). <br /><br />For instance, a length may be reported as<br />\[<br />\widehat{l} = 10.1\,\mathrm{cm},\; \frac{u(\widehat{l})}{\widehat{l}} = 3\% \;.<br />\]<br />An uncertain number for \( \widehat{l} \) in this case can be defined by writing<br /><pre class="brush: python;">l = 10.1 * ureal(1.0,0.03)<br /></pre>which corresponds to the relation <br />\[<br />l = \widehat{l} \cdot e_\mathrm{rel} \;,<br />\]<br />with \( e_\mathrm{rel} \) an unknown relative error (estimated as unity).<br /><h3><br />Expanded uncertainty</h3>Whether relative or absolute, uncertainties are more often reported "expanded", than as standard uncertainties. <br /><br />An expanded uncertainty is just a standard uncertainty, or a relative standard uncertainty, multiplied by a <i>coverage factor</i>. Formally, the expanded uncertainty is represented by a capital \( U \) and the coverage factor by \( k_p \), so for length we might write<br />\[<br />U(\widehat{l}) = k_p \cdot u(\widehat{l}) \; ,<br />\]<br />\( p \) is a coverage probability, or level of confidence, associated with the expanded uncertainty.<br /><br />Expanded uncertainties should <b>never</b> be used to define uncertain numbers. They should be converted to a standard uncertainty, by dividing by the coverage factor. The coverage factor should be used to find a number of degrees of freedom. The standard uncertainty and the degrees of freedom can then be used to define an uncertain number. <br /><br />For example, suppose a length has been reported as \( \widehat{l} = 10.1 \) with an expanded uncertainty of \( k_p = 2.1 \) at a level of confidence of \(p = 95\%\). Let's say we already have the following lines of GTC code<br /><pre class="brush: python;">U_l = 0.63<br />k_95 = 2.1<br /></pre>The standard uncertainty and degrees of freedom can be easily found (note, '<span style="font-family: "Courier New",Courier,monospace;">k_to_dof</span>' is a function in the GTC '<span style="font-family: "Courier New",Courier,monospace;">reporting</span>' module) <br /><pre class="brush: python;">u_l = U_l / k_95<br />df_l = reporting.k_to_dof( k_95, 95 )<br /></pre>and an uncertain number for \( \widehat{l} \) can be defined with these values<br /><pre class="brush: python;">l = 10.1 - ureal(0,u_l,df_l)<br /></pre><br /><br /><br /><br /><br />Il Cuoconoreply@blogger.com0