wave-shape fidelity index

Share icon

Hi team of IGOR

I tried to calculate wave-shape fidelity index. It is a simple formula as shown below

\Z18\Z14\[0F =\Z20\F'Symbol'\[1S\F]0\X1\Z20\B\B  i\M\X1\Z20\S\S\Z09 N\M f((Wb-Wt)\S2\M/n)


Where Wb is the wave of baseline group and Wt is the wave of treatment group. I have successfully generated the algorithm as shown below
///Wave-shape fidelity index
Function WSFI(baselinewave, treatmentwave)
Wave baselinewave, treatmentwave
duplicate baselinewave resultwave
wave resultwave
wavestats/q resultwave
variable numberPoints = V_npnts
resultwave = (((baselinewave-treatmentwave)^2)/numberPoints)
variable WSFIvalue = sum(resultwave)
print WSFIvalue, resultwave
KillWaves/A/Z  
End


However, when I tried to introduce a selective range (bracket) base on X axis (left to right), I didn't get what I want.
I tried wavestats/R= (left, right) ...but failed to pick the range

Here is the code
///Wave-shape fidelity index
Function WSFI(baselinewave, treatmentwave, Xleft, Xright)
variable Xleft
variable Xright
Wave baselinewave, treatmentwave
duplicate baselinewave resultwave
wave resultwave
wavestats/q/R=(Xleft, Xright) resultwave
variable numberPoints = V_npnts
resultwave = (((baselinewave-treatmentwave)^2)/numberPoints)
variable WSFIvalue = sum(resultwave)
print WSFIvalue, resultwave
KillWaves/A/Z  
End

The selection failed to select the range but the numberpoint works. I don't know how to correct it. Any suggestions?
Thanks in advance.
Consider these two functions.



The code you have programmed is the first function.

What is n? Is it the index position in the wave? If so, you have to rethink your coding to give the second form.

Otherwise, you will want to select the x range in Wb and Wt with syntax such as Wb(xl,xr) or Wt(xl,xr).

--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAHuntsville
equation.png (11.84 KB)

January 27, 2014 at 06:22 pm - Permalink

wavestats/q/R=(Xleft, Xright) resultwave
variable numberPoints = V_npnts
resultwave = (((baselinewave-treatmentwave)^2)/numberPoints)


Note that wavestats/R has no effect on the later wave assignment, for which there is no selection. It only calculates statistic parameters for the range given by the /R flag. You probably want Duplicate/R, but I guess this also needs to be applied to baselinewave and treatmentwave.

January 28, 2014 at 12:03 am - Permalink

ChrLie wrote:
wavestats/q/R=(Xleft, Xright) resultwave
variable numberPoints = V_npnts
resultwave = (((baselinewave-treatmentwave)^2)/numberPoints)


Note that wavestats/R has no effect on the later wave assignment, for which there is no selection. It only calculates statistic parameters for the range given by the /R flag. You probably want Duplicate/R, but I guess this also needs to be applied to baselinewave and treatmentwave.


Thanks ChrLie, your suggestion worked. here is the code i modified based on your suggestion.
Function WSFI(baselinewave, treatmentwave,left, right)
variable left,right
Wave baselinewave
Wave treatmentwave
duplicate/r=(left,right) baselinewave baselinewaveSet
duplicate/r=(left,right) treatmentwave treatmentwaveSet
duplicate/r=(left,right) baselinewave resultwave
wave resultwave
wavestats/q/r=(left, right) resultwave
variable numberPoints = V_npnts
resultwave = (((baselinewaveSet-treatmentwaveSet)^2)/numberPoints)
variable WSFIvalue = sum(resultwave)
print WSFIvalue, resultwave
KillWaves/A/Z  
End

January 28, 2014 at 11:07 am - Permalink

jjweimer wrote:
Consider these two functions.



The code you have programmed is the first function.

What is n? Is it the index position in the wave? If so, you have to rethink your coding to give the second form.

Otherwise, you will want to select the x range in Wb and Wt with syntax such as Wb(xl,xr) or Wt(xl,xr).

--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAHuntsville


I see your point. Sorry that I did not make it clear. The n is the value of y-axis based on labmda. it is the former formula. The N is the total of points.
Thank you for the input.

January 28, 2014 at 11:00 am - Permalink

Perhaps this will work as well ...

Function WSFI(bw, tw, xl, xr)
    wave bw, tw
    variable xl, xr
   
    variable WSFIvalue
   
    duplicate/R=(xl, xr)/FREE bw, rw, bdw
    duplicate/R=(xl, xr)/FREE tw, tdw
    rw = (bdw - tdw)^2 / numpnts(rw)
    WSFIvalue = sum(rw)
   
    return WSFIvalue
end


Less overhead is here to create and then kill waves. You can do a print WSFI(...) and get the value printed when you want (rather than all of the time). If you really want to save the resultwave, then perhaps pass it by reference from the outside and return it modified.

--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAHuntsville

January 28, 2014 at 02:21 pm - Permalink

...even less overhead
Function WSFI(bw, tw, xl, xr)
    wave bw, tw
    variable xl, xr
 
    duplicate/R=(xl, xr)/FREE bw, bdw
    duplicate/R=(xl, xr)/FREE tw, tdw
    MatrixOP/FREE aa=Sum(magsqr(bdw - tdw))/numPoints(tdw)
    return aa[0]
End


The code above would be more efficient when the range [xl,xr] includes many points.

January 28, 2014 at 02:56 pm - Permalink

jjweimer wrote:
Perhaps this will work as well ...

Function WSFI(bw, tw, xl, xr)
    wave bw, tw
    variable xl, xr
   
    variable WSFIvalue
   
    duplicate/R=(xl, xr)/FREE bw, rw, bdw
    duplicate/R=(xl, xr)/FREE tw, tdw
    rw = (bdw - tdw)^2 / numpnts(rw)
    WSFIvalue = sum(rw)
   
    return WSFIvalue
end


Less overhead is here to create and then kill waves. You can do a print WSFI(...) and get the value printed when you want (rather than all of the time). If you really want to save the resultwave, then perhaps pass it by reference from the outside and return it modified.

--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAHuntsville


Thank you for the beautiful code (neat and simple). It worked.
Much to learn.



January 28, 2014 at 03:05 pm - Permalink

Igor wrote:
...even less overhead
Function WSFI(bw, tw, xl, xr)
    wave bw, tw
    variable xl, xr
 
    duplicate/R=(xl, xr)/FREE bw, bdw
    duplicate/R=(xl, xr)/FREE tw, tdw
    MatrixOP/FREE aa=Sum(magsqr(bdw - tdw))/numPoints(tdw)
    return aa[0]
End


The code above would be more efficient when the range [xl,xr] includes many points.


What to say!!! I like it!
Thanks for the code!

January 28, 2014 at 03:10 pm - Permalink

Igor wrote:
...even less overhead...


It seems I am going to have to learn MatrixOP in addition perhaps to python somewhere along the line.

And my students think, they can just learn to hack around in MathCad or Mathematica and be done with all their programming work for the rest of their careers. :-)

--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAHuntsville

January 28, 2014 at 03:37 pm - Permalink