diff --git a/Packages/PFFFT/Sources/PFFFT/FFTDoubleImpl.swift b/Packages/PFFFT/Sources/PFFFT/FFTDouble.swift similarity index 50% rename from Packages/PFFFT/Sources/PFFFT/FFTDoubleImpl.swift rename to Packages/PFFFT/Sources/PFFFT/FFTDouble.swift index 6d872f9..7b2699e 100644 --- a/Packages/PFFFT/Sources/PFFFT/FFTDoubleImpl.swift +++ b/Packages/PFFFT/Sources/PFFFT/FFTDouble.swift @@ -1,10 +1,16 @@ internal import PFFFTLib -public final class FFTDoubleImpl: FFTImplProtocol { +public final class FFTDouble: PFFFTProtocol { let ptr: OpaquePointer let n: Int let type: FFTType + static let sharedCache = SetupCache() + + public static func setup(for n: Int, type: FFTType) throws -> FFTDouble { + try sharedCache.get(for: n, type: type) + } + public init(n: Int, type: FFTType) throws { guard let ptr = pffftd_new_setup(Int32(n), pffft_transform_t(type)) else { throw FFTError.invalidSize } self.ptr = ptr @@ -12,37 +18,56 @@ public final class FFTDoubleImpl: FFTImplProtocol { self.type = type } - public func fft(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?, sign: FFTSign) { + + func transform(_ input: borrowing Buffer, _ output: borrowing Buffer, _ work: borrowing Buffer?, _ dir: pffft_direction_t) { + checkFftBufferCounts(n: n, type: type, input: input, output: output, work: work) + + let workAddress: UnsafeMutablePointer! = switch work { + case let .some(b): b.baseAddress + case .none: nil + } + pffftd_transform_ordered(ptr, input.baseAddress, output.baseAddress, workAddress, dir) + } + + func transformUnordered(_ input: borrowing Buffer, _ output: borrowing Buffer, _ work: borrowing Buffer?, _ dir: pffft_direction_t) { checkFftBufferCounts(n: n, type: type, input: input, output: output, work: work) let workAddress: UnsafeMutablePointer! = switch work { case let .some(b): b.baseAddress case .none: nil } - pffftd_transform_ordered(ptr, input.baseAddress, output.baseAddress, workAddress, pffft_direction_t(sign)) + pffftd_transform(ptr, input.baseAddress, output.baseAddress, workAddress, dir) } - public func fftUnordered(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?, sign: FFTSign) { - checkFftBufferCounts(n: n, type: type, input: input, output: output, work: work) - let workAddress: UnsafeMutablePointer! = switch work { - case let .some(b): b.baseAddress - case .none: nil - } - pffftd_transform(ptr, input.baseAddress, output.baseAddress, workAddress, pffft_direction_t(sign)) + public func forward(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?) { + transform(input, output, work, PFFFT_FORWARD) } - public func zReorder(input: borrowing Buffer, output: borrowing Buffer, sign: FFTSign) { + public func inverse(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?) { + transform(input, output, work, PFFFT_BACKWARD) + } + + public func forwardUnordered(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?) { + transformUnordered(input, output, work, PFFFT_FORWARD) + } + + public func inverseUnordered(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?) { + transformUnordered(input, output, work, PFFFT_BACKWARD) + } + + + public func reorderSpectrum(input: borrowing Buffer, output: borrowing Buffer) { checkFftBufferCounts(n: n, type: type, input: input, output: output, work: nil) - pffftd_zreorder(ptr, input.baseAddress, output.baseAddress, pffft_direction_t(sign)) + pffftd_zreorder(ptr, input.baseAddress, output.baseAddress, PFFFT_FORWARD) } - public func zConvolveAccumulate(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Double) { + public func convolveAccumulate(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Double) { checkConvolveBufferCounts(n: n, type: type, a: dftA, b: dftB, ab: dftAB) pffftd_zconvolve_accumulate(ptr, dftA.baseAddress, dftB.baseAddress, dftAB.baseAddress, scaling) } - public func zConvolve(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Double) { + public func convolve(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Double) { checkConvolveBufferCounts(n: n, type: type, a: dftA, b: dftB, ab: dftAB) pffftd_zconvolve_no_accu(ptr, dftA.baseAddress, dftB.baseAddress, dftAB.baseAddress, scaling) } @@ -65,5 +90,5 @@ public final class FFTDoubleImpl: FFTImplProtocol { } extension Double: FFTElement { - public typealias FFTImpl = FFTDoubleImpl + public typealias FFTImpl = FFTDouble } diff --git a/Packages/PFFFT/Sources/PFFFT/FFTFloatImpl.swift b/Packages/PFFFT/Sources/PFFFT/FFTFloat.swift similarity index 73% rename from Packages/PFFFT/Sources/PFFFT/FFTFloatImpl.swift rename to Packages/PFFFT/Sources/PFFFT/FFTFloat.swift index c61008f..9b0d05a 100644 --- a/Packages/PFFFT/Sources/PFFFT/FFTFloatImpl.swift +++ b/Packages/PFFFT/Sources/PFFFT/FFTFloat.swift @@ -1,10 +1,16 @@ internal import PFFFTLib -public final class FFTFloatImpl: FFTImplProtocol { +public final class FFTFloat: PFFFTProtocol { let ptr: OpaquePointer let n: Int let type: FFTType + static let sharedCache = SetupCache() + + public static func setup(for n: Int, type: FFTType) throws -> FFTFloat { + try sharedCache.get(for: n, type: type) + } + public init(n: Int, type: FFTType) throws { guard let ptr = pffft_new_setup(Int32(n), pffft_transform_t(type)) else { throw FFTError.invalidSize } self.ptr = ptr @@ -12,7 +18,7 @@ public final class FFTFloatImpl: FFTImplProtocol { self.type = type } - public func fft(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?, sign: FFTSign) { + public func forward(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?, sign: FFTSign) { checkFftBufferCounts(n: n, type: type, input: input, output: output, work: work) let workAddress: UnsafeMutablePointer! = switch work { @@ -32,17 +38,17 @@ public final class FFTFloatImpl: FFTImplProtocol { pffft_transform(ptr, input.baseAddress, output.baseAddress, workAddress, pffft_direction_t(sign)) } - public func zReorder(input: borrowing Buffer, output: borrowing Buffer, sign: FFTSign) { + public func reorderSpectrum(input: borrowing Buffer, output: borrowing Buffer) { checkFftBufferCounts(n: n, type: type, input: input, output: output, work: nil) - pffft_zreorder(ptr, input.baseAddress, output.baseAddress, pffft_direction_t(sign)) + pffft_zreorder(ptr, input.baseAddress, output.baseAddress, PFFFT_FORWARD) } - public func zConvolveAccumulate(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Float) { + public func convolveAccumulate(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Float) { checkConvolveBufferCounts(n: n, type: type, a: dftA, b: dftB, ab: dftAB) pffft_zconvolve_accumulate(ptr, dftA.baseAddress, dftB.baseAddress, dftAB.baseAddress, scaling) } - public func zConvolve(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Float) { + public func convolve(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Float) { checkConvolveBufferCounts(n: n, type: type, a: dftA, b: dftB, ab: dftAB) pffft_zconvolve_no_accu(ptr, dftA.baseAddress, dftB.baseAddress, dftAB.baseAddress, scaling) } @@ -65,5 +71,5 @@ public final class FFTFloatImpl: FFTImplProtocol { } extension Float: FFTElement { - public typealias FFTImpl = FFTFloatImpl + public typealias FFTImpl = FFTFloat } diff --git a/Packages/PFFFT/Sources/PFFFT/PFFFT.swift b/Packages/PFFFT/Sources/PFFFT/PFFFT.swift index 002a7f6..007f6d5 100644 --- a/Packages/PFFFT/Sources/PFFFT/PFFFT.swift +++ b/Packages/PFFFT/Sources/PFFFT/PFFFT.swift @@ -17,12 +17,22 @@ public enum FFTError: Error { } public protocol FFTElement { - associatedtype FFTImpl: FFTImplProtocol + associatedtype FFTImpl: PFFFTProtocol } -public protocol FFTImplProtocol: ~Copyable { +public protocol PFFFTProtocol { associatedtype Element + /// Get an FFT interface for the given size and FFT Type. + /// + /// This call is backed by a global cache to avoid repeated setup costs. + /// - Parameters: + /// - n: The size of the FFT. + /// - type: The type of FFT. + /// - Returns: An FFT interface. + /// - Throws: `FFTError.invalidSize` if the size is invalid. + static func setup(for n: Int, type: FFTType) throws -> Self + /// Initialize the FFT implementation with the given size and type. /// - Parameters: /// - n: The size of the FFT. @@ -30,7 +40,7 @@ public protocol FFTImplProtocol: ~Copyable { /// - Throws: `FFTError.invalidSize` if the size is invalid. init(n: Int, type: FFTType) throws - /// Perform a forward or backward FFT on the input buffer. + /// Perform a forward FFT on the input buffer. /// /// The input and output buffers may be the same. /// The data is stores in order as expected (interleaved complex components ordered by frequency). @@ -50,9 +60,11 @@ public protocol FFTImplProtocol: ~Copyable { /// - output: The output buffer. /// - work: An optional work buffer. Must have capacity of at least `n` for real FFTs and `2 * n` for complex FFTs. /// - sign: The direction of the FFT. - func fft(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?, sign: FFTSign) + func forward(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?) - /// Perform a forward or backward FFT on the input buffer, with implementation defined order. + func inverse(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?) + + /// Perform a forward FFT on the input buffer, with implementation defined order. /// /// This function behaves similarly to `fft` however the z-domain data is stored in most efficient ordering, /// which is suitable for transforming back with this function, or for convolution. @@ -61,9 +73,11 @@ public protocol FFTImplProtocol: ~Copyable { /// - output: The output buffer. /// - work: An optional work buffer. Must have capacity of at least `n` for real FFTs and `2 * n` for complex FFTs. /// - sign: The direction of the FFT. - func fftUnordered(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?, sign: FFTSign) + func forwardUnordered(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?) - func zReorder(input: borrowing Buffer, output: borrowing Buffer, sign: FFTSign) + func inverseUnordered(input: borrowing Buffer, output: borrowing Buffer, work: borrowing Buffer?) + + func reorder(input: borrowing Buffer, output: borrowing Buffer, sign: FFTSign) /// Perform a convolution of two complex signals in the frequency domain. /// @@ -74,7 +88,7 @@ public protocol FFTImplProtocol: ~Copyable { /// - dftB: The second input buffer of frequency domain data. /// - dftAB: The output buffer of frequency domain data. /// - scaling: The scaling factor to apply to the result. - func zConvolveAccumulate(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Element) + func convolveAccumulate(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Element) /// Perform a convolution of two complex signals in the frequency domain. /// @@ -85,7 +99,7 @@ public protocol FFTImplProtocol: ~Copyable { /// - dftB: The second input buffer of frequency domain data. /// - dftAB: The output buffer of frequency domain data. /// - scaling: The scaling factor to apply to the result. - func zConvolve(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Element) + func convolve(dftA: borrowing Buffer, dftB: borrowing Buffer, dftAB: borrowing Buffer, scaling: Element) /// Returns the minimum FFT size for the given type. /// diff --git a/Sources/EcgSynKit/EcgSynKit.swift b/Sources/EcgSynKit/EcgSynKit.swift index 81e2ab7..96e2860 100644 --- a/Sources/EcgSynKit/EcgSynKit.swift +++ b/Sources/EcgSynKit/EcgSynKit.swift @@ -63,10 +63,25 @@ func stdev(_ data: [Double]) -> Double { return sqrt(data.lazy.map { ($0 - mean) * ($0 - mean) }.reduce(0.0, +) / (n - 1)) } -let fftlib = SetupCache() +struct RRProcess : ~Copyable { + let nrr: Int -func rrprocess(params: Parameters, nrr: Int) -> [Double] { - let w1 = 2.0 * .pi * params.flo + let input: Buffer + let output: Buffer + let work: Buffer + let fft: FFTDouble + + init(nrr: Int) { + self.nrr = nrr + self.input = Buffer(capacity: nrr + 2) + self.output = Buffer(capacity: nrr) + self.work = Buffer(capacity: nrr) + self.fft = try! FFTDouble.setup(for: nrr, type: .real) + } + + func generate(params: Parameters) -> [Double] { + + let w1 = 2.0 * .pi * params.flo let w2 = 2.0 * .pi * params.fhi let c1 = 2.0 * .pi * params.flostd let c2 = 2.0 * .pi * params.fhistd @@ -80,10 +95,6 @@ func rrprocess(params: Parameters, nrr: Int) -> [Double] { let sf = Double(params.sfInternal) let df = sf / Double(nrr) - - let input = Buffer(capacity: (nrr + 2)) - let output = Buffer(capacity: nrr) - input.withUnsafeMutableBufferPointer { swc in for i in 0 ..< nrr / 2 + 1 { let w = df * Double(i) * 2.0 * .pi @@ -103,7 +114,7 @@ func rrprocess(params: Parameters, nrr: Int) -> [Double] { swc[1] = swc[nrr] } - try! fftlib.get(for: nrr, type: .real).fft(input: input, output: output, work: nil, sign: .backward) + self.fft.forward(input: input, output: output, work: nil, sign: .backward) var rr = output.withUnsafeMutableBufferPointer { outptr in return outptr.map { $0 * 1.0 / Double(nrr) } @@ -116,6 +127,8 @@ func rrprocess(params: Parameters, nrr: Int) -> [Double] { rr[i] = rr[i] * ratio + rrmean } return rr + + } }