Applications of Sparse L1 Pursuits to Precision Reference

[rhb]

The  properties of sparse L1 pursuits were initially presented in a flurry of papers by David Donoho and his former student Emmanuel Candes.  Candes had done an experiment in which he attempted to fit random signals with a combination of Fourier series and impulse spikes.  The results were far better than expected and led both to an intense effort spanning most of 2004.

The best summary of the importance of their work is stated in the introduction to this paper:

https://statweb.stanford.edu/~donoho/Reports/2004/l1l0EquivCorrected.pdf

Unless you enjoy mathematics a lot, I suggest just reading the introduction.  Why it works is not as important as that it does work and that software for applying the technique is readily available.

What Donoho proves in the paper cited above is that a least summed absolute error (L1) solution to an ordinary system of linear equations Ax=y is exactly equivalent to the optimal (L0) solution which is NP-hard and can only be obtained by exhaustive search.  Previously, the only known solution to such problems was to try every possible combination of the columns of A.  For significant problems this is impossible.  On the fastest computer available it would take until all the stars in the universe burned out.

This has profound implications for a great many applications.  In the context of metrology, it permits solving problems which account for all the known and characterized system errors.  In the case of errors due to quantum effects which can only be described statistically, it provides tight bounds on the irreducible error.

So if we parameterize the system output  of a voltage reference and the associated buffer amplifier with an expression such as:

expected_value = vref_initial_value + ref_aging + vref_1/f_noise + vref_temp + amp_initial_offset + amp_aging + amp_temp + r1_initial + r1_aging + r1_temp + r2_initial + r2_aging + r2_temp + r3_initial + r3_aging +r3_temp

Given reasonable approximation of the values in the summation as arbitrary functions of measured parameters and unknown coefficients, we can compute the functions for a large number  of possible coefficients.over a wide range.  A problem with 10,000 choices for each term is solvable in a matter of a few minutes by means of linear programming using the simplex method.  Faster algorithms based on the properties of regular polytopes in N dimensional space exist.  But the simplex method, invented by Dantzig  to solve optimization problems in operations research does an excellent job and high quality FOSS  software is readily available.

If one computes 10,000 possible instances for each of the terns above as functions of measured values such as temperature, age, etc, there are 10,000**16 possible combinations.  An L0 solution is computationally intractable.  It can't be done.  But an L1 solution, if it exists, has been proven to be identical if the columns of A are uncorrelated in all combinations and the x vector in x of Ax-y are sparse, that is, most of the coefficients of x are zero.

To put this in concrete terms, suppose we know from testing that the value of a resistor is f(a,b,c)  where:
a=initial value
b= change due to ambient temperature
c=change due to current flowing through the resistor

but that we do not know  a, b or c.  All we know is that the change due to ambient temperature  is some function g(b) and the change due to heating by the current is h(c).  If we know the form of h(c) as say a polynomial h(c) = k0 + k1*c + k2*((c-d)**2) from measuring a number of resistors where k0, k1 and k2  are all vary with each device, then we can evaluate a large number of values of h(c) for a range of values of k0, k1 and k2. All of these are added to the A matrix as possible answers.  The "dictionary" in the terms of sparse pursuits.

So what one does in practice is create models of the functional forms of the various terms.  One then builds a massive dictionary of all the possible values for each term as a function of the known values.  Unknown values are simply varied over a range of expected values and the expression evaluated.  In the example cited, with 10,000 instances for each term the A matrix would consist of 160,000 columns and as many rows as one had measurements.  When i was in grad school 30 years ago and are computers consisted of an 11/780 and a MicroVAX II solving such a problem was not possible even if we had known how.  Now it's no harder than computing a long FFT on a desktop machine.

[texaspyro]

That made my noggin throb...   

[TiN]

It's all fun and dandy, but do we really need third thread on essentially same topic? 

[ramon]

Yes, we do need a third thread! I really want to know if we can predict noise with that sparse L1 evil thing.

[zhtoor]

Yes, we do need a third thread! I really want to know if we can predict noise with that sparse L1 evil thing.

it is not about predicting noise, it is about predicting the drift pattern based on some sophisticated model and adaptive predictors (like kalman filters) if at all.

regards.

-zia

[IconicPCB]

... and the rope gets longer...

[Conrad Hoffman]

There are three kinds of people, those that understand advanced math, and those that don't.
 

[Theboel]

There are three kinds of people, those that understand advanced math, and those that don't.
 

and the third are ?

[ramon]

Zia, the words were not randomly chosen. Why do you make such difference between drift and noise?

[zhtoor]

Zia, the words were not randomly chosen. Why do you make such difference between drift and noise?

look at LTZ1000 datasheet, on the first page, look at "FEATURES", and i quote:-

1.2uVP-P Noise
2uV/rtkHr Long-Term Stability
Very Low Hysteresis
0.05ppm/°C Drift
Temperature Stabilized
400°C/W Thermal Resistance for LTZ1000A Reduces Insulation Requirements
Specified for –55°C to 125°C Temperature Range
Offered in TO-99 package

best regards.

-zia

[rhb]

Sparse L1 pursuits are far more general than the particular application of the OSHW voltage reference which is an inverse problem.  It is the basis of the solution that won the Netflix prize for best predictor of other movies people would like.  In mathematical jargon that is called "matrix completion".  Other applications in include error correction, machine learning and compressive sensing.

I started this thread specifically to address some ideas for improving high precision measurements. I think it unhelpful to clutter the OSHW Vref thread with details which are mostly germane to sub-ppm cases.  My intent with the OSHW project is a cost-performance design effort.  This thread is is best possible performance. National lab level work.

I also wanted to separate out the mathematical details so I didn't find myself writing the same thing multiple times to explain a point.

Consider measurement of a standard resistor in an oil bath with precise control of the temperature and  measurement current using a state of the art voltmeter, the details of which will not be treated except in the abstract.  Initially I shall presume that the meter is perfect.  Such things do not exist, but I want to keep the initial discussion as simple as possible.  The complications  of real meters will be dealt with later.

No matter how precisely we control the oil bath, current flowing through the resistor will produce self heating and unavoidable thermal gradients in the measurement path.  The current under consideration will be limited such that it does not cause permanent damage to the resistor.

When the measurement is started, the resistance will change as the current heats the device.  This is a simple curve which asymptoticaly approaches a steady state value.  However, we should have to wait forever to reach that value. The mathematical form is an instance of Fourier's heat equation, so the asymptotic behavior is well studied.  Measurements during the start of the run will allow very precisely determining the value at any instant in the course of the measurement run.

For the purpose of this discussion, I shall stipulate a wirewound resistor.  With that construction the resistor also has an inductance and an interwinding capacitance.  That might seem irrelevant to a DC measurement, but quantum effects come into play.  In particular, the current source will have a number of quantum level noise components, 1/f being the most troublesome.  The resistor itself is also a source of noise.  As the latter noise source is broadband, the parasitic inductance and capacitance of the resistor form frequency dependent filters that will alter the value of the noise at certain frequencies.  These effects will be ignored for now.

To make a precise measurement of the resistance using a precision current source we need to account for a number of small, but vexing problems.  We have thermal noise,  i/f noise and more.  At present I know about the form of the solution but not the equations of the problem. I shall limit myself here to just the thermal noise and the 1/f noise.  The first two layers of the onion we wish to peel apart.

If we take a long series of measurements which start before power is applied to the resistor we will observe the effect of heating of the resistor and we will observe all the various quantum level noise effects produced by the physical devices in the circuit,

If we Fourier transform this series of measurements with our ideal meter we shall see the 1/f noise of the current source, the thermal noise of the resistor which is at a slightly higher temperature than the oil bath as it is a heat source and  the Fourier components of the resistance change induced by the current flowing through the resistor.  The latter will change  over time as the resistor heats up, but I shall ignore that complication for now.  It can be dealt with at the cost of making the system of equations larger.

At present I shall assume that the 1/f noise is the least well determined error term and that the self heating is adequately described by the heat equation so that we need only determine the asymptotic value and the time constant.  The equation is precisely known and we can calculate the value during the course of the measurement for a large number of possible values at infinite time and a range of time constants.

The thermal noise of the resistor is also known to a degree, but we do not know the exact temperature of the resistor, just the temperature of the oil bath. Again, we can compute the power spectral density of the thermal noise for a range of temperatures covering the expected values based upon the temperature of the oil bath.

This leaves the 1/f noise of the current source.  This is the most difficult term in our problem.  Obviously the noise does not become infinite after an infinite period of time.  So clearly, the correct functional form is more complex than 1/(f**a), though that may suffice for measurements over practical time periods.  A variety of more complex forms have been proposed.  Once again we can compute a large number of estimates of the power spectrum density for a range of coefficients and estimation functions.

With the exception of the self heating error, the other components are essentially random. We can make some general statements about the expected values of the power spectrum such as the mean and standard deviation but no more.

So our mathematical problem consists of D(f), the Fourier transform of the measurements, and the dictionary of error term models for the resistor and the current source.

D(f) = R(0) + F(f) + S(f) + G(f)

Where:

R(0) is the DC component of the transform

F(f) is the 1/f noise

S(f) is the transform of the self heating

G(f) is the approximately Gaussian thermal noise

Thus we have Ax=y with a rather large A matrix and an x vector which should contain 4 non-zero elements which give us the DC resistance, the coefficients of the 1/f noise, the actual device temperature and the time constant and resistance change of the resistor due to self heating.  The next layer of the onion is to account for the voltmeter at a similar level of detail by adding the appropriate equations to the system.  To measure the drift the problem outlined here would need to be be modified to to use R(f) instead of R(0).  As a practical matter solving for R(0) once per measurement run and then using those results to solve for the drift equation is more tractable to set up and solve.

To the best of my knowledge no one has attempted anything comparable.  Should anyone know of work at this level of refinement I should be most grateful for references and contact information.

[MisterDiodes]

RHB:

So lets get back to a practical application of "Chasing a few little ones" - And how it works on Reg's Cheap Calibrator for example, or even on a good Vref.

Without dropping author names, your credential history or equations, walk me through a one paragraph explanation of why you think this might work on a system that also includes an unknown chaotic history based real worlld situations:

Let's say you find a Vref out of the blue that you've supplied a couple years prior and maybe at some point a few of it's "L1's have been pursued".  In other words you think you might be able to predict it's actual Absolute Value of Voltage to some level of uncertainty at some point in the future.

BUT - There is no known history of the device in the few years it's been in the field.  You have no idea the time-integral effect of load, humidity, temperature effects, voltage bias points, etc, etc, etc.  Was the device shipped across the Arctic freezing on an airplane or did it sit in a shipping container for months, or was it in constant use the whole time in a lab.  All of that is unknown.

A really good one that throws the college interns claiming to be able to predict a future value is having them look at the RATE of CHANGE of temperature on a Vref or resistors, and you can also look at the Jerk effects as well (rate of rate of change).  See if that's predictable and how that affects long term drift.  That's an effect that we look at in real world situations too, all of which have long-lasting effects on variables that you could have sworn you had a handle on.

The point is: you are holding in your hands a device on which you have NO idea of it's past use history or environmental exposure:  How does an L1 pursuit technique predict it's current value, especially when you have absolutely no clue about it's past history?

[zhtoor]

RHB:

So lets get back to a practical application of "Chasing a few little ones" - And how it works on Reg's Cheap Calibrator for example, or even on a good Vref.

Without dropping author names, your credential history or equations, walk me through a one paragraph explanation of why you think this might work on a system that also includes an unknown chaotic history based real worlld situations:

Let's say you find a Vref out of the blue that you've supplied a couple years prior and maybe at some point a few of it's "L1's have been pursued".  In other words you think you might be able to predict it's actual Absolute Value of Voltage to some level of uncertainty at some point in the future.

BUT - There is no known history of the device in the few years it's been in the field.  You have no idea the time-integral effect of load, humidity, temperature effects, voltage bias points, etc, etc, etc.  Was the device shipped across the Arctic freezing on an airplane or did it sit in a shipping container for months, or was it in constant use the whole time in a lab.  All of that is unknown.

A really good one that throws the college interns claiming to be able to predict a future value is having them look at the RATE of CHANGE of temperature on a Vref or resistors, and you can also look at the Jerk effects as well (rate of rate of change).  See if that's predictable and how that affects long term drift.  That's an effect that we look at in real world situations too, all of which have long-lasting effects on variables that you could have sworn you had a handle on.

The point is: you are holding in your hands a device on which you have NO idea of it's past use history or environmental exposure:  How does an L1 pursuit technique predict it's current value, especially when you have absolutely no clue about it's past history?

lets simplify things a little bit to make this discussion more fruitful:-

assume:-

1. a JJA is available right next to your setup. (for periodic / aperiodic transfers)
2. you are running a setup which includes some kind of a predictor (L1 persuits, kalman filters, sensor fusion, neural networks etc.)
3. now you have a unit under test (UUT) which is being monitored and its predictor updated at pre-deteremined times.
4. a particular UUT is characterized using this method and dominant environmental parameters established from model analysis.

and then:-

build a UUT with the dominant parameter logging / recording built-in to enable UUT transport / transfers.

best regards.

-zia

[rhb]

If one has no knowledge about past device history, the only thing one can do is ascertain the effects of noise on the measurement in a more sophisticated fashion than merely averaging values over a number of power line cycles.  Even the most modest consideration should make that obvious.

This is NOT the OSHW Vref project thread.  It is about developing  system models for precision references and using those to improve the accuracy of the measurements.  My example was intended to illustrate the concept, not be a comprehensive description of the physics of  a resistor or the mechanics of solving the problem.   To the best of my knowledge, this is a national lab level subject which has had little or no investigation. 

As you seem to resent mathematical analysis, perhaps you should sit this one out.  I believe the jargon is, "This is a league game".  This thread is about the application of recent developments in applied mathematics to state of the art problems in metrology. 

[MisterDiodes]

I guess I'm still confused when you started with a $50 Vref with 4ppm drift using cheap resistors and "...Being...ignorant", and then the discussion moved to describe a $150 calibrator for 5.5 digit DMMs...and now you're doing an NIST project.  OK, I'm just trying to keep up. I thought the threads were somehow related.

Don't get me wrong - I'm always interested in a REALISTIC and ROBUST mathematical approach to make GOOD measures and equipment BETTER.  We use those techniques all the time!.   The reason we like to have fun with the summer college interns is hopefully they go back to university with eyes opened a bit - again we try to get them to put down the books for a sec and observe the real world.  Because half the time their profs haven't got a flippin' clue of how or why or when to apply the math.  It's bad enough we get EE kids that have no idea of which end of the soldering iron to hold (They can simulate it though!)...but that's a different story. 

OK, OK - I promise my mind is open, but I'll just pull up a chair and shut up.

Sell it to me!

[rhb]

I should like to pose a few questions.

The procedure I intend to follow  is to develop a system model, make some measurements, setup Ax=y, solve it and then examine how well the model matches future measurements.

I'd like to do this with existing data to start.  TiN has made a lot of data available.  I've not yet looked at more than a small fraction of it.  So I should like him to suggest the best choice.

I want to address 1/f noise and thermal noise in the first iteration.  The data requirements for this are measurement runs of moderate length with temperature data.  More environmental data may be useful, but this is the minimum.  It does not appear that it matters whether the measurements are resistors or voltage references.

For resistors we have the 1/f noise of the ADC voltage reference, the  voltage reference and ADC input buffer amps and the current source and its associated circuitry.  Are there other components which generate 1/f noise?  We  also have thermal noise from a number of resistors besides the DUT.

For voltage references the noise sources are similar except that we have 1/f noise from the DUT rather than the current source.

Because I propose to address the 1/f noise in the frequency domain, long measurement runs are important. My present candidate would be the cyrogenic test in February just looking at the data when the device was normalized to the liquid, but I shall consider any other data that meets the requirements for setting up the problem.

Whatever problem is chosen, I need an equivalent circuit model which includes all the noise sources.  AS I am not familiar with the details of the instruments I'd much appreciate of someone who is would define the equivalent circuit model,  I shall have to write several programs to generate the GMPL input to the glpsol program fin GLPK and to do that I have to work out how to properly linearize the problem.  So it will take a while to do this.  The input problem files will be many megabytes, so writing them by hand is not a practical undertaking.

I also need  references to the professional and OEM literature.  Typically for a serious project I will collect and read 50-100 papers.  I usually try to find a few good recent papers and then go back through the references in those.  I shall be spending sometime at the nearest university library on the 10th, so if anyone can suggest a few good starting points related to 1/f noise in voltage references, zener diodes and op amps that would be very helpful.

I'd like to keep a clean separation between this and the OSHW project.  This is a science project.  The OSHW is an engineering project.  If this effort is successful then it can be applied to the OSHW project during the development of the UI and calibration component of that device.

[MisterDiodes]

You can look at a PWW resistor which has virtually no 1/f noise to begin with...That might narrow down the variables you're looking at on a first pass, and allow you to concentrate on other effects.

EDIT:  Use the same approach with the amplifier for your first tests:  Use a chopper amp to eliminate 1/f noise contributed by the amp.

This might be something to look at in your L1 Pursuit / OR1MP experiments as well, from an overall efficiency point of view:  At what point does spending slightly more money (overall) on good components leave a smaller carbon footprint than using less stable components that require a lot of CPU cycles (= Energy = Cost) forever trying to compensate for more noise and drift?

[rhb]

Mathematically, I don't think the complexity matters.  Mostly I want to avoid spending a huge amount of time researching the noise characteristics of numerous  devices.  However, I do want to be complete to first order.  The most important aspect is to use existing measurements for the first few trials.  So error terms for which the associated variables were not measured will be left out.

Not sure what you mean by "OR1MP".  Could you please explain?  One only need  solve the inverse problem when a new calibration is done.  Otherwise it's a very simple evaluation of Vdc = F(temperature, age, etc).  Something very trivial for an MSP430 to do for a few picowatts.  The RTC and display will consume vastly more power.

There is no way to compensate for the noise.  It's important to include it though.  If this effort succeeds it  will dominate the error bars.  The point of including it in the system model is to remove the effects of noise from the aging function estimate and the estimate of the true value.

I found a paper on resistor noise which looks to be fairly complete.  I'd be grateful for a review with particular attention to omissions or errors in the equations as well as links to other papers.

https://dcc.ligo.org/public/0002/T0900200/001/current_noise.pdf

This thread is really about improving the calibration process for electrical measurements and references, particularly transfer standards. I'll leave the national lab staff to apply the techniques to primary references if the method succeeds with transfer standards.  If it meets my expectations for aging and temperature corrections I'll take a crack  at the hysteresis effects which plague traveling references. That may not be solvable beyond a very  limited degree.

For the first trials I'm looking for  a nice smooth, repeatable temperature curve separated from noise curves which are uncorrelated from run to run.

Errors in the DMM will be largely uncorrelated with the reference being measured.  So by reading two or more references during a run the low frequency portion of the DMM errors can be separated and corrected.

The mathematics of doing this are just basic algebra.  The mathematics of the solution are ugly, but there's no real reason to go into that in detail unless edge cases arise.

[MisterDiodes]

OR1MP = Orthogonal Rank One Matrix Pursuit (or just Rank One Matrix Pursuit) which is not the same as your Sparse L1 Pursuit I know, but produces -very- similar results, at least on the places where I've seen it used.  As you would know.  That's what I thought of based on your initial description, and I got that stuck in my head.

That paper you reference above even mentions on Page 8 that PWW and BMF resistor have little to no excess (1/f) noise.  I'm not sure but I don't think they have many PWW resistors represented on their graphs, but double check.  A few of the model numbers I'm not familiar with.  One of the mistakes in that experiment of course is the use of an aluminum shield box, which does just about nothing for low freq EMI - which will find it's way into the resistor test.

The man who wrote the book on PWW is here on the forum, Edwin Pettis, and I suggest you check his past articles - and the fact that PWW are devoid of 1/f this is also something that's been known for many decades in the older literature (back to the 50's and before), a lot of which is not online.

The same mechanisms that produce VC and 1/f are present in all other resistors because (very simplified) a PWW uses a round conductor, and when it's properly annealed and formed that offers very few areas where electron traps can form, and the atomic bonds are all of very "strong" type.  Any time you try to use a planar resistive element you'll run into many 90° corners of the conductive structure, all of which can form electron traps which add to -chaotic- electron flow (noise) - and flat, non-smooth planes are areas where tend to find a lot of weak atomic bonds (allows atoms to swivel around and align in an E-field), tending to produce VC effects.  In a very simplified nutshell.  We run into these problems all the time on film / diffused resistors built onto a crystal lattice such as GaAs, Si, Sapphire, etc.

But Edwin is the expert in this department, and he can add more.

 

[zhtoor]

hello,

please find attached test bench #1:-

1. there are N voltage references operating in tandem.
2. the N voltage references are connected to switch groups SA(1..N) and SB(1..N).
3. any of VREF(1..N) can be selected under automated control of SA(1..N) to produce average of the selected refs on line VA.
4. any of VREF(1..N) can be selected under automated control of SB(1..N) to produce average of the selected refs on line VB.
5. line VA and VB is measured differentially and analyzed.
6. it is assumed in the start that switches SA(1..N), SB(1..N) and resistors RA(1..N), RB(1..N) are noise/thermal emf free.
7. the maximum difference between VREF(1..N) is under K ppm.

to solve:-
1. which subset of the VREF(1..N) is an outlier and by what margin?
2. which subset(s) of the VREF(1..N) form a group in such a way that their combined drift is a minimum?
3. what are the individual drift patterns vs. the most stable group?

best regards.

-zia

[rhb]

@MisterDiodes  I'd never heard of OR1MP.  I'll need to look into it.  There were a number of sparse L1 pursuits in use by a wide variety of names before the work of Candes and Donoho. 

@zhtoor  Geez!  Did you need to bring such a *big* baseball bat? ;-) 

I figured out how to remove the measurement system reference error.  That's correlated between the measurements  You have posed a seriously difficult problem.  Do you know of a solution or is this an open problem?

I've been contemplating this for N = 3 without success. I haven't been able to determine even whether it can be solved.

I'll work on it some more.

[rhb]

Professor zhtoor :-)

If we measure the output of N references over a sufficiently long period of time and crosscorrelate the N references with the N-1 other references, the random noise will have zero correlation.  Therefore the references with the largest crosscorrelation coefficient will have the most similar drift.  This result is independent of the thermal noise of the switches and resistors.  So it applies even if those are non-zero.

The selection of the best set of references may be done by sorting the references by crosscorrelation coefficient and selecting the references above the larger coefficient knee in a sorted plot of the set of all the crosscorrelations. 

1) The outliers may be found by comparing the crosscorrelation of each device with the pair of devices with the highest correlation coefficients. 

2) Compute the Allan deviation of the mean of all P subsets of N, 1< P <= N,  the drift is a minimum for the set of size P whose mean has the most constant Allan deviation 

3) Compute the Allan deviation of the difference of the individual references and the mean of the most stable group

QED
Reg

Actually easier than I thought.  But rather intimidating on first reading.  A *very* good question.

Edit: NB: The preceding analysis has nothing to do with sparse L1 pursuits.  It is a classical analysis based on the work of Norbert Wiener.  The best reference on the topic is "Random Data" by Bendat and Piersol.   I can think of no monograph on applied mathematics which is more useful for practical work.  I have a substantial fraction of the classic works on time series analysis.  Bendat and Piersol are superior by an order of magnitude.

Edit 2:  I should point out that the above assumed that the series were all of the same length.  If they are not of the same length, the statements still apply, but correctly computing the crosscorrelation coefficients becomes more nuanced.  This is particularly the case for long series for which using the FFT is desirable.  One of these nuances handed my head to me on a platter 20 years ago.  It was a very anxious week or so before I figured out what was going on.

[Svgeesus]

But an L1 solution, if it exists, has been proven to be identical if the columns of A are uncorrelated in all combinations and the x vector in x of Ax-y are sparse, that is, most of the coefficients of x are zero.

To put this in concrete terms, suppose we know from testing that the value of a resistor is f(a,b,c)  where:
a=initial value
b= change due to ambient temperature
c=change due to current flowing through the resistor

The constraint is that the columns are uncorrelated, but your example has two heating terms. Doesn't the heating due to current affect the local, ambient temperature?

[rhb]

But an L1 solution, if it exists, has been proven to be identical if the columns of A are uncorrelated in all combinations and the x vector in x of Ax-y are sparse, that is, most of the coefficients of x are zero.

To put this in concrete terms, suppose we know from testing that the value of a resistor is f(a,b,c)  where:
a=initial value
b= change due to ambient temperature
c=change due to current flowing through the resistor

The constraint is that the columns are uncorrelated, but your example has two heating terms. Doesn't the heating due to current affect the local, ambient temperature?

A very good question.  Solutions of the heat equation, a*F(kt) are only very slightly correlated for different values of k.  The ambient temperature term is not correlated with the self heating term. However, measurement during warmup is essential to solving for c.  Once the self heating approaches steady state they cannot be separated.

 I spent a good bit of time investigating the correlation of solutions of the heat equation  long before I learned about sparse L1 pursuits. I studied sparse L1 pursuits because I had been taught you could not solve the problems I had been solving.

[zhtoor]

Professor zhtoor :-)

If we measure the output of N references over a sufficiently long period of time and crosscorrelate the N references with the N-1 other references, the random noise will have zero correlation.  Therefore the references with the largest crosscorrelation coefficient will have the most similar drift.  This result is independent of the thermal noise of the switches and resistors.  So it applies even if those are non-zero.

The selection of the best set of references may be done by sorting the references by crosscorrelation coefficient and selecting the references above the larger coefficient knee in a sorted plot of the set of all the crosscorrelations. 

1) The outliers may be found by comparing the crosscorrelation of each device with the pair of devices with the highest correlation coefficients. 

2) Compute the Allan deviation of the mean of all P subsets of N, 1< P <= N,  the drift is a minimum for the set of size P whose mean has the most constant Allan deviation 

3) Compute the Allan deviation of the difference of the individual references and the mean of the most stable group

QED
Reg

Actually easier than I thought.  But rather intimidating on first reading.  A *very* good question.

Edit: NB: The preceding analysis has nothing to do with sparse L1 pursuits.  It is a classical analysis based on the work of Norbert Wiener.  The best reference on the topic is "Random Data" by Bendat and Piersol.   I can think of no monograph on applied mathematics which is more useful for practical work.  I have a substantial fraction of the classic works on time series analysis.  Bendat and Piersol are superior by an order of magnitude.

Edit 2:  I should point out that the above assumed that the series were all of the same length.  If they are not of the same length, the statements still apply, but correctly computing the crosscorrelation coefficients becomes more nuanced.  This is particularly the case for long series for which using the FFT is desirable.  One of these nuances handed my head to me on a platter 20 years ago.  It was a very anxious week or so before I figured out what was going on.

thanks. i am no professor, actually far from being one.

what would be the an efficient (optimal?) way of achieving the above?

actually, this is an attempt to generalize the 12-ball weighing problem for voltage reference characterization.
the real problem is the "weighing procedure", ie; how long to measure, where are the thresholds?

best regards.

-zia

[rhb]

Well, I felt as if I'd been called out in class ;-)

Given:

Aj, Bj, Cj..... Nj  at each round of measurements, calculate all the pair wise products and add to the accumulated correlation coefficient.

Phi(A,B) += Aj*Bj is the running estimate of the crosscorrelation coefficient of A & B.  If I were doing this I would put all the correlation coefficients in a file as ASCII text:

number pair

sort -k 1n <file >sorted

and then plot sorted using gnuplot. In general one gets a plot of the cumulative distribution function.  As most data are quasi-Gaussian this takes on a sigmoidal shape with the outliers at the right and left.

The optimal measurement length is as long as possible, as often as possible.  However, daily over a few years is adequate.  The lower the drift, the longer the series needed to characterize it.

I should like to note that the measurement drift of the DVM is correlated across devices.  That introduces an error in the correlation coefficient estimate.  I'm still contemplating the particulars, but a sparse L1 pursuit should be able to separate the DVM drift from the reference drift if there are multiple references.  The part I've not determined is how the number of references in excess of two affects the accuracy of the DVM drift estimate.  I'm also a bit unclear if it's possible using Wiener's L2 methods.  Likely one can, but the procedure may be rather involved.

If there is interest from TiN, cellularmitosis, Andreas and others who are measuring and comparing multiple references I'll write a bespoke program to do the calculations and generate files for use with gnuplot.

I know little about the Allan deviation as I have never coded it, just looked at plots of clocks.  So I'll need to read the definition closely and work out the arithmetic required.  The article I checked suggested it was far simpler than I thought.

If you are at all interested in problems such as posed by Zia, you really should get a copy of "Random Data" by Bendat and Piersol.  It is very readable and describes in clear language how to solve most problems.  There are a few things that sparse L1 pursuits can do that can't be done with classical L2 methods, but pursuits are not a replacement for classical methods. Allan Piersol passed away, so the 4th ed is the last.  Octave/MATLAB will handle all the numerical chores.

[zhtoor]

I should like to note that the measurement drift of the DVM is correlated across devices.  That introduces an error in the correlation coefficient estimate.

this might help.

http://rubiola.org/pdf-articles/archives/2010-arxiv-1003.0113v1-xspectrum.pdf

best regards.

-zia

[zhtoor]

hello,

refer to test bench #1 as posited above (reply #19):-

1. for a given k <= floor(N/2) select a random subset of k VREFS(1..N)
2. select the switch group A to this subset.
3. select the switch group B to a random subset of remaining VREFS. (N maybe odd)
4. do a measurement / analysis on lines VA, VB.
5. repeat 1..4

now:-

1. can we detect an outlier and in how many measurement cycles ?
2. can we draw a curve of the number of measurements required to achieve (1) vs. N/k, what form would it be ?

best regards.

-zia

[rhb]

I'm going to restate the problem:

Given N references and two DVMs, how can the worst references be found in the least time?

divide the references in two groups of approximately equal size

compute the Allan deviation of the two groups

the group with the largest deviation contains an outlier

divide the largest deviation group into two groups and repeat until the device with the largest drift is found

This is search by bisection which is guaranteed to converge in P steps for 2**P <= N.  The length of time it takes depends upon how large the drift of the outlier is.  Each step in the search must run long enough to produce a robust estimate of the variance.  In term of wall clock time to find the best and worst references measuring them individually and computing the cross spectra is the best approach.  Finding the next worst reference would require completing the search again with the worst reference removed from the pool.

Which brings up the problem of multiplexing the input to a DVM.  What's the best way to do this?  FET switches/multiplexers?  Relays?  What type of relays?

I'd like to read N references and N*(N-1) differentials with my 3478A.  In my case N=3, at least at present.  I tried a search of the forum with google, but didn't find anything beyond buying commercial units and it's hard even with commercial units.

As an aside, "Numerical Recipes" includes a number of "optimal" searches.  But they are not guaranteed to be faster than bisection and can have search times which are significantly worse.  I have never understood why anyone would use them unless they know with certainty  a priori that the worst cases cannot arise. I was in grad school when the first edition came out and wildly enthusiastic.  I bought the 2nd, but I haven't looked at either in many years and can't think of a reason I would.  I have *never* used anything but bisection in my professional work.  It is guaranteed to be reliable if you code it correctly.  But it was several years between the first publication of a bisection algorithm and a correct bisection algorithm.  It's maddenlngly easy to get the test edge cases wrong.

[zhtoor]

Which brings up the problem of multiplexing the input to a DVM.  What's the best way to do this?  FET switches/multiplexers?  Relays?  What type of relays?

have a look at this:-

https://www.eevblog.com/forum/metrology/low-thermal-emf-scanner-and-ordinary-non-latching-relays/msg1064900/#msg1064900

best regards.

-zia

[zhtoor]

I'm going to restate the problem:

Given N references and two DVMs, how can the worst references be found in the least time?

how about restating the problem as a compressed sensing problem?

best regards.

-zia

[rhb]

I'm going to restate the problem:

Given N references and two DVMs, how can the worst references be found in the least time?

how about restating the problem as a compressed sensing problem?

Sparse L1 pursuits are a method of solution that allows solving problems traditional L2 methods can't handle.  Amazing as it it is, it doesn't perform miracles.

If we knew with high certainty the functional form of the drift curve and we had measurements from initial startup, it *might* be possible to solve for the individual drift curves from a measurement of the sum.  But it seems unlikely.  Given the above information it would take a week working full time to determine if it was possible at all.  And then even with the answer, which unit is which?

I was just pricing a 3457A w/ 44492A.  The price made me ask why am I doing this?  Building a multiplexer seems more attractive. It's creation rather than consumption.

[rhb]

Which brings up the problem of multiplexing the input to a DVM.  What's the best way to do this?  FET switches/multiplexers?  Relays?  What type of relays?

have a look at this:-

https://www.eevblog.com/forum/metrology/low-thermal-emf-scanner-and-ordinary-non-latching-relays/msg1064900/#msg1064900

best regards.

-zia

That's very interesting.  It appears that it would be pretty easy to correct. Also rather annoying that google didn't find it.

[zhtoor]

I was just pricing a 3457A w/ 44492A.  The price made me ask why am I doing this?  Building a multiplexer seems more attractive. It's creation rather than consumption.

hello,

the hp3456a may be a better and a cheaper choice. a number of pro's on this board swear by it.

best regards.

-zia

[rhb]

The 3457A has a a couple of different multiplexer modules the 44491A & 44492A.  The 3456A doesn't.  But the cheapest I could find was a 3457A w/ 44492A for $500 which is pushing what I want to spend on this.  A 44492A alone was typically about $250-300.

Based on the information in this thread:

https://www.eevblog.com/forum/metrology/low-thermal-emf-scanner-and-ordinary-non-latching-relays/msg1064900/#msg1064900

I shouldn't have much trouble building a  10-12 channel scanner and simply correcting for the thermal effects.  The biggest obstacle is my lack of PCB design tool experience.  I tried doing some simple stuff a few years ago and got really frustrated trying to sort out how to do what I wanted.  I never completed a board design.  This is not something I would build without a full ground plane and guarding.  So unless someone volunteers to design a board, it won't happen soon.

I put a small loop on a scope probe and was able to track down SMPS noise from my PC which I was able to squelch with a piece of hardware cloth, but I'm now trying to figure out the path for the SMPS noise from the Instek MSO-2204EA.  Somehow it is coupling into the probe on the power feed to the basic-lm399 boards.  The power lead is now shielded PTFE twisted pair with the shield connected only at the PSU end.  I have an assortment of clamp on ferrite chokes on order.

Screening the PC PSU was so successful that I'm going to pull the light fixture down and examine the screen bonding and the integrity of the mains ground lead.  It was difficult working on a ladder over my head.

[Kleinstein]

The 3457 is not that bad, though a bid noisy and has odd ranges ( 3 V and 30 V).


PCB design can be quite time consuming, especially with more than 1 layer.

For a scanner it depends on the requirements.  Mechanical relays usually have better isolation, but may not be 100% reproducible in the sub µV range. If some leakage can be tolerated CMOS switches are a really easy choice. The DG508 can be a simple way that even includes some protection.  JFET switching can get tricky to get the drive levels - one may need some kind of guard amplifier to get the right level to drive the gates.  With already a scanner card in the 3457 one may not need an extra scanner.

[rhb]

R-on of ~175 ohms  seems to make the DG508B unattractive for measuring voltage references unless there's a way to null it out.   My 4.5 digit Tenma DMM reads 0.50 uA  connected in series with the leads from an LM399 to my 3478A which shows a  300 uV drop  from the shunt reisstor in the Tenma.

The tempco of the DG508B from the datasheet is about  0.76 ohms/C which translates to 0.38 uV/C and a nominal offset voltage drop of 85 uV feeding the 3478A.  Of greater concern was the lack of any indication of aging performance,

I've been chasing down EMI in my work area much of the day, so my brain is elsewhere.  Do you know of a circuit for using the DG508/9A for ppm measurements?  The relay behavior that @zhtoor pointed to seems to me easier to model and correct.

My big interest in metrology at the moment is sparse L1 pursuits.  And I have the continual annoyance of the UIs on my Instek scopes as motivation to get back to learning to make a Zynq dance to my tune.  I'm getting very tempted to get a 2nd 3478A and build a 12 input relay scanner if someone will help out with the scanner PCB.  I feel very confident I can measure and correct the thermal EMF shown in the thread @zhtoor cited.

I want to *do* something, not buy something.  I invested a lot of time learning about sparse L1 pursuits.  I'd like to use it to do something original.  Personal conflicts with my supervisor at Austin led to my not getting my doctorate. I don't care about the certificate, but I am unhappy about not being allowed to complete the project.  It was good work, but was done a few years later by someone else.

I hope no one will take umbrage at my mentioning this here,  but I connected a $2.99 HF DMM to one of the LM399s today and trimmed it to match the 3478A.  Amazingly, it was within 0.5% when I tested it.  Rather  fiddly to adjust, but I got it to match the first 3 digits of the 3478A. .  I have no idea how old it is other than they have sold 2 newer models since I bought this one.

[cellularmitosis]

Hey Reg,

I would definitely encourage the DIY scanner route over the $500 used HP gear route.  $500 buys a lot of latching relays, and I think doktor pyta just recently reported that his DIY design had a thermal error of something like 10nV.  Pretty amazing for a budget of only a few bucks per channel.

(I'm actually working on a little DIY LM399 scanner myself, but it might be a few weeks before I have time to work out a board design)

If you were constrained by not being able to design a PCB, you could pick up one of the cards which goes into those scanners, and use an Arduino to drive it directly.  For example, searching for "03497-66509" (on ebay) shows a number of listings for under $50.

Either way, a fun little project for sure!

edit: link to doktor pyta's result: https://www.eevblog.com/forum/metrology/measuring-16-signals-with-one-multimeter/msg1463391/#msg1463391

edit2: looks like he is using G6KU-2F-Y relays, which are under $5 https://www.digikey.com/product-detail/en/omron-electronics-inc-emc-div/G6KU-2F-Y-TR-DC5/Z3706CT-ND/3908609

[rhb]

Thanks for the tip on the boards.  They're reasonably priced and the critical analog part is done to the best professional standards.  I just agreed to buy an HP 8569A.  So I think it's probably time to slow down buying equipment

I'd still like to try correcting the thermal errors in non-latchng  relays.

[rhb]

This is my candidate  system model:

initial voltage
1/f noise
thermal Johnson noise
thermal errors
aging drift

This is my candidate equation:

Vref(f) = V_initial + 1/f_noise(f) + Johnson_noise(f) + thermal_errors(f) + aging_drift(f) + mains_noise(f)

The 1/f noise has curvature at the LF end of the spectrum and the Johnson noise should be flat until it reaches the anti-alias filter Fc,.  The thermal errors  correlate with temperature.  The aging drift has curvature early in the series which flattens out over time, so I'll have to window the data into segments to derive that term.  The mains noise is narrowband, so it should be easy to model.

The early aging data is critical to accurately modeling aging.  Once the drift curve starts to flatten out it becomes increasingly difficult to separate from the 1/f noise.

PLC != 0 applies a sinc(f) filter to the data.  This has the consequence of making it difficult to separate the Johnson noise from the 1/f noise because the full amplitude frequency range is very narrow.  It remains to be seen how severe an impact it has on the separation.

[rhb]

After reading various descriptions of the 3458A I think it wise to skip over the multislope operation for the initial model.  My current model candidate is as follows:

1) Generate the noise component limited to one 1/f noise source and one thermal noise source

          compute Vnoise = Gaussian(Tmean, Tsigma) + A*Poisson/f in the frequency domain from DC to 3 KHz  to create a one hour record using     
         four independent Mersenne Twister PRNG sequences for the 1/f and thermal noise, one for the real part and one for the imaginary part.   
         The frequency sampling is 278 uHz.

2) inverse transform to time and add Vdc and mains noise

3) integrate for PLC=10 and subsample to 6 Samples/S

This basically assumes a "perfect" 3458A.  So all the noise is attributed to the reference which is assumed to have no drift over the course of one hour.

The Nyquist is 3 Hz, but the first notch is at 6 Hz, so everything above 3 Hz is aliased.  This seems to me the easiest way to reproduce the effect of the integration in the 3458A.

The experiment is to see how much data I need to recover Vdc to  < 1 ppm if I account for the noise explicitly in an L1 inverse problem.  This is not a sparse L1 pursuit, that implies multiple models for the 1/f noise and trying to find the one that fits best. Or more than one noise source with different characteristics.  It's the first step towards doing that and it seems that no one has done this. based on what I found searching IEEE Explore.

What I have in mind doing was not practical when the 3458A was designed.  I was in grad school in Austin at the time.  We could barely do L2 problems of this size on an 11/780.   I don't even want to think about how long a 22 million sample FFT would have taken.  Probably a couple of days.   Now it's a trivial problem.

[copied from separate thread which I should not have started]

This is the first step in testing the model in the previous post.  Once I have verified this can be solved I'll add the other errors to the model.

I had been concerned about aliasing, but as thermal noise is uncorrelated  one gets an interesting filter effect.  The higher order aliases cause the noise to increasingly cancel as the frequency approaches the Nyquist of 3 Hz.

I'll need to evaluate the expected amplitude of the thermal noise and the 1/f noise and study the cross correlation to get an idea of the limitations on separating Vdc from the noise.  The driving concept is that we can estimate the magnitude of the error at DC from the shape and amplitude of the noise components at higher frequencies and use that to correct the estimate of Vdc.

[RandallMcRee]

I created a java program to model voltage references a while ago.

https://github.com/RMcRee/ModelVoltageReference

It models pink noise and temperature drift. It does use Mersenne Twister, FYI.

Maybe it will be useful for you. The idea there was to recover individual reference behaviour by measuring the differences between each reference to one another. Presuming a known starting point, of course. Used LSQ to impute the actual voltages.

Randall

[rhb]

What were your results?  Did you take it any farther?  Tell me more.

[RandallMcRee]

So, yes, if all that happens to your vrefs is related to temperature you can track a set of references together.

What discourages me is that there is so much unaccounted for: long term drift, popcorn noise, thermal lag. This is not to mention the corresponding problems in the measurement system itself. Although this is presumably mitigated somewhat by measuring voltage differences.

The problem with my setup is that, once lost, the LSQ cannot recover the actual state of the whole system without redo-ing the initial (assumed) precise, calibrated measurement. That defeats the purpose since, if you had the ability to calibrate at any time, why bother with the LSQ ensemble? LSQ always comes up with a solution and you have no way of knowing if things have gone awry due to one of the unknown, not-accounted-for factors above.

All in all, I think it may be possible for you to achieve some of your goals with your more sophisticated techniques. I can't say my experiment proves it one way or another. But perhaps it does give you some things to think about.

Sort of on a side note. Beginners might wonder why the industry produces so many LT6655 type voltage references that have comparatively low initial accuracy and 2ppm/C temp change. Most ADCs for example have built-in references with those characteristics or worse. The answer is that they are suitable for ratiometric measurements. For metrology the bar is higher. Just thought I would throw that in since I have not seen this written down anywhere. Now, I'm a beginner, too, so if any of this is wrong or misleading, wait a few hours for the corrections to pour in!

[rhb]

Thanks.  I don't have a solution for the lack of an initial reference.  I'm afraid that is essential.  All the rest I *think* can be handled.  That's why I have to try this out.

I've been working on the effect of thermal (aka Johnson) noise in an integrating ADC.  I'm solving for the amplitude of the aliased noise after integration over a fixed window.  As it is Gaussian, zero mean noise, summing the aliases cancels by 1/sqrt(n).  So it seems to me that the amplitude of the thermal noise after aliasing and integration is:

J(f) = sinc(f) - sinc(f-pi)
J(f) += sqrt(2)*(sinc(f-pi) - sinc(f+pi))
J(f) += sqrt(3)*(sinc(f+pi) - sinc(f-2*pi))
J(f) += sqrt(4)*(sinc(f-2*pi)-sinc(f+2*pi))
.
.
.


so you have the sum of one noise process, two processes, three, etc. In addition to the noise level of each alias getting smaller, the difference of the terms is getting smaller.

the model is a sinc(f) convoived with a spike series in f. The terms converge rapidly to zero as j increases. So at any frequency you have all of the positive and negative frequency aliases summed together.  I've never dealt with anything like this before.  I'd be grateful for input from someone else with a DSP background.  I don't have a Stanford PhD across the hall to talk to anymore. The wall only responds when I walk into it.  Talking to it does no good and walking into it hurts.

I've attached a plot of the alias envelope amplitudes.  The aliases of concern  are to the left of Pi/2.

[edit; typo in plot.  fixed expression ]

[rhb]

Yesterday I couldn't figure out an algorithm for calculating the amplitude of thermal noise from an integrating ADC.  This morning it was easy.

It turns out I was quite wrong about convergence.  It is very slow.

The plots have the frequency normalized to the sample rate. So 0.5 is Nyquist.

Whereas 1/f noise decreases with frequency, the thermal noise increases quite a bit.  Much more than I expected.  The opposite slopes  bode well for separating the thermal and 1/f noise from the estimated DC value.  The statistical meaning of Vdc after accounting for the noise is yet to be determined.

The plots should be self explanatory.  Nyquist is in the middle.  Not quite sure what to call the right hand half from 0.5 to 1.0.  That is not in the spectrum of the data.

C code is attached if any want to play.

[rhb]

In light of  @RandallMcRee's efforts, I think I should note that what I'm doing right now can just as well be done with least squares. The Marquardt-Levenberg implementation in gnuplot is particularly good. A solver has to be able to beat that for me to consider using it.  However, gnuplot is configured to do function approximations rather than as a general Ax=y solver, so it is more convenient to use the machinery I already have in place for L1 pursuits.  Because the thermal noise has Gaussian properties, L2  will fit it just as well as L1.

The sparse L1 pursuit doesn't offer a significant  advantage until you need to select one of a large number of choices.  That presents an underdetermined problem which L2 cannot handle.  If the residue after removal of the thermal noise shows significant deviation from a simple 1/f relationship, then a sparse pursuit becomes important.  But I'll need a bunch of data to resolve that.

At the moment I am just building and testing the machinery.

The takeaway from this post is L1 is *really* good at some things L2 can't do.  But for a *lot* of things L2 is just as good and several orders of magnitude cheaper to compute.

[rhb]

How many angels can dance on the head of a pin?

The thermal noise amplitude spectrum shows the standard deviation of the thermal noise as a function of frequency after accounting for aliasing.  The other two plots indicate how high a frequency bleeds through the ADC and some notion of the numerical error of the modeling.  Actually quantifying the numerical error requires finding a solution to the series expansion.

However, you can't do anything about it.  It's just a precise statement of the uncertainty.  The same applies to the 1/f spectrum.  You can't correct for it except by collecting a large number of measurements and invoking the central limit theorem.

The sole benefit I can see at the moment is that having rigorous statements of the expected variance (i.e. the power spectrum) makes it possible to identify the presence of other noise sources in the measurement. Whether that is of significant practical benefit remains to be seen.

The ultimate goal of this project is to be able to accurately predict the aging curve for a voltage reference from data collected during the first 1000 or so hours of operation.  Properly quantifying the uncertainty created by the irreducible noise is a small, but essential step in that direction.

[zhtoor]

Hello @rhb,

find attached a sample dataset of 10,000 samples, sampling frequency is 3.89Hz, data starts at line #10.
a scatter plot is also attached.

questions:
1. how well (and if?) can the next 1000 samples be predicted from the previous 1000?
2. this scatter plot does look sparse to me, how can we further "sparsify" it?
3. how would you specify the L1 pursuit problem for this dataset?

ref:
further details available at: https://www.eevblog.com/forum/metrology/vlfn-characterization-of-reference-devices-topologies/msg1534340/#msg1534340

best regards.

-zia

[rhb]

The series is not long enough.  An eyeball inspection indicates there is no significant trend which can be predicted.

There might be a small time dependent trend, but the quantization obscures it.

The first order test to apply would be to examine the autocorrelation aka power spectrum and compare that to a spike convolved with the window function imposed by the series lenth.