While working in Athena today I noticed something peculiar while looking at the derivative plots of mu(E) data. When I plot the current group (red button) the plot can be spiky. However, when I plot as marked group (purple button) everything is much smoother and more rounded. Smoothed data will do the same thing. It always happens whether I have mu(E) or norm(E) selected. I have never noticed a difference in how any of the other plots current when switching between marked and active plotting. I didn't see anything in the documentation about how the derivatives are calculated so I'm turning to all of you. Have any of you noticed this? Does anyone know if there is a difference in how the derivatives are calculated in each plotting function? Andy Korinda Notestein Lab Northwestern University
On Tuesday, February 22, 2011 05:03:02 pm Andrew Korinda wrote:
While working in Athena today I noticed something peculiar while looking at the derivative plots of mu(E) data. When I plot the current group (red button) the plot can be spiky. However, when I plot as marked group (purple button) everything is much smoother and more rounded. Smoothed data will do the same thing. It always happens whether I have mu(E) or norm(E) selected. I have never noticed a difference in how any of the other plots current when switching between marked and active plotting. I didn't see anything in the documentation about how the derivatives are calculated so I'm turning to all of you. Have any of you noticed this? Does anyone know if there is a difference in how the derivatives are calculated in each plotting function?
Hi Andy, You have uncovered a very interesting mistake that I have never noticed before. First some background info. There is a preference that defined how many smoothing iterations are made when plotting a derivative spectrum. If you edit the preferences (Settings menu -> Edit preference), then click down to the plot -> smoothderiv parameter, you will see that it is set to 3. (Well, by default it should be, it could be different on your computer.) The plot->smoothderiv parameter defines the munber of smoothing iterations before plotting a derivative spectrum. If you set this to 0, then your data will not be smoothed. The default value of three means to run your data through the smoothing filter 3 times before plotting. When plotting with the purple button, the derivative spectrum is computed then smoothed. When plotting with the red button the mu(E) spectrum is smoothed then the derivative is computed. In an example that I just looked at, I can see a difference similar to what you describe due to the order of operation. The behavior of the red button does not seem right to me. I don't recall a reason why I made these behave differently. It seems to be a mistake that no one has noticed and reported to me until now. Good eye! I'll put it on my list of things to fix. B P.S. The order of operations is clearly documented in the Ifeffit buffer, which you can find in the Edit menu. The numerical derivative is computed as the numerical differential of the mu(E) spectrum divided by the numerical differential of the energy array. Smoothing happens as explained above. Justw atch the lines that get printed to the buffer when you click the plot buttons. -- Bruce Ravel ------------------------------------ bravel@bnl.gov National Institute of Standards and Technology Synchrotron Methods Group at NSLS --- Beamlines U7A, X24A, X23A2 Building 535A Upton NY, 11973 My homepage: http://xafs.org/BruceRavel EXAFS software: http://cars9.uchicago.edu/~ravel/software/exafs/
The numerical derivative is computed as the numerical differential of the mu(E) spectrum divided by the numerical differential of the energy array. Smoothing happens as explained above. Justw atch the lines that get printed to the buffer when you click the plot buttons. - B. Ravel
I wonder if that's the right way to go for data that, for instance, might have gaps or changes in step size. What I do in my code is run a cubic spline through the data, with the abscissa as is, then take the derivative of that spline. Thus, if the data is smooth enough to be locally approximated by a cubic, the derivative will come out smooth regardless of how it's tabulated. I'm not sure that's the case with Bruce's algorithm. Is the 'numerical differential' the difference E(i+1)-E(i) or (1/2)*(E(i+1)-E(i-1)) or some form that's intended to approximate dE(i)/di (i= point index)? It's not clear to me whether it's any more correct to do smoothing before or after differentiating. You can take the red pill or the blue pill :-) mam
On Tuesday, February 22, 2011 06:19:16 pm Matthew Marcus wrote:
The numerical derivative is computed as the numerical differential of the mu(E) spectrum divided by the numerical differential of the energy array. Smoothing happens as explained above. Justw atch the lines that get printed to the buffer when you click the plot buttons. - B. Ravel
I wonder if that's the right way to go for data that, for instance, might have gaps or changes in step size. What I do in my code is run a cubic spline through the data, with the abscissa as is, then take the derivative of that spline. Thus, if the data is smooth enough to be locally approximated by a cubic, the derivative will come out smooth regardless of how it's tabulated. I'm not sure that's the case with Bruce's algorithm. Is the 'numerical differential' the difference E(i+1)-E(i) or (1/2)*(E(i+1)-E(i-1)) or some form that's intended to approximate dE(i)/di (i= point index)?
It's not clear to me whether it's any more correct to do smoothing before or after differentiating. You can take the red pill or the blue pill :-) mam
Well ... red pill or purple pill. They both have side effects that 4 out of 5 doctors advise against :) So, I do what's easy using Ifeffit as Ifeffit is implemented. I would have to code dive (or ask Matt) what the deriv() function actually does. A big gap in the data might result in a weird derivative spectrum, but changes in step size isn't so much of a problem because I ask Ifeffit to compute d(mu)/d(E). I agree that "smooth then differentiate" is not obviously better than "differentiate then smooth". What is certainly wrong is to do this inconsistently. B -- Bruce Ravel ------------------------------------ bravel@bnl.gov National Institute of Standards and Technology Synchrotron Methods Group at NSLS --- Beamlines U7A, X24A, X23A2 Building 535A Upton NY, 11973 My homepage: http://xafs.org/BruceRavel EXAFS software: http://cars9.uchicago.edu/~ravel/software/exafs/
Hi Bruce, Matthew,
Oh, deriv() (and smooth()) are not nearly as sophisticated as what
Matthew suggests. They are simply applied to a single array,
sequentially, without regard to "the x axis" (as, in fact, it does
not know what x might be):
Deriv would translate (sorry for the python) to
def deriv(x):
npts = len(x)
out = array(npts)
out[0] = x[1] - x[0]
out[npts-1] = x[npts-1] - x[npts-2]
for i in range(1, npts-2): # (1, 2, 3, ... npts-2)
out[i] = (x[i+1] - x[i])/2.0
return out
and Smooth:
def smooth(x):
npts = len(x)
out = array(npts)
out[0] = 0.75*x[0] + 0.25 * x[1]
out[npts-1] = 0.75*x[npts-1] + 0.25 * x[npts-2]
for i in range(1, npts-2): # (1, 2, 3, ... npts-2)
out[i] = 0.5*x[i] + 0.25*x[i+1] + 0.25*x[i-1]
return out
Like Matthew suggests, interpolating onto a finer, uniform x (energy)
grid and then taking the derivative might be a better approach, at
least for data which needs smoothing.
--Matt
On Wed, Feb 23, 2011 at 8:46 AM, Bruce Ravel
On Tuesday, February 22, 2011 06:19:16 pm Matthew Marcus wrote:
The numerical derivative is computed as the numerical differential of the mu(E) spectrum divided by the numerical differential of the energy array. ÃÂ Smoothing happens as explained above. ÃÂ Justw atch the lines that get printed to the buffer when you click the plot buttons. - B. Ravel
I wonder if that's the right way to go for data that, for instance, might have gaps or changes in step size. ÃÂ What I do in my code is run a cubic spline through the data, with the abscissa as is, then take the derivative of that spline. ÃÂ Thus, if the data is smooth enough to be locally approximated by a cubic, the derivative will come out smooth regardless of how it's tabulated. ÃÂ I'm not sure that's the case with Bruce's algorithm. ÃÂ Is the 'numerical differential' the difference E(i+1)-E(i) or (1/2)*(E(i+1)-E(i-1)) or some form that's intended to approximate dE(i)/di (i= point index)?
It's not clear to me whether it's any more correct to do smoothing before or after differentiating. ÃÂ You can take the red pill or the blue pill :-) ÃÂ ÃÂ mam
Well ... red pill or purple pill. ÃÂ They both have side effects that 4 out of 5 doctors advise against :)
So, I do what's easy using Ifeffit as Ifeffit is implemented. ÃÂ I would have to code dive (or ask Matt) what the deriv() function actually does.
A big gap in the data might result in a weird derivative spectrum, but changes in step size isn't so much of a problem because I ask Ifeffit to compute d(mu)/d(E).
I agree that "smooth then differentiate" is not obviously better than "differentiate then smooth". ÃÂ What is certainly wrong is to do this inconsistently.
B
--
ÃÂ Bruce Ravel ÃÂ ------------------------------------ bravel@bnl.gov
ÃÂ National Institute of Standards and Technology ÃÂ Synchrotron Methods Group at NSLS --- Beamlines U7A, X24A, X23A2 ÃÂ Building 535A ÃÂ Upton NY, 11973
ÃÂ My homepage: ÃÂ ÃÂ http://xafs.org/BruceRavel ÃÂ EXAFS software: http://cars9.uchicago.edu/~ravel/software/exafs/ _______________________________________________ Ifeffit mailing list Ifeffit@millenia.cars.aps.anl.gov http://millenia.cars.aps.anl.gov/mailman/listinfo/ifeffit
Should the line
out[i] = (x[i+1] - x[i])/2.0
be
out[i] = (x[i+1] - x[i-1])/2.0 ?
That would account for the 2.0. That is, in fact, exactly what I proposed was meant by 'numerical differential'.
I don't actually resample. I do a cubic spline interpolation onto the given grid and pull the derivatives out that way. Of course,
you have to test for
grids that aren't strictly monotonic. The smooth() function as given would definitely cause problems with non-uniform tabulation
unless resampled.
My program (2-column Editor on the ALS Beamline 10.3.2 website) offers resampling and tensioned spline fitting but does not smooth
the derivative
before or after - the user has to take control of that using the tensioned spline fit. I thought about smoothing but never got
around to figuring out what makes
a decent smooth algorithm for general data. A possibility would be a slow version of Savitzky-Golay in which you do the polynomial
fits explicitly instead of
being able to use pre-computed coefficients as in the actual S-G method. Similarly, a VSFT (Very Slow Fourier Transform) in which
you fit the data to
a sum of sines and cosines (linear fit - the frequencies are fixed), filter using your favorite low-pass filter, then
inverse-transform, might be good. Or, there's
always reading the literature to see what actual numerical mathematicians have done :-)
There are two obvious ways to get non-uniform tabulation in XAS data. One is that it's taken on a varying grid, either
constant-k-step or varying E step. The
other is that pieces have been chopped out due to glitches, Bragg peaks or some other problem. Also, the 2-column Editor program
was designed to be general,
not restricted to XAS-like data.
mam
----- Original Message -----
From: "Matt Newville"
On Tuesday, February 22, 2011 06:19:16 pm Matthew Marcus wrote:
The numerical derivative is computed as the numerical differential of the mu(E) spectrum divided by the numerical differential of the energy array. Ã, Smoothing happens as explained above. Ã, Justw atch the lines that get printed to the buffer when you click the plot buttons. - B. Ravel
I wonder if that's the right way to go for data that, for instance, might have gaps or changes in step size. Ã, What I do in my code is run a cubic spline through the data, with the abscissa as is, then take the derivative of that spline. Ã, Thus, if the data is smooth enough to be locally approximated by a cubic, the derivative will come out smooth regardless of how it's tabulated. Ã, I'm not sure that's the case with Bruce's algorithm. Ã, Is the 'numerical differential' the difference E(i+1)-E(i) or (1/2)*(E(i+1)-E(i-1)) or some form that's intended to approximate dE(i)/di (i= point index)?
It's not clear to me whether it's any more correct to do smoothing before or after differentiating. Ã, You can take the red pill or the blue pill :-) Ã, Ã, mam
Well ... red pill or purple pill. Ã, They both have side effects that 4 out of 5 doctors advise against :)
So, I do what's easy using Ifeffit as Ifeffit is implemented. Ã, I would have to code dive (or ask Matt) what the deriv() function actually does.
A big gap in the data might result in a weird derivative spectrum, but changes in step size isn't so much of a problem because I ask Ifeffit to compute d(mu)/d(E).
I agree that "smooth then differentiate" is not obviously better than "differentiate then smooth". Ã, What is certainly wrong is to do this inconsistently.
B
--
Ã, Bruce Ravel Ã, ------------------------------------ bravel@bnl.gov
Ã, National Institute of Standards and Technology Ã, Synchrotron Methods Group at NSLS --- Beamlines U7A, X24A, X23A2 Ã, Building 535A Ã, Upton NY, 11973
Ã, My homepage: Ã, Ã, http://xafs.org/BruceRavel Ã, EXAFS software: http://cars9.uchicago.edu/~ravel/software/exafs/ _______________________________________________ Ifeffit mailing list Ifeffit@millenia.cars.aps.anl.gov http://millenia.cars.aps.anl.gov/mailman/listinfo/ifeffit
_______________________________________________ Ifeffit mailing list Ifeffit@millenia.cars.aps.anl.gov http://millenia.cars.aps.anl.gov/mailman/listinfo/ifeffit
Hi Matthew,
On Tue, Feb 22, 2011 at 9:42 PM, Matthew Marcus
Should the line out[i] = (x[i+1] - x[i])/2.0 be out[i] = (x[i+1] - x[i-1])/2.0 ?
That would account for the 2.0. That is, in fact, exactly what I proposed was meant by 'numerical differential'.
Yes, sorry for the typo, it should b def deriv(x): npts = len(x) out = array(npts) out[0] = x[1] - x[0] out[npts-1] = x[npts-1] - x[npts-2] for i in range(1, npts-2): # (1, 2, 3, ... npts-2) out[i] = (x[i+1] - x[i-1])/2.0 return out
I don't actually resample. I do a cubic spline interpolation onto the given grid and pull the derivatives out that way. Of course, you have to test for grids that aren't strictly monotonic. The smooth() function as given would definitely cause problems with non-uniform tabulation unless resampled. My program (2-column Editor on the ALS Beamline 10.3.2 website) offers resampling and tensioned spline fitting but does not smooth the derivative before or after - the user has to take control of that using the tensioned spline fit. I thought about smoothing but never got around to figuring out what makes a decent smooth algorithm for general data. A possibility would be a slow version of Savitzky-Golay in which you do the polynomial fits explicitly instead of being able to use pre-computed coefficients as in the actual S-G method. Similarly, a VSFT (Very Slow Fourier Transform) in which you fit the data to a sum of sines and cosines (linear fit - the frequencies are fixed), filter using your favorite low-pass filter, then inverse-transform, might be good. Or, there's always reading the literature to see what actual numerical mathematicians have done :-)
Smooth() is definitely over-simplistic and will work poorly for non-uniform grids. Using a Savitzky-Golay smoothing function is probably a good idea. Having the abilty to use more sophisticated numerical algorithms on demand would be very nice.... --Matt
participants (4)
-
Andrew Korinda -
Bruce Ravel -
Matt Newville -
Matthew Marcus