Wave arithmetic for mismatched XY wave pairs

Average rating
(0 votes)

Performs addition, subtraction, multiplication or division on two XY pairs of waves. Determines the overlap X range and creates an XY output pair with names determined by the string outwavename; the X wave of the output pair has "_X" appended the the specified name. The output X wave contains all X values from both input pairs that fall within the common range. To do the arithmetic operations, a Y value is interpolated from each input Y wave at each of the output X values. The arithmetic is done using this pair of interpolated Y values. No attempt is made to eliminate duplicated X values.

// Constants to use for the "operation" input
Constant XYArith_Add=0
Constant XYArith_Subtract=1
Constant XYArith_Multiply=2
Constatn XYArith_Divide=3
 
Function UnmatchedXYPairArithmetic(WAVE axwave, WAVE aywave, WAVE bxwave, WAVE bywave, Variable operation, String outwavename)
 
	Sort axwave, axwave, aywave
	Sort bxwave, bxwave, bywave
 
	WaveStats/Q axwave
	Variable amin = V_min
	Variable aminloc = V_minloc
	Variable amax = V_max
	Variable amaxloc = V_maxloc
	WaveStats/Q bxwave
	Variable bmin = V_min
	Variable bminloc = V_minloc
	Variable bmax = V_max
	Variable bmaxloc = V_maxloc
 
	if ((amin > bmax) || (amax < bmin))
		DoAlert 0, "X waves don't overlap"
	endif
	Variable beginx = amin < bmin ? bmin : amin
	Variable endx = amax < bmax ? amax : bmax
	Variable beginp = binarysearch(axwave, beginx)+1
	Variable endp = binarysearch(axwave, endx)
	Duplicate/FREE/R=[beginp, endp] axwave,srcax
	Duplicate/FREE/R=[beginp, endp] aywave,srcay
	beginp = binarysearch(bxwave, beginx)+1
	endp = binarysearch(bxwave, endx)
	Duplicate/FREE/R=[beginp, endp] bxwave,srcbx
	Duplicate/FREE/R=[beginp, endp] bywave,srcby
 
	KillWaves/Z $(outwavename+"_X")
	Concatenate/NP/O {srcax, srcbx}, $(outwavename+"_X")/WAVE=xw
	Duplicate/O xw, $outwavename/WAVE=yw
 
	switch(operation)
		case 0:
			yw = aywave[binarysearchinterp(axwave, xw[p])] + bywave[binarysearchinterp(bxwave, xw[p])]
			break;
		case 1:
			yw = aywave[binarysearchinterp(axwave, xw[p])] - bywave[binarysearchinterp(bxwave, xw[p])]
			break;
		case 2:
			yw = aywave[binarysearchinterp(axwave, xw[p])] * bywave[binarysearchinterp(bxwave, xw[p])]
			break;
		case 3:
			yw = aywave[binarysearchinterp(axwave, xw[p])] / bywave[binarysearchinterp(bxwave, xw[p])]
			break;
	endswitch
end
 
// IsMonotonic() returns true if the wave has dwave(x)/dx >= 0 for all x.
Static Function IsMonotonic(WAVE wv)
 
	Differentiate/METH=1/EP=1 wv/D=difwave
	WAVE difwave
	Variable wmin = WaveMin(difwave)
	Variable wmax = WaveMax(difwave)
 
	KillWaves difwave
	return (wmin < 0 && wmax <= 0) || (wmin >= 0 && wmax > 0)
End

Back to top