steve@nuchat.sccsi.com (Steve Nuchia) (02/14/91)
-2) This (comp.dsp) seems to be the best place. Redirection welcome. sci.math added as an afterthought just in case somebody there finds this interesting or happens to know the answer. -1) I am not a regular reader of comp.dsp (one can only read so many things) so I would appreciate mailed responses. I do subscribe to sci.math, but haven't been able to keep up recently. 0) I'm doing consulting for a corporate client -- I don't expect a lot of free help here. Pointers to accessible literature will be gratefully followed up. 1) I'm a senior EE student with 7 years professional software (ahem) engineering and consulting experience. I've had basicly only the cartoon version of digital filtering: the second semester linear systems course covering discrete time Fourier theory and such. 2) I'm trying to analyse the following situation: A package of measureing instruments (we can assume the singular without loss of generality) traverses some physical environment at an unregulated speed, with its position indirectly measured. It is traditional in this industry to ignore the slop in the mechanical system and use the indicated position for downstream processing directly. Let us accept this simplification for this discussion, as it is orthogonal to the main point. The instruments are sampled at intervals. Some instruments are designed for fixed-frequency sampling, others are sampled on a temporal basis for historical reasons, and still others are sampled on a very roughly spacial basis (first opportunity following floor(pos/interval) changes, usually). The underlying phenomena are assumed to be spacially band-limited. Output samples of an estimate of the underlying phenomena are required at regular spacial intervals, in real time (but a constant lag is OK -- phase linear filters and all that) I want to understand the theory well enough to make a sound recomendation on the design of digital filters (to be implemented on a general purpose CPU in real time) to compute the estimate samples from the available sample data. 3) Is this a well-studied problem? If so, a reference should be all I need. 4) If not, I would appreciate any advice on how to procede. The following no-doubt-amateurish remarks are my thoughts to date. If they are not completely bogus, perhaps someone could give some advice on how to tighten it up? I have to admit to posting this before taking the math as far as I could, but I'm trying to avoid having to do all those transforms :-( Anyway, I'm thinking of proceeding like this: Normalize all cases to depth-referenced samples on a fine (high spacial frequency) comb. The synthetic sampling spacial frequency would be selected based on the output spacial frequency -- say 2x? 1x? 4x? Actual samples will be scattered randomly about in any given filtering window. Let us model this as multiplication of a hypothetical ideal sampling of the input at the synthetic frequency with a (Poisson-distributed, I suppose. Not important yet) random unit impulse comb. If we consider the random comb (ie, the jitter filter) and the low-pass filter as acting on the ideal input we get a mess, since its DC gain, if nothing else, will be all over the map. By convolving the jitter comb itself with the LPF function we can get a number measuring (at least in some sense) that mess. In particular, that number will be the DC gain of the composed filter. There will be other distortions not capturable in a single scalar of course. Assuming that the jitter comb points are sufficiently dense and uniform it would seem that normalizing the output wrt the DC gain estimate would be a reasonable thing to do. An alternate approach would be to do an LPF on the temporal axis and use interpolated samples of its output as spacial samples. This would seem to be essentially the same analysis except that the properties (in the filtering domain) of the jitter comb would be much more stable (less dependent on the instrument velocity). The problem with this approach is that the *spacial* cutoff frequency of the filtering process would then depend on the velocity. Actually that last may be the desired behaviour, but let us assume it is not. (A colleague is working on the physics of both the phenomena and the transucer behaviour to answer that question.) One could dynamically adjust the implicit scale of the FIR function, withing limits, to at least partially compensate for changes in instrument velocity. Does cutoff frequency tweaking by stride modulation in a FIR lookup table work OK, or would it be better for instance to "autorange" among a set of filter tables? Or is there a good fast way of doing the convolution without using a lookup table? Am I even asking the right questions? Should I be thinking IIR? Should I just hook the plotter up to a canned demo data stream? [:-)] Using a sample-hold approach might seem simpler, and it is in fact essentially what the existing system is using (along with really gross causal filters). The problem is that if one wants to use higher order filters and properly account for jitter, the number of points in the convolution climbs rapidly, and we only have so much CPU to spend. The comb-compensation mechanism proposed above uses all samples, allows for high-order (low quantization noise) underlying filter functions, and runs in time proportional to the number of available samples in the window. 5) Finally, if this isn't well studied, is it worth writing up? :-) Anybody with a stronger math background want to co-author? Any Rice people on these lists want to straighten me out over some expense-account suds? Thanks, all. -- Steve Nuchia South Coast Computing Services (713) 964-2462 "Innocence is a splendid thing, only it has the misfortune not to keep very well and to be easily misled." --- Immanuel Kant, Groundwork of the Metaphysic of Morals
jbuck@galileo.berkeley.edu (Joe Buck) (02/15/91)
If you have available the exact x values where the samples are taken,
and your signal is band-limited, you can reconstruct exactly the
continuous-time signal. Let x[n] be the space coordinate of the nth
sample; let y[n] be the measured value of the signal at that time.
First, consider a function that is zero for all x except at the
x[n] values, and has a Dirac delta function of height y[n] at
x = x[n]. Pass this signal through an ideal low-pass filter corresponding
to the band limit, and you have the continuous-time signal. This
reconstruction works if the highest frequency in the data, f, is
less than 1/2T, where T is the maximum spacing between samples.
For this to work right at the boundary f = 1/2T you must have no
noise. I think you could tolerate noise if you oversampled by a
good deal, so the LPF would take it out.
You can then resample the continuous-time signal at your desired
rate.
I recommend that you simulate this procedure before actually relying
on it.
--
Joe Buck
jbuck@galileo.berkeley.edu {uunet,ucbvax}!galileo.berkeley.edu!jbuck steve@nuchat.sccsi.com (Steve Nuchia) (02/17/91)
In <11145@pasteur.Berkeley.EDU> jbuck@galileo.berkeley.edu (Joe Buck) writes: >If you have available the exact x values where the samples are taken, >and your signal is band-limited, you can reconstruct exactly the >continuous-time signal. Let x[n] be the space coordinate of the nth Uhm, no, I don't think so... We proved in class that the sampling theorem holds when the pulse spacing is irregular *but periodic*. I am working with aperiodic sample spacing in which the (average) sampling frequency is not constant. A simple counterexample serves to show that one cannot recover the underlying signal without compensating for sampling frequency changes: DC. Update: the name for the class of problems to which mine belongs seems to be "Multirate Filtering". I've found a book on the subject and will be wrapping my head around it. >You can then resample the continuous-time signal at your desired rate. Part 2 of my problem is to compute the samples directly, rather than go through an intermediate. Particularly not an analog intermediate. I suppose that may be obvious, but the quoted line makes me fear I failed to make myself clear. I also misspelled "spatial". Sigh. -- Steve Nuchia South Coast Computing Services (713) 964-2462 "Innocence is a splendid thing, only it has the misfortune not to keep very well and to be easily misled." --- Immanuel Kant, Groundwork of the Metaphysic of Morals
jfa0522@hertz.njit.edu (john f andrews ece) (02/18/91)
In article <1991Feb16.160801.7117@nuchat.sccsi.com> steve@nuchat.sccsi.com (Steve Nuchia) writes: >In <11145@pasteur.Berkeley.EDU> jbuck@galileo.berkeley.edu (Joe Buck) writes: >>If you have available the exact x values where the samples are taken, ...stuff deleted >>You can then resample the continuous-time signal at your desired rate. > >Part 2 of my problem is to compute the samples directly, rather than >go through an intermediate. Particularly not an analog intermediate. Just for clarity I would like to point out that resampling does not require an analog intermediate but simply resampling of the resulting digital output... BTW, can anyone recall the discussion on applying similar interpolation to a sampled EKG some time ago? I believe someone wanted to reconstruct a complete EKG waveform from sampled data, and was actually testing these digital techniques and reporting results here? Perhaps I am incorrect but I beleive many of these ideas were tossed about then ... ----------------------------------------------------------------------------- john f andrews SYSOP The Biomedical Engineering BBS 24 hrs 300/1200/2400 (201) 596-5679 ----------------------------------------------------------------------------- INTERNET jfa0522@hertz.njit.edu LabRat@faraday.njit.edu CIS 73710,2600 -----------------------------------------------------------------------------
jcs@math.ufl.edu (Joshua C Sasmor) (02/19/91)
In article <2341@njitgw.njit.edu> jfa0522@hertz.njit.edu (john f andrews ece) writes: >In some article steve@nuchat.sccsi.com (Steve Nuchia) writes: >BTW, can anyone recall the discussion on applying similar interpolation >to a sampled EKG some time ago? I believe someone wanted to reconstruct a >complete EKG waveform from sampled data, and was actually testing these >digital techniques and reporting results here? Perhaps I am incorrect but >I beleive many of these ideas were tossed about then ... This past summer I was working as a research engineer's assistant at a pacemaker company (which shall remain nameless) and we were converting a digitized EKG from one sample rate to another. The original was sampled at 1000 samples/sec and the new digital scope needed to have it's rates as powers of 2. So we did a direct linear interpolation on the 1000 sample/sec tape and reconstructed an EKG at 1024 samples/sec. It was a trivial programming job once the algorithm was worked out properly. Four people with engineering PhDs came to me for math help :) (I will graduate with a BS in math in may.) There was only one unusual result. We converted 1000 to 1024 and then back to 1000 and did an error check. At the nth power of 2 there were approximately 2n terms in error, evenly distributed around the power of 2 itself. This was not a particularly large error (10^-16) but the error doubled with every power as well. Perhaps someone here could explain that one for me. The engineers simply said it was negligible and ignored it. ******************************************************************************* Joshua C. Sasmor -- famous mathematician-to-be * jcs@lab4.math.ufl.edu ************************************************* brassman@maple.circa.ufl.edu "The highest form of pure thought is in * x9999bis@maple.circa.ufl.edu mathematics" -- Plato * *******************************************************************************
mcmahan@netcom.COM (Dave Mc Mahan) (02/20/91)
In a previous article, jfa0522@hertz.njit.edu (john f andrews ece) writes: >In article <1991Feb16.160801.7117@nuchat.sccsi.com> steve@nuchat.sccsi.com (Steve Nuchia) writes: >>Part 2 of my problem is to compute the samples directly, rather than >>go through an intermediate. Particularly not an analog intermediate. > >Just for clarity I would like to point out that resampling does not require >an analog intermediate but simply resampling of the resulting digital >output... > >BTW, can anyone recall the discussion on applying similar interpolation >to a sampled EKG some time ago? I believe someone wanted to reconstruct a >complete EKG waveform from sampled data, and was actually testing these >digital techniques and reporting results here? Perhaps I am incorrect but >I beleive many of these ideas were tossed about then ... I believe you were talking about me. I did manage to get some FINE looking waveforms from my re-constituted data. I haven't been closely following this thread. I'll back track a bit and catch up, then I'll give you my 2 cents. The algorithm runs quite quickly and works strictly in the time domain. I can process one data point every 1.5 milliseconds on a 68000 running at 10 MHz and can generate 4 points (including the original point) to feed to the next stage of my process. With a bit of tweaking, I think I can get the CPU time down to about 1 millisecond per input data point. >john f andrews SYSOP The Biomedical Engineering BBS >INTERNET jfa0522@hertz.njit.edu LabRat@faraday.njit.edu CIS 73710,2600 -dave -- Dave McMahan mcmahan@netcom.com {apple,amdahl,claris}!netcom!mcmahan
mcmahan@netcom.COM (Dave Mc Mahan) (02/23/91)
In a previous article, steve@nuchat.sccsi.com (Steve Nuchia) writes: >2) I'm trying to analyse the following situation: > > A package of measureing instruments (we can assume the > singular without loss of generality) traverses some > physical environment at an unregulated speed, with > its position indirectly measured. > > The instruments are sampled at intervals. Some instruments > are designed for fixed-frequency sampling, others are sampled > on a temporal basis for historical reasons, and still others > are sampled on a very roughly spacial basis (first opportunity > following floor(pos/interval) changes, usually). > > The underlying phenomena are assumed to be spacially band-limited. > > Output samples of an estimate of the underlying phenomena are > required at regular spacial intervals, in real time (but a > constant lag is OK -- phase linear filters and all that) > > I want to understand the theory well enough to make a sound > recomendation on the design of digital filters (to be implemented > on a general purpose CPU in real time) to compute the estimate samples > from the available sample data. > >Anyway, I'm thinking of proceeding like this: Normalize all cases >to depth-referenced samples on a fine (high spacial frequency) comb. >The synthetic sampling spacial frequency would be selected based on >the output spacial frequency -- say 2x? 1x? 4x? First off, I assume that you KNOW the time at which samples were taken with a 'large' (relative to the overall system) degree of accuracy. Your problem is that you just can't control WHEN the samples are taken. The data you take (in terms of the measured position) is also of sufficient accuracy that you don't have to worry about positional uncertainty. For the sake of argument you 'know' with infinite precision both the time at which a sample was taken and the position of the sample. Next, I assume you wish to synthesize some time/position related output data of what happened between samples. This may mean that you wish to guess-timate the position at any time between two samples, or you may wish to calculate velocity at some point in time (or space) of the thing begin sampled. Assuming that you can guess-timate the position at any desired time, I assume you can extrapolate any other desired information such as velocity or acceleration. I also assume that if some of your measurements are taken as velocity, others are acceleration measurements, and still others are displacement measurements, you can convert (with enough math) to a consistent base for the problem. From here on, I will assume you can convert to or have taken all measurements in terms of displacment. Finally, we shall assume that you have enough time between measurements and enough CPU horsepower to do a reasonable amount of calculations and are therefore able to obtain the desired answers quickly enough to be useful. I'm not sure how many measurements per second you are taking, but I assume it is not more than about 500 or 1000 samples per second, probably much less than this. For the sake of argument, I'm going to assume that you have something on the order of a 80386 (or larger) available to do your desired math, and all the I/O hardware and instrumentation required to take the samples, timestamp them, feed them to the CPU, and then distribute the results in whatever form you need. I also assume that you don't need answers with 16 digits of precision and that your data doesn't span several orders of magnitude in range. If any of my assumptions or problem re-statements above are wrong, please let me know. I'm having a bit of trouble trying to visualize just exactly what you are trying to accomplish. I think I understand, but I am not certain. I would say that it is quite possible for you to get good results with your synthetic sampling interval approach. It strikes me that the spacing between points that you wish to synthesize should be as small as the shortest time between sample points. Since this approach will lead to most of your synthetic samples (the ones you create to just 'fill' sample points you don't know between the ones you do know) being zero, you will be able to write a smart/fast enough program to ignore multiplication by zero and save lots of CPU time. One problem with this is that if your samples are more-or-less evenly spaced but are not exact multiples of each other, you could end up generating lots of synthetic samples that will later need to be discarded. An example of this is if you have samples that are sometimes .9 seconds apart, sometimes 1 second apart, and sometimes 1.05 samples apart, you would either need to wave your hands and just assume that all are 1 second apart (generating the least synthetic samples but losing the most amount of information), assume that all samples should occur at .1 second intervals (thus saving more information on about your data but forcing you to make the decision as to whether 1.05 second samples get binned up to 1.1 seconds or down to 1.0 seconds) or assume that all samples should occur at .05 second intervals (thus extracting the most information from your data at the expense of more CPU time required to crunch numbers). Your synthetic sample interval period will also be dictated by your output data requirement. In the above example, what should you use if you need to know the position estimate at every .025 seconds? I guess the answer is, "It depends." (It depends on your other system requirements, which I don't really know). In this case, you can probably get by with just using the standard 'sync' function method of reconstructing data between known sample points. This method is outlined in the book, "Digital Signal Processing" by Oppenheim & Schafer on page 29. This is esentially creating a low pass filter and passing your data thru it by convolution. This can be done in the time domain and is fairly straightforward. You need to pick a finite number of filter points that will give you good results. Without knowing the nature of your data, I haven't a clue as to how many that is. I would suggest if you do pick this approach that you use an appropriate windowing function to smooth the sync filter function before you pass your data thru it. I have found that with a finite number of points in the sync function, the output results look MUCH better when the sync function has been windowed. One other approach you may wish to think about instead of using filtering is that of curve fitting. It can be done quite quickly and is pretty good and interpolation if your data can be made to fit onto a smoothly varying curve. Curve fitting is somewhat of an art, and there are many tricks to getting good results. I know one or two just by having observed a few masters at the art, but it generally boils down to knowing something about the nature of your data. If (for example) you can expect all of your data to fit a logrithmic curve, there are tricks to play that will give you optimal results for that type of curve. This works quite well with things like temperature sensors that exhibit general responses that fit on a log curve but need to be corrected based on a few know sample points taken at calibration time. The tricks you pick all depend on your expectation of what your data will look like. Anyway, curve fitting to interpolate between data points may be what you need, but I can't say for sure without more info. >One could dynamically adjust the implicit scale of the >FIR function, withing limits, to at least partially compensate for >changes in instrument velocity. You lost me on this one. Are you saying that you are trying to take measurements that aren't time-invariant? This could be big problems. It seems to me that if your instrument passes point A at one velocity and you take a measurement, you better have the ability to pass point A at a different velocity later and measure the same position. If you don't then you will need to record the velocity (or temperature, or direction, or whatever affects your position measuring) and filter that out before you ever begin. >Am I even asking the right questions? Should I be thinking IIR? >Should I just hook the plotter up to a canned demo data stream? [:-)] You should definately test drive any approach you come up with against a known (and repeatable!) data set so you can test different approaches. Having some way to play back the same data into different systems and comparing it with known-good results always gives you ideas of how you can make your system better. I'm sure I have missed a few (or more) points you think I should know. Please post to comp.dsp or e-mail me directly and I'll give you my input. >-- >Steve Nuchia South Coast Computing Services (713) 964-2462 -dave -- Dave McMahan mcmahan@netcom.com {apple,amdahl,claris}!netcom!mcmahan
gah@hood.hood.caltech.edu (Glen Herrmannsfeldt) (02/28/91)
The Nyquist theorem does not require equal spaced points. Equal spacing is convenient, and good in terms of S/N ratio. If your samples are infinitely narrow (points in time), then, theory says that if you multiply each point by a dirac delta function centered at that time, and filter it through the appropriate band limiting filter, you will get the same result as the original signal through such band limiting filter. If the points are not equally spaced, then errors in the position of those points are magnified in the final result. If you have velocity or acceleration, you use first or second derivative of the dirac delta function, I believe. (I have never actually tried this.) As a side thought, if you know all the derivatives of a function at a single point in time, you can find the value of the function for all times. The problem is you have to know them exactly.
jumper@spf.trw.com (Greg Jumper) (03/01/91)
> As a side thought, if you know all the derivatives of a function at a > single point in time, you can find the value of the function for all times. This statement is only true if the function is so-called "real-analytic," which is more restrictive than even "smooth" (C-infinity). Admittedly, functions which are not real-analytic are "pathological," particularly in a "real-world" setting. Ironically (since the subject is sampling theory), it turns out that the examples of functions whose Taylor series representations do not converge, even though all their derivatives exist everywhere, are constructed using Fourier methods -- essentially by taking advantage of "Gibb's phenomenon" to produce non-convergence. (Even "almost everywhere", if I remember correctly.) An interesting counter-example! Greg Jumper TRW Data Systems Center jumper@spf.trw.com
deghare@daisy.waterloo.edu (Dave Hare) (03/16/91)
In article <27CD919E.5D23@deneva.sdd.trw.com> jumper@spf.trw.com (Greg Jumper) writes: [The attribution for the original posting seems to have been lost in the follow-up-- DH] >> As a side thought, if you know all the derivatives of a function at a >> single point in time, you can find the value of the function for all times. > >This statement is only true if the function is so-called "real-analytic," >which is more restrictive than even "smooth" (C-infinity). Admittedly, >functions which are not real-analytic are "pathological," particularly in a >"real-world" setting. > >Ironically (since the subject is sampling theory), it turns out that the >examples of functions whose Taylor series representations do not converge, >even though all their derivatives exist everywhere, are constructed using >Fourier methods -- essentially by taking advantage of "Gibb's phenomenon" to >produce non-convergence. (Even "almost everywhere", if I remember correctly.) > >An interesting counter-example! As is / exp(-1/x^2) x <> 0 f(x) = < \ 0 x = 0 The original statement, even as qualified, is not precisely true, as it assumes an infinite radius of convergence for the Taylor series. *D*