SF: Add instantaneous pair method for apfrequency operation

The method number is 3.

Packages/MIES/MIES_SweepFormula.ipf

@@ -52,6 +52,7 @@ static Constant SF_MAX_NUMPOINTS_FOR_MARKERS = 1000
 static Constant SF_APFREQUENCY_FULL          = 0x0
 static Constant SF_APFREQUENCY_INSTANTANEOUS = 0x1
 static Constant SF_APFREQUENCY_APCOUNT       = 0x2
+static Constant SF_APFREQUENCY_INSTANTANEOUS_PAIR = 0x3
 
 static StrConstant SF_OP_MINUS = "-"
 static StrConstant SF_OP_PLUS = "+"
@@ -4295,7 +4296,8 @@ End
 // apfrequency(data, [frequency calculation method], [spike detection crossing level])
 static Function/WAVE SF_OperationApFrequency(variable jsonId, string jsonPath, string graph)
 
-	variable numArgs, i
+	variable i, numArgs, keepX
+	string xUnit, xLabel
 
 	numArgs = SFH_GetNumberOfArguments(jsonID, jsonPath)
 	SFH_ASSERT(numArgs <=3, "ApFrequency has 3 arguments at most.")
@@ -4314,15 +4316,32 @@ static Function/WAVE SF_OperationApFrequency(variable jsonId, string jsonPath, s
 		WAVE method = SFH_GetArgumentSingle(jsonID, jsonPath, graph, SF_OP_APFREQUENCY, 1, checkExist=1)
 		SFH_ASSERT(DimSize(method, ROWS) == 1, "Too many input values for parameter method")
 		SFH_ASSERT(IsNumericWave(method), "method parameter must be numeric.")
-		SFH_ASSERT(method[0] == SF_APFREQUENCY_FULL || method[0] == SF_APFREQUENCY_INSTANTANEOUS ||  method[0] == SF_APFREQUENCY_APCOUNT, "method parameter is invalid")
+		SFH_ASSERT(method[0] == SF_APFREQUENCY_FULL || method[0] == SF_APFREQUENCY_INSTANTANEOUS ||  method[0] == SF_APFREQUENCY_APCOUNT || method[0] == SF_APFREQUENCY_INSTANTANEOUS_PAIR, "method parameter is invalid")
 	else
 		Make/FREE method = {SF_APFREQUENCY_FULL}
 	endif
 
 	WAVE/WAVE output = SFH_CreateSFRefWave(graph, SF_OP_APFREQUENCY, DimSize(input, ROWS))
 	output = SF_OperationApFrequencyImpl(input[p], level[0], method[0])
 
-	SFH_TransferFormulaDataWaveNoteAndMeta(input, output, SF_OP_APFREQUENCY, SF_DATATYPE_APFREQUENCY)
+	if(method[0] == SF_APFREQUENCY_INSTANTANEOUS_PAIR)
+		if(DimSize(input, ROWS))
+			WAVE/Z data = input[0]
+			if(WaveExists(data))
+				xUnit = WaveUnits(data, ROWS)
+				if(!CmpStr(xUnit, "ms"))
+					xLabel = "kHz"
+				elseif(IsEmpty(xUnit))
+					xLabel = "1/x"
+				else
+					xLabel = "1/" + xUnit
+				endif
+				JWN_SetStringInWaveNote(output, SF_META_XAXISLABEL, xLabel)
+			endif
+		endif
+		keepX = 1
+	endif
+	SFH_TransferFormulaDataWaveNoteAndMeta(input, output, SF_OP_APFREQUENCY, SF_DATATYPE_APFREQUENCY, keepX=keepX)
 
 	return SFH_GetOutputForExecutor(output, graph, SF_OP_APFREQUENCY)
 End
@@ -4334,24 +4353,21 @@ static Function/WAVE SF_OperationApFrequencyImpl(WAVE data, variable level, vari
 	WAVE levels = FindLevelWrapper(data, level, FINDLEVEL_EDGE_INCREASING, FINDLEVEL_MODE_MULTI)
 	numSets = DimSize(levels, ROWS)
 	Make/FREE/N=(numSets) levelPerSet = str2num(GetDimLabel(levels, ROWS, p))
+	WAVE/Z pairTimes = $""
 
 	// @todo we assume that the x-axis of data has a ms scale for FULL/INSTANTANEOUS
 	switch(method)
 		case SF_APFREQUENCY_FULL:
 			Make/N=(numSets)/D/FREE outD = levelPerSet[p] / (DimDelta(data, ROWS) * DimSize(data, ROWS) * MILLI_TO_ONE)
 			break
 		case SF_APFREQUENCY_INSTANTANEOUS:
-			Make/N=(numSets)/D/FREE outD
-
-			for(i = 0; i < numSets; i += 1)
-				if(levelPerSet[i] <= 1)
-					outD[i] = 0
-				else
-					Make/FREE/D/N=(levelPerSet[i] - 1) distances
-					distances[0, levelPerSet[i] - 2] = levels[i][p + 1] - levels[i][p]
-					outD[i] = 1.0 / (Mean(distances) * MILLI_TO_ONE)
-				endif
-			endfor
+			WAVE outD = SF_ApFrequencyInstantaneous(levels, levelPerSet)
+			break
+		case SF_APFREQUENCY_INSTANTANEOUS_PAIR:
+			[outD, pairTimes] = SF_ApFrequencyInstantaneousPairs(levels, levelPerSet)
+			if(WaveExists(outD))
+				JWN_SetWaveInWaveNote(outD, SF_META_XVALUES, pairTimes)
+			endif
 			break
 		case SF_APFREQUENCY_APCOUNT:
 			Make/N=(numSets)/D/FREE outD = levelPerSet[p]
@@ -4361,6 +4377,46 @@ static Function/WAVE SF_OperationApFrequencyImpl(WAVE data, variable level, vari
 	return outD
 End
 
+static Function/WAVE SF_ApFrequencyInstantaneous(WAVE levels, WAVE levelPerSet)
+
+	variable i
+	variable numSets = DimSize(levels, ROWS)
+
+	Make/N=(numSets)/D/FREE result
+
+	for(i = 0; i < numSets; i += 1)
+		if(levelPerSet[i] <= 1)
+			result[i] = 0
+		else
+			Make/FREE/D/N=(levelPerSet[i] - 1) distances
+			distances[0, levelPerSet[i] - 2] = levels[i][p + 1] - levels[i][p]
+			result[i] = 1.0 / (Mean(distances) * MILLI_TO_ONE)
+		endif
+	endfor
+
+	return result
+End
+
+static Function [WAVE/D result, WAVE/D pairTimes] SF_ApFrequencyInstantaneousPairs(WAVE levels, WAVE levelPerSet)
+
+	SFH_ASSERT(DimSize(levels, ROWS) == 1, "Expected a single set")
+
+	if(levelPerSet[0] == 0)
+		return [$"", $""]
+	endif
+
+	Make/FREE/D/N=(max(levelPerSet[0] - 1, 1)) result, pairTimes
+	if(levelPerSet[0] == 1)
+		result[0] = 0
+		pairTimes[0] = levels[0][0]
+	else
+		result = 1.0 / ((levels[0][p + 1] - levels[0][p]) * MILLI_TO_ONE)
+		pairTimes = levels[0][p]
+	endif
+
+	return [result, pairTimes]
+End
+
 // `store(name, ...)`
 static Function/WAVE SF_OperationStore(variable jsonId, string jsonPath, string graph)