31 March 2013

RF power measurement

An example of how uncertainties combine

This post discusses one aspect of an uncertainty calculation for a radio frequency (RF) power measurement. It builds on the earlier post Formulating uncertainty (3) in which a simple voltmeter model was described.

We consider a measurement setup that is often used in calibration laboratories. A signal is applied to the input to a power splitter (see lazy about differentiation) and a pair of power sensors are connected to the two output ports.

The sensors convert the RF signal into a voltage that is measured with a voltmeter (not shown).

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
P_\mathrm{sens} = \frac{(V_\mathrm{off}-V_\mathrm{on})(V_\mathrm{off}+V_\mathrm{on})} {R_\mathrm{std}}
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.1

However, the ratio of port powers has to be calculated to benefit from using a splitter. So we are interested in
R_\mathrm{sens} = \frac{P_\mathrm{sens \cdot 1}}{P_\mathrm{sens \cdot 2}} \;.
Assuming different sensor electronics (and hence different resistors) at each port
 R_\mathrm{sens} = \frac{R_\mathrm{std\cdot 2}}{R_\mathrm{std\cdot 1}}
\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})}
So, what is the uncertainty in the power ratio?

Reusing a model

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.
class Sensor(object):
    def __init__(self,R,u_R,tag=None):
        """Initialise a new sensor object"""
        label = "R"
        if tag is not None: label += "_%s" % tag
        self.R = ureal(R,u_R,label=label)

    def power(self,v_off,v_on):
        """Return an estimate of power"""
        return (v_off-v_on)*(v_off+v_on)/self.R
Using this class and the Meter class defined here, we can proceed as follows.

First define objects for a voltmeter and two sensors
m = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7,tag='m1')
s1 = Sensor(R=200,u_R=2E-4,tag='1')
s2 = Sensor(R=200,u_R=2E-4,tag='2')
Then convert raw voltage readings into uncertain numbers, obtain power estimates at each port and finally calculate the power ratio
V1, V2, = m.voltage(2.6,tag='1'), m.voltage(2.4,tag='2')
V3, V4, = m.voltage(2.55,tag='3'), m.voltage(2.41,tag='4')

P1 = s1.power(V1,V2)
P2 = s2.power(V3,V4)

X = P1/P2

The results

First, the uncertainty budget of the power ratio
print summary(X)
for cpt in rp.budget(X,trim=0):
    print "%s: %G" % cpt
The output is
1.4400922, u=4.9E-06, df=inf
e_ran_3: 2.69706E-06
e_ran_4: 2.40904E-06
e_ran_1: 1.947E-06
e_ran_2: 1.65899E-06
R_1: 1.44009E-06
R_2: 1.44009E-06
e_zero_m1: 4.64546E-09
e_gain_m1: 5.0822E-21
Random errors associated with the voltage measurements dominate this estimate.

A different picture emerges if we look at just one of the power estimates
print summary(P1)
for cpt in rp.budget(P1,trim=0):
    print "%s: %G" % cpt
The output is dominated by uncertainty in the residual meter gain error
0.005000000, u=3.2E-08, df=inf
e_gain_m1: 3E-08
e_ran_1: 6.76E-09
e_ran_2: 5.76E-09
R_1: 5E-09
e_zero_m1: 2E-09
If instead we look at a single voltage reading
print summary(V1)
for cpt in rp.budget(V1,trim=0):
    print "%s: %G" % cpt
We see that the random errors were dominated by both the residual meter gain and offset uncertainties.
2.6000000, u=7.9E-06, df=inf
e_gain_m1: 7.8E-06
e_zero_m1: 1E-06
e_ran_1: 2.6E-07

In conclusion

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 Meter class definition.

This post has also illustrated how to combine classes of objects in data processing: a Meter object was used here 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 Sensor objects. The use of classes makes data processing code easier to read, understand and maintain.

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.

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

16 March 2013

Lazy about differentiation

Calculating derivatives

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) .

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.

I decided to check their assumptions, which eventually led me to calculate some partial derivatives with GTC.


The splitter and divider networks

The power splitter has two resistors, each of impedance \( Z_0 \).

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.

The power divider has three resistors, each of impedance \( \frac{Z_0}{3} \).

There are also three ports. The input is shown here on the left and the two outputs are on the right.

Ideal behaviour and sensitivity to errors

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.

That is an ideal configuration, and both the power splitter and power divider behave the same way.

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.

In other words, the power splitter would be less sensitive to variation of the impedance at the output ports.

This seems reasonable, because the impedances in the splitter are three times bigger of those in the divider. However, I decided to check.

Checking the assumption

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 \).

For the power splitter, the equation is1
Z_\mathrm{in} = \frac{2 Z_0 (Z_0 + Z)}{3Z_0 + Z}
and for the power divider, it is
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}
To find the sensitivity of \( Z_\mathrm{in} \) to changes, these equations must be differentiated with respect to \( Z \).

GTC can do these calculations quite easily (I'm a bit out of practice with differentiation).

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 rp.sensitivity() obtains the first partial derivative of \( Z_\mathrm{in} \) with respect to \( Z \).

It is worth noting, that when Z 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 Z was a constant, so rp.sensitivity() would have returned a value of zero.

Compare the sensitivity of the input impedance to departures from Z0.

The Application Note, claims that the 2-resistor splitter is less  
sensitive at the input port to changes in the output port impedance. 
It isn't!!!

Z0 = 50         
Z0_3 = Z0/3.0

# This is the impedance of the output port. 
# The second parameter can be anything other than zero. 
Z = ureal(Z0,1)

# Splitter expression for input impedance
Z_in = 2*Z0*(Z + Z0)/(Z + 3*Z0)

print value(Z_in)
print "sensitivity of Z_in to Z = %G" % rp.sensitivity(Z_in,Z)

# Divider expression for input impedance
num = 4*Z0_3 * (Z0_3 + Z)
den = 4*Z0_3 + (Z0_3 + Z)
Z_in = Z0_3 + num/den

print value(Z_in)
print "sensitivity of Z_in to Z = %G" % rp.sensitivity(Z_in,Z)

And the results are ... 

Input impedance Z0 = 50
sensitivity of Z_in to Z = 0.25

Input impedance Z0 = 50
sensitivity of Z_in to Z = 0.25
In both cases, the value of derivative of \( Z_\mathrm{in} \) with respect to \( Z \) is \( \frac{1}{4} \). 

So the Application Note was wrong!

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.

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.

[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.

4 March 2013

Formulating problems (3)

Object-oriented programming and systematic errors

In the previous post (March 3, 2013), I showed how to use GTC when both random and systematic errors contribute to the overall measurement uncertainty.

The key is that an estimate of each measurement error in the model must be associated with one, and only one, uncertain number.

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.

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.

A simple voltmeter

An equation for the meter reading \(x_\mathrm{dis} \) in terms of the applied voltage \( V \) and instrument errors is
x_\mathrm{dis}  = (1 + e_\mathrm{gain} + e_\mathrm{ran}) V + e_\mathrm{zero}
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.

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
\widehat{V} = \frac{x_\mathrm{dis} - \widehat{e}_\mathrm{zero}}{1 + \widehat{e}_\mathrm{gain} + \widehat{e}_\mathrm{ran}}
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.

I'll now define a Python class to implement the required uncertain-number data processing.

A voltmeter class

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).

Here is a class definition for the voltmeter

class Meter(object):
    def __init__(self,u_e_gain,u_e_zero,u_e_ran=0,tag=None): 
        """Initialise a new voltmeter object"""

        self.u_e_ran = u_e_ran
        tag = "_%s" % tag if tag is not None else ""
        label = 'e_gain%s' % tag
        self.e_gain = ureal(0,u_e_gain,label=label)
        label = 'e_zero%s' % tag
        self.e_zero = ureal(0,u_e_zero,label=label)

        self.u_e_ran = u_e_ran

    def voltage(self,x,tag=None):
        """Return an uncertain number for the voltage"""

        tag = "_%s" % tag if tag is not None else ""
        if self.u_e_ran != 0:
            label = 'e_ran%s' % tag
            e_ran = ureal(0,self.u_e_ran,label=label)
            e_ran = 0

        return (x - self.e_zero)/(1 + self.e_gain + e_ran)
When objects of the Meter class are created, their data is initialised by the __init__ 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.

The voltage function uses the object-data to return an uncertain number for an applied voltage estimate. When voltage 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.

Resistance measurement with one meter

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} \)

To do this, I create an object of the Meter class by supplying numbers for the standard uncertainty parameters: u_e_gain, u_e_zero and u_e_ran:
m = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7)
Then use the meter object to obtain a pair of uncertain numbers for the voltage readings
Vs = m.voltage(5.0100,tag='Vs')
Vx = m.voltage(4.9885,tag='Vx')

The code to calculate the unknown resistance and display the results (apart from the Meter class definition) now looks like this
m = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7)

Rs = 1000 + ureal(0,1E-3,label='Rs')

Vs = m.voltage(5.0100,tag='Vs')
Vx = m.voltage(4.9885,tag='Vx')

Rx = Rs * Vx/Vs

print "Rx = %G, u(Rx) = %G" % (value(Rx), uncertainty(Rx))
for cpt in reporting.budget(Rx,trim=0):
    print "%s: %G" % cpt
we obtain (as before)
Rx = 995.709, u(Rx) = 0.00100562
Rs: 0.000995709
e_ran_Vs: 9.95709E-05
e_ran_Vx: 9.95709E-05
e_zero: 8.5657E-07
e_gain: 0
This is easier to read, write and understand than the code I posted earlier.

Resistance measurement with two meters

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.

So with just a few simple changes to the code we have
m1 = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7,tag='m1')
m2 = Meter(u_e_gain=3E-6,u_e_zero=1E-6,u_e_ran=1E-7,tag='m2')

Vs = m1.voltage(5.0100,tag='Vs')
Vx = m2.voltage(4.9885,tag='Vx')

Rs = 1000 + ureal(0,1E-3,label='Rs')

Rx = Rs * Vx/Vs
print "Rx = %G, u(Rx) = %G" % (value(Rx), uncertainty(Rx))
for cpt in reporting.budget(Rx,trim=0):
    print "%s: %G" % cpt
which gives the same results as we obtained last time
Rx = 995.709, u(Rx) = 0.0043516
e_gain_m1: 0.00298713
e_gain_m2: 0.00298713
Rs: 0.000995709
e_zero_m2: 0.000199601
e_zero_m1: 0.000198744
e_ran_Vs: 9.95709E-05
e_ran_Vx: 9.95709E-05

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.

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.

3 March 2013

Formulating problems (2)

Systematic errors

In the previous post (Feb 10, 2013), 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.

A simple voltmeter

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
x_\mathrm{dis}  = (1 + e_\mathrm{gain}) V + e_\mathrm{zero}
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
x_\mathrm{dis} = (1 + e_\mathrm{gain} + e_\mathrm{ran}) V + e_\mathrm{zero}
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} \).

How can this be handled with GTC?


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.

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.

A voltage ratio measurement

We will consider a simple measurement in which a resistance value \( R_\mathrm{x} \) is estimated by measuring a voltage ratio.

The circuit is shown in this figure,

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.

With \( R_\mathrm{s} \) known, \( R_\mathrm{x} \) can be calculated
R_\mathrm{x} = R_\mathrm{s}\frac{V_\mathrm{x}}{V_\mathrm{s}}
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}\).

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.

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.

I will define uncertain numbers for the estimates of \( R_\mathrm{s} \), \( V_\mathrm{s} \) and \( V_\mathrm{x} \) and then work out
R_\mathrm{x} = R_\mathrm{s} \frac{V_\mathrm{x}}{V_\mathrm{s}}
Here is the standard resistance
R_s = 1000 + ureal(0,1E-3,label='R_s')
and here are the voltmeter errors
e_zero = ureal(0,1E-6,label='e_zero')
e_gain = ureal(1,3E-6,label='e_gain')
Uncertain numbers representing the individual voltage measurements are then (note an uncertain number for the meter noise is created independently for each reading)
e_ran_1 = ureal(1,1E-7,label='e_ran_1')
V_s = (5.0100 - e_zero) / (e_gain * e_ran_1) 

e_ran_2 = ureal(1,1E-7,label='e_ran_2')
V_x = (4.9885 - e_zero) / (e_gain * e_ran_2) 
And the unknown resistance is
R_x = R_s * V_x/V_s
The value, the uncertainty and the uncertainty budget of the result are displayed by
print "R_x = %G, u(R_x) = %G" % (value(R_x), uncertainty(R_x))
for cpt in reporting.budget(R_x,trim=0):
    print "%s: %G" % cpt
We obtain the following results:
R_x = 995.709, u(R_x) = 0.00100562
R_s: 0.000995709
e_ran_1: 9.95709E-05
e_ran_2: 9.95709E-05
e_zero: 8.5657E-07
e_gain: 0
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.

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.

Different meters will have different residual gain and zero setting errors. So the calculation must now define uncertain numbers for each meter.

The code now looks like this
R_s = 1000 + ureal(0,1E-3,label='R_s')

# Meter 1
e_zero_1 = ureal(0,1E-6,label='e_zero_1')
e_gain_1 = ureal(1,3E-6,label='e_gain_1')

# Meter 2
e_zero_2 = ureal(0,1E-6,label='e_zero_2')
e_gain_2 = ureal(1,3E-6,label='e_gain_2')

# Readings
e_ran_1 = ureal(1,1E-7,label='e_ran_1')
V_s = (5.0100 - e_zero_1) / (e_gain_1 * e_ran_1) 

e_ran_2 = ureal(1,1E-7,label='e_ran_2')
V_x = (4.9885 - e_zero_2) / (e_gain_2 * e_ran_2) 

# Result
R_x = R_s * V_x/V_s

# Display results
print "R_x = %G, u(R_x) = %G" % (value(R_x), uncertainty(R_x))
for cpt in reporting.budget(R_x,trim=0):
    print "%s: %G" % cpt
We obtain the following results:
R_x = 995.709, u(R_x) = 0.0043516
e_gain_1: 0.00298713
e_gain_2: 0.00298713
R_s: 0.000995709
e_zero_2: 0.000199601
e_zero_1: 0.000198744
e_ran_1: 9.95709E-05
e_ran_2: 9.95709E-05
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.

In summary

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.

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.

The next post will revisit this example. I will show then how object-oriented programming techniques can be used to implement the voltmeter model.

10 February 2013

Formulating problems (1)

Developing an uncertainty calculation with GTC

The traditional approach to uncertainty calculation described in the GUM requires a mathematical model of the measurement
Y = f(X_1, X_2,\cdots)
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.

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.

An example should make this clearer.

Consider a measurement of current through a standard resistor (also called current shunt) made using a voltmeter, as shown

If the resistance \( R \) of the current shunt is known, \( I \) can be calculated using Ohm's law
I = \frac{V}{R}
Four voltage readings have been recorded (all in mV)
100.6512 \,, \; 100.6401 \,, \; 100.6420 \,, \; 100.6101 \,,
and the resistor has been calibrated
R = 0.009011 \, {\Omega} \quad \frac{u(R)}{R} = 0.0003
So the current can be calculated using uncertain numbers for the voltage and resistance estimates
x_obs = [100.6512E-3, 100.6401E-3, 100.6420E-3, 100.6101E-3]
V = type_a.estimate(x_obs,label='V_type_a')

e_R_cal = ureal(1,0.0003,label='e_R_cal')
R = 0.009001 * e_R_cal

I = V/R
I.label = 'I'
print summary(I)
The result (in amperes) is
I: 11.1805, u=0.0035, df=465.6

Enhancing the model: temperature dependence

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.

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
e_R_T = ureal(1,4E-6*type_b.uniform(3),label='e_R_T')
R = 0.009001 * e_R_cal * e_R_T
Repeating the calculation
I = V/R
I.label = 'I'
print summary(I)
we obtain the same result
I: 11.1805, u=0.0035, df=465.6
Apparently, the temperature does not have a big influence. To investigate further, we can look at the uncertainty budget for \( I \)
for cpt in reporting.budget(I,trim=0):
    print "%s: %G" % cpt
This displays the magnitude of the component of uncertainty in \( I \) due to each of the three influence quantities
e_R_cal: 0.00335416
V_type_a: 0.000990883
e_R_T: 7.74609E-05

Enhancing the model: voltmeter errors

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
x_\mathrm{dis}  = (1 + e_\mathrm{gain}) V + e_\mathrm{zero}
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\) ).

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 \).

The measurement model can be extended by adding two more uncertain numbers
e_gain = ureal(0,25E-6,label='e_V_gain')
e_off = ureal(0,500E-6,label='e_V_off')
V = (type_a.estimate(v_obs,label='V_type_a') - e_off ) /(1 + e_gain)
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).

Recalculating the current now
I = V/R
I.label = 'I'
print summary(I)
for cpt in reporting.budget(I,trim=0):
    print "%s: %G" % cpt
we obtain
I: 11.181, u=0.056, df=inf

e_V_off: 0.0555494
e_R_cal: 0.00335416
V_type_a: 0.000990883
e_V_gain: 0.000279513
e_R_T: 7.74609E-05
In this case, the voltmeter offset uncertainty makes a significant contribution to the current measurement.

Putting it all together

The GTC code now looks like this:
e_R_cal = ureal(1,0.0003,label='e_R_cal')
e_R_T = ureal(1,4E-6*type_b.uniform(3),label='e_R_T')

R = 0.009001 * e_R_cal * e_R_T

x_obs = [100.6512E-3, 100.6401E-3, 100.6420E-3, 100.6101E-3]

e_gain = ureal(0,25E-6,label='e_V_gain')
e_off = ureal(0,500E-6,label='e_V_off')

V = (type_a.estimate(x_obs,label='V_type_a') - e_off ) /(1 + e_gain)

I = V/R
I.label = 'I'

print summary(I)

for cpt in reporting.budget(I,trim=0):
    print "%s: %G" % cpt
Try to imagine writing down the measurement model for this problem as a simple function \( f() \). It is not straight forward.

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.

That makes GTC a very convenient tool for doing simple uncertainty calculations.

20 January 2013

Reporting results

Obtaining numerical results from uncertain numbers

Every uncertain real number has the attributes:
  • value
  • uncertainty
  • degrees of freedom
Value is the estimate of a quantity of interest. Uncertainty is the standard uncertainty of this estimate and the degrees of freedom is, well, the degrees of freedom associated with the uncertainty.

The numerical values of these attributes are used to calculate an expanded uncertainty for a quantity of interest.

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
y = ureal(1.1,2.5,6)
print value(y)
print uncertainty(y)
print dof(y)
as does the more succinct form using object attributes
print y.x
print y.u
print y.df

Expanded uncertainty

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).

For the uncertain number above (and for a 95% level of confidence)
k_95 = reporting.k_factor(y.df,p=95)
U_95 = k_95 * y.u
Alternatively, the lower and upper bounds of an uncertainty interval can be calculated directly
Y_lb, Y_ub = reporting.uncertainty_interval(y,p=95)
The numbers obtained here are -5.017 and 7.217

Uncertainty budgets

An uncertainty budget breaks down the various contributions to the combined uncertainty of a result.

For example, the following three uncertain numbers are associated with measurements of voltage, current and phase in an electrical circuit
V = ureal(4.999,3.2E-3,label='V')        # volt
I = ureal(19.661E-3,9.5E-6,label='I')    # amp
phi = ureal(1.04446,7.5E-4,label='phi')  # radian
Based on these measurements, an uncertain number for the resistance \( R \) is
R = V * cos(phi) / I
We obtain a budget for this by writing
for cpt in reporting.budget(R):
    print "%s: %G" % cpt
which generates a list of names and components of uncertainty in the same units as the measurement result and in order of decreasing magnitude
phi: 0.164885
V: 0.0817649
I: 0.0617189
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.

Correlation coefficients

Measurement results may become correlated if they share influence quantities, or are calculated from common data sets.

Continuing with the example above, an estimate of reactance is
X = V * sin(phi) / I
Now, the estimates associated with R and X will not be independent, because they were calculated using the same set of data. Their correlation coefficient can be calculated by
(the value obtained is 0.058)

13 January 2013

Defining uncertain numbers

Representing estimates

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.

The first task then, when using GTC, is to define uncertain numbers.

Absolute uncertainty

When an estimate and standard uncertainty are given in the same units, the definition of an uncertain number is straight forward.

Suppose, a length is reported
\widehat{l} = 10.1\,\mathrm{cm},\; u(\widehat{l}) = 0.3\,\mathrm{cm}\;,
an uncertain number for \( \widehat{l} \) may be defined by writing
l = ureal(10.1,0.3) 
although it would be better to write
l = 10.1 - ureal(0,0.3)
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
l = \widehat{l} - e \;,
with \( e \) being the (unknown) measurement error.

Relative uncertainty

Sometimes an uncertainty is expressed relative to the magnitude of an estimate (i.e., as a standard relative uncertainty).

For instance, a length may be reported as
\widehat{l} = 10.1\,\mathrm{cm},\; \frac{u(\widehat{l})}{\widehat{l}} = 3\% \;.
An uncertain number for \( \widehat{l} \) in this case can be defined by writing
l = 10.1 * ureal(1.0,0.03)
which corresponds to the relation
l = \widehat{l} \cdot e_\mathrm{rel} \;,
with \( e_\mathrm{rel} \) an unknown relative error (estimated as unity).

Expanded uncertainty

Whether relative or absolute, uncertainties are more often reported "expanded", than as standard uncertainties.

An expanded uncertainty is just a standard uncertainty, or a relative standard uncertainty, multiplied by a coverage factor. Formally, the expanded uncertainty is represented by a capital \( U \) and the coverage factor by \( k_p \), so for length we might write
U(\widehat{l}) = k_p \cdot u(\widehat{l}) \; ,
\( p \) is a coverage probability, or level of confidence, associated with the expanded uncertainty.

Expanded uncertainties should never 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.

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
U_l = 0.63
k_95 = 2.1
The standard uncertainty and degrees of freedom can be easily found (note, 'k_to_dof' is a function in the GTC 'reporting' module)
u_l = U_l / k_95
df_l = reporting.k_to_dof( k_95, 95 )
and an uncertain number for \( \widehat{l} \) can be defined with these values
l = 10.1 - ureal(0,u_l,df_l)