before protocol removal
This commit is contained in:
@@ -14,7 +14,6 @@ let package = Package(
|
||||
],
|
||||
dependencies: [
|
||||
.package(url: "https://github.com/apple/swift-algorithms", from: "1.2.0"),
|
||||
.package(path: "./Packages/KissFFT/"),
|
||||
.package(path: "./Packages/PFFFT/"),
|
||||
],
|
||||
targets: [
|
||||
|
||||
@@ -6,12 +6,14 @@ import PackageDescription
|
||||
let package = Package(
|
||||
name: "PFFFT",
|
||||
products: [
|
||||
// Products define the executables and libraries a package produces, making them visible to other packages.
|
||||
.library(
|
||||
name: "PFFFT",
|
||||
targets: ["PFFFT", "PFFFTLib"]
|
||||
),
|
||||
],
|
||||
dependencies: [
|
||||
.package(url: "https://github.com/apple/swift-numerics", from: "1.0.0"),
|
||||
],
|
||||
targets: [
|
||||
.target(
|
||||
name: "PFFFTLib",
|
||||
@@ -25,7 +27,7 @@ let package = Package(
|
||||
),
|
||||
.target(
|
||||
name: "PFFFT",
|
||||
dependencies: ["PFFFTLib"]
|
||||
dependencies: ["PFFFTLib", .product(name: "Numerics", package: "swift-numerics")]
|
||||
),
|
||||
.testTarget(
|
||||
name: "PFFFTTests",
|
||||
|
||||
@@ -1,27 +1,27 @@
|
||||
let bufferAlignment = 32
|
||||
|
||||
@frozen
|
||||
public struct Buffer<Element>: ~Copyable {
|
||||
let buffer: UnsafeMutableBufferPointer<Element>
|
||||
public struct Buffer<T>: ~Copyable {
|
||||
let buffer: UnsafeMutableBufferPointer<T>
|
||||
var count: Int { buffer.count }
|
||||
var baseAddress: UnsafeMutablePointer<Element> { buffer.baseAddress! }
|
||||
var baseAddress: UnsafeMutablePointer<T> { buffer.baseAddress! }
|
||||
|
||||
public init(capacity: Int) {
|
||||
buffer = UnsafeMutableRawBufferPointer.allocate(
|
||||
byteCount: MemoryLayout<Element>.stride * capacity,
|
||||
byteCount: MemoryLayout<T>.stride * capacity,
|
||||
alignment: bufferAlignment
|
||||
).bindMemory(to: Element.self)
|
||||
).bindMemory(to: T.self)
|
||||
}
|
||||
|
||||
deinit {
|
||||
buffer.deallocate()
|
||||
}
|
||||
|
||||
public func withUnsafeMutableBufferPointer<R>(_ body: (UnsafeMutableBufferPointer<Element>) throws -> R) rethrows -> R {
|
||||
public func withUnsafeMutableBufferPointer<R>(_ body: (UnsafeMutableBufferPointer<T>) throws -> R) rethrows -> R {
|
||||
try body(buffer)
|
||||
}
|
||||
|
||||
public func withUnsafeBufferPointer<R>(_ body: (UnsafeBufferPointer<Element>) throws -> R) rethrows -> R {
|
||||
public func withUnsafeBufferPointer<R>(_ body: (UnsafeBufferPointer<T>) throws -> R) rethrows -> R {
|
||||
try body(UnsafeBufferPointer(buffer))
|
||||
}
|
||||
|
||||
|
||||
276
Packages/PFFFT/Sources/PFFFT/FFT.swift
Normal file
276
Packages/PFFFT/Sources/PFFFT/FFT.swift
Normal file
@@ -0,0 +1,276 @@
|
||||
internal import PFFFTLib
|
||||
import ComplexModule
|
||||
import RealModule
|
||||
|
||||
@frozen
|
||||
public enum FFTType {
|
||||
case real
|
||||
case complex
|
||||
}
|
||||
|
||||
@frozen
|
||||
public enum FFTSign: Int {
|
||||
case forward = -1
|
||||
case backward = 1
|
||||
}
|
||||
|
||||
public enum FFTError: Error {
|
||||
case invalidSize
|
||||
}
|
||||
|
||||
public protocol FFTElemental {
|
||||
associatedtype ScalarType: FFTScalar
|
||||
associatedtype ComplexType = Complex<ScalarType>
|
||||
|
||||
static func setupPfft(_ n: Int, _ type: FFTType) throws -> OpaquePointer
|
||||
static func pffftMinFftSize(_ type: FFTType) -> Int
|
||||
static func pffftIsValidSize(_ n: Int, _ type: FFTType) -> Bool
|
||||
static func pffftNearestValidSize(_ n: Int, _ type: FFTType, _ higher: Bool) -> Int
|
||||
}
|
||||
|
||||
public protocol FFTScalar: Real {
|
||||
static func pffftTransformOrdered(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Self>, _ output: UnsafeMutablePointer<Self>, _ work: UnsafeMutablePointer<Self>?, _ dir: FFTSign)
|
||||
static func pffftTransform(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Self>, _ output: UnsafeMutablePointer<Self>, _ work: UnsafeMutablePointer<Self>?, _ dir: FFTSign)
|
||||
static func pffftZreorder(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Self>, _ output: UnsafeMutablePointer<Self>, _ dir: FFTSign)
|
||||
static func pffftZconvolveAccumulate(_ ptr: OpaquePointer, _ dftA: UnsafeMutablePointer<Self>, _ dftB: UnsafeMutablePointer<Self>, _ dftAB: UnsafeMutablePointer<Self>, _ scaling: Self)
|
||||
static func pffftZconvolveNoAccu(_ ptr: OpaquePointer, _ dftA: UnsafeMutablePointer<Self>, _ dftB: UnsafeMutablePointer<Self>, _ dftAB: UnsafeMutablePointer<Self>, _ scaling: Self)
|
||||
}
|
||||
|
||||
extension Float: FFTScalar {
|
||||
public static func pffftTransformOrdered(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Self>, _ output: UnsafeMutablePointer<Self>, _ work: UnsafeMutablePointer<Self>?, _ dir: FFTSign) {
|
||||
pffft_transform_ordered(ptr, input, output, work, pffft_direction_t(dir))
|
||||
}
|
||||
|
||||
public static func pffftTransform(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Self>, _ output: UnsafeMutablePointer<Self>, _ work: UnsafeMutablePointer<Self>?, _ dir: FFTSign) {
|
||||
pffft_transform(ptr, input, output, work, pffft_direction_t(dir))
|
||||
}
|
||||
|
||||
public static func pffftZreorder(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Self>, _ output: UnsafeMutablePointer<Self>, _ dir: FFTSign) {
|
||||
pffft_zreorder(ptr, input, output, pffft_direction_t(dir))
|
||||
}
|
||||
|
||||
public static func pffftZconvolveAccumulate(_ ptr: OpaquePointer, _ dftA: UnsafeMutablePointer<Self>, _ dftB: UnsafeMutablePointer<Self>, _ dftAB: UnsafeMutablePointer<Self>, _ scaling: Self) {
|
||||
pffft_zconvolve_accumulate(ptr, dftA, dftB, dftAB, scaling)
|
||||
}
|
||||
|
||||
public static func pffftZconvolveNoAccu(_ ptr: OpaquePointer, _ dftA: UnsafeMutablePointer<Self>, _ dftB: UnsafeMutablePointer<Self>, _ dftAB: UnsafeMutablePointer<Self>, _ scaling: Self) {
|
||||
pffft_zconvolve_no_accu(ptr, dftA, dftB, dftAB, scaling)
|
||||
}
|
||||
}
|
||||
|
||||
extension Double: FFTScalar {
|
||||
public static func pffftTransformOrdered(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Double>, _ output: UnsafeMutablePointer<Double>, _ work: UnsafeMutablePointer<Double>?, _ dir: FFTSign) {
|
||||
pffftd_transform_ordered(ptr, input, output, work, pffft_direction_t(dir))
|
||||
}
|
||||
|
||||
public static func pffftTransform(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Double>, _ output: UnsafeMutablePointer<Double>, _ work: UnsafeMutablePointer<Double>?, _ dir: FFTSign) {
|
||||
pffftd_transform(ptr, input, output, work, pffft_direction_t(dir))
|
||||
}
|
||||
|
||||
public static func pffftZreorder(_ ptr: OpaquePointer, _ input: UnsafeMutablePointer<Double>, _ output: UnsafeMutablePointer<Double>, _ dir: FFTSign) {
|
||||
pffftd_zreorder(ptr, input, output, pffft_direction_t(dir))
|
||||
}
|
||||
|
||||
public static func pffftZconvolveAccumulate(_ ptr: OpaquePointer, _ dftA: UnsafeMutablePointer<Double>, _ dftB: UnsafeMutablePointer<Double>, _ dftAB: UnsafeMutablePointer<Double>, _ scaling: Double) {
|
||||
pffftd_zconvolve_accumulate(ptr, dftA, dftB, dftAB, scaling)
|
||||
}
|
||||
|
||||
public static func pffftZconvolveNoAccu(_ ptr: OpaquePointer, _ dftA: UnsafeMutablePointer<Double>, _ dftB: UnsafeMutablePointer<Double>, _ dftAB: UnsafeMutablePointer<Double>, _ scaling: Double) {
|
||||
pffftd_zconvolve_no_accu(ptr, dftA, dftB, dftAB, scaling)
|
||||
}
|
||||
}
|
||||
|
||||
extension Complex: FFTElemental where RealType: FFTElemental & FFTScalar {
|
||||
public typealias ScalarType = RealType
|
||||
|
||||
public static func setupPfft(_ n: Int, _: FFTType) throws -> OpaquePointer {
|
||||
return try ScalarType.setupPfft(n, .complex)
|
||||
}
|
||||
|
||||
public static func pffftMinFftSize(_: FFTType) -> Int {
|
||||
return ScalarType.pffftMinFftSize(.complex)
|
||||
}
|
||||
|
||||
public static func pffftIsValidSize(_ n: Int, _: FFTType) -> Bool {
|
||||
return ScalarType.pffftIsValidSize(n, .complex)
|
||||
}
|
||||
|
||||
public static func pffftNearestValidSize(_ n: Int, _: FFTType, _ higher: Bool) -> Int {
|
||||
return ScalarType.pffftNearestValidSize(n, .complex, higher)
|
||||
}
|
||||
}
|
||||
|
||||
extension Double: FFTElemental {
|
||||
public typealias ScalarType = Double
|
||||
|
||||
public static func setupPfft(_ n: Int, _ type: FFTType) throws -> OpaquePointer {
|
||||
guard let ptr = pffftd_new_setup(Int32(n), pffft_transform_t(type)) else { throw FFTError.invalidSize }
|
||||
return ptr
|
||||
}
|
||||
|
||||
public static func pffftMinFftSize(_ type: FFTType) -> Int {
|
||||
Int(pffftd_min_fft_size(pffft_transform_t(type)))
|
||||
}
|
||||
|
||||
public static func pffftIsValidSize(_ n: Int, _ type: FFTType) -> Bool {
|
||||
pffftd_is_valid_size(Int32(n), pffft_transform_t(type)) != 0
|
||||
}
|
||||
|
||||
public static func pffftNearestValidSize(_ n: Int, _ type: FFTType, _ higher: Bool) -> Int {
|
||||
Int(pffftd_nearest_transform_size(Int32(n), pffft_transform_t(type), higher ? 1 : 0))
|
||||
}
|
||||
}
|
||||
|
||||
extension Float: FFTElemental {
|
||||
public typealias ScalarType = Float
|
||||
|
||||
public static func setupPfft(_ n: Int, _ type: FFTType) throws -> OpaquePointer {
|
||||
guard let ptr = pffft_new_setup(Int32(n), pffft_transform_t(type)) else { throw FFTError.invalidSize }
|
||||
return ptr
|
||||
}
|
||||
|
||||
public static func pffftMinFftSize(_ type: FFTType) -> Int {
|
||||
Int(pffft_min_fft_size(pffft_transform_t(type)))
|
||||
}
|
||||
|
||||
public static func pffftIsValidSize(_ n: Int, _ type: FFTType) -> Bool {
|
||||
pffft_is_valid_size(Int32(n), pffft_transform_t(type)) != 0
|
||||
}
|
||||
|
||||
public static func pffftNearestValidSize(_ n: Int, _ type: FFTType, _ higher: Bool) -> Int {
|
||||
Int(pffft_nearest_transform_size(Int32(n), pffft_transform_t(type), higher ? 1 : 0))
|
||||
}
|
||||
}
|
||||
|
||||
public final class FFT<T: FFTElemental> {
|
||||
public typealias ComplexType = T.ComplexType
|
||||
public typealias ScalarType = T.ScalarType
|
||||
|
||||
let ptr: OpaquePointer
|
||||
let n: Int
|
||||
let work: Buffer<ScalarType>?
|
||||
|
||||
/// Initialize the FFT implementation with the given size and type.
|
||||
/// - Parameters:
|
||||
/// - n: The size of the FFT.
|
||||
/// - type: The type of FFT.
|
||||
/// - Throws: `FFTError.invalidSize` if the size is invalid.
|
||||
public init(n: Int) throws {
|
||||
ptr = try T.setupPfft(n, .real)
|
||||
self.n = n
|
||||
let capacity = T.self == Complex<ScalarType>.self ? 2 * n : n
|
||||
work = n > 4096 ? Buffer<ScalarType>(capacity: capacity) : nil
|
||||
}
|
||||
|
||||
public func makeSignalBuffer(extra _: Int = 0) -> Buffer<T> {
|
||||
Buffer(capacity: n)
|
||||
}
|
||||
|
||||
public func makeSpectrumBuffer(extra _: Int = 0) -> Buffer<ComplexType> {
|
||||
Buffer(capacity: n)
|
||||
}
|
||||
|
||||
@inline(__always)
|
||||
func toAddress(_ work: borrowing Buffer<ScalarType>?) -> UnsafeMutablePointer<ScalarType>? {
|
||||
switch work {
|
||||
case let .some(b): return b.baseAddress
|
||||
case .none: return nil
|
||||
}
|
||||
}
|
||||
|
||||
@inline(__always)
|
||||
func rebind<I>(_ buffer: borrowing Buffer<I>) -> UnsafeMutablePointer<ScalarType>! {
|
||||
UnsafeMutableRawBufferPointer(buffer.buffer).bindMemory(to: ScalarType.self).baseAddress
|
||||
}
|
||||
|
||||
@inline(__always)
|
||||
func checkFftBufferCounts(signal: borrowing Buffer<T>, spectrum: borrowing Buffer<ComplexType>) {
|
||||
guard signal.count >= n else {
|
||||
fatalError("signal buffer too small")
|
||||
}
|
||||
guard spectrum.count >= n else {
|
||||
fatalError("spectrum buffer too small")
|
||||
}
|
||||
}
|
||||
|
||||
@inline(__always)
|
||||
func checkFftInternalLayoutBufferCounts(signal: borrowing Buffer<T>, spectrum: borrowing Buffer<ScalarType>) {
|
||||
guard signal.count >= n else {
|
||||
fatalError("signal buffer too small")
|
||||
}
|
||||
guard spectrum.count >= (T.self == Complex<ScalarType>.self ? 2 * n : n) else {
|
||||
fatalError("spectrum buffer too small")
|
||||
}
|
||||
}
|
||||
|
||||
@inline(__always)
|
||||
func checkConvolveBufferCounts(a: borrowing Buffer<ScalarType>, b: borrowing Buffer<ScalarType>, ab: borrowing Buffer<ScalarType>) {
|
||||
let minCount = T.self == Complex<ScalarType>.self ? 2 * n : n
|
||||
|
||||
guard a.count >= minCount else {
|
||||
fatalError("a buffer too small")
|
||||
}
|
||||
guard b.count >= minCount else {
|
||||
fatalError("b buffer too small")
|
||||
}
|
||||
guard ab.count >= minCount else {
|
||||
fatalError("ab buffer too small")
|
||||
}
|
||||
}
|
||||
|
||||
public func forward(signal: borrowing Buffer<T>, spectrum: borrowing Buffer<ComplexType>) {
|
||||
checkFftBufferCounts(signal: signal, spectrum: spectrum)
|
||||
ScalarType.pffftTransformOrdered(ptr, rebind(signal), rebind(spectrum), toAddress(work), .forward)
|
||||
}
|
||||
|
||||
public func inverse(spectrum: borrowing Buffer<ComplexType>, signal: borrowing Buffer<T>) {
|
||||
checkFftBufferCounts(signal: signal, spectrum: spectrum)
|
||||
ScalarType.pffftTransformOrdered(ptr, rebind(spectrum), rebind(signal), toAddress(work), .backward)
|
||||
}
|
||||
|
||||
public func forwardToInternalLayout(signal: borrowing Buffer<T>, spectrum: borrowing Buffer<ScalarType>) {
|
||||
checkFftInternalLayoutBufferCounts(signal: signal, spectrum: spectrum)
|
||||
ScalarType.pffftTransform(ptr, rebind(signal), spectrum.baseAddress, toAddress(work), .forward)
|
||||
}
|
||||
|
||||
public func inverseFromInternalLayout(spectrum: borrowing Buffer<ScalarType>, signal: borrowing Buffer<T>) {
|
||||
checkFftInternalLayoutBufferCounts(signal: signal, spectrum: spectrum)
|
||||
ScalarType.pffftTransform(ptr, spectrum.baseAddress, rebind(signal), toAddress(work), .backward)
|
||||
}
|
||||
|
||||
public func reorder(spectrum: borrowing Buffer<ScalarType>, output: borrowing Buffer<ComplexType>) {
|
||||
guard spectrum.count >= (T.self == Complex<ScalarType>.self ? 2 * n : n) else {
|
||||
fatalError("signal buffer too small")
|
||||
}
|
||||
guard output.count >= n else {
|
||||
fatalError("outputbuffer too small")
|
||||
}
|
||||
ScalarType.pffftZreorder(ptr, spectrum.baseAddress, rebind(output), .forward)
|
||||
}
|
||||
|
||||
public func convolve(dftA: borrowing Buffer<ScalarType>, dftB: borrowing Buffer<ScalarType>, dftAB: borrowing Buffer<ScalarType>, scaling: ScalarType) {
|
||||
checkConvolveBufferCounts(a: dftA, b: dftB, ab: dftAB)
|
||||
ScalarType.pffftZconvolveNoAccu(ptr, dftA.baseAddress, dftB.baseAddress, dftAB.baseAddress, scaling)
|
||||
}
|
||||
|
||||
public func convolveAccumulate(dftA: borrowing Buffer<ScalarType>, dftB: borrowing Buffer<ScalarType>, dftAB: borrowing Buffer<ScalarType>, scaling: ScalarType) {
|
||||
checkConvolveBufferCounts(a: dftA, b: dftB, ab: dftAB)
|
||||
ScalarType.pffftZconvolveAccumulate(ptr, dftA.baseAddress, dftB.baseAddress, dftAB.baseAddress, scaling)
|
||||
}
|
||||
|
||||
public static func minFftSize() -> Int {
|
||||
T.pffftMinFftSize(.real)
|
||||
}
|
||||
|
||||
public static func isValidSize(_ n: Int) -> Bool {
|
||||
T.pffftIsValidSize(n, .real)
|
||||
}
|
||||
|
||||
public static func nearestValidSize(_ n: Int, higher: Bool) -> Int {
|
||||
T.pffftNearestValidSize(n, .real, higher)
|
||||
}
|
||||
|
||||
deinit {
|
||||
pffftd_destroy_setup(ptr)
|
||||
}
|
||||
}
|
||||
@@ -1,94 +0,0 @@
|
||||
internal import PFFFTLib
|
||||
|
||||
public final class FFTDouble: PFFFTProtocol {
|
||||
let ptr: OpaquePointer
|
||||
let n: Int
|
||||
let type: FFTType
|
||||
|
||||
static let sharedCache = SetupCache<Element>()
|
||||
|
||||
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
|
||||
self.n = n
|
||||
self.type = type
|
||||
}
|
||||
|
||||
|
||||
func transform(_ input: borrowing Buffer<Element>, _ output: borrowing Buffer<Element>, _ work: borrowing Buffer<Element>?, _ dir: pffft_direction_t) {
|
||||
checkFftBufferCounts(n: n, type: type, input: input, output: output, work: work)
|
||||
|
||||
let workAddress: UnsafeMutablePointer<Element>! = 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<Double>, _ output: borrowing Buffer<Double>, _ work: borrowing Buffer<Double>?, _ dir: pffft_direction_t) {
|
||||
checkFftBufferCounts(n: n, type: type, input: input, output: output, work: work)
|
||||
|
||||
let workAddress: UnsafeMutablePointer<Double>! = switch work {
|
||||
case let .some(b): b.baseAddress
|
||||
case .none: nil
|
||||
}
|
||||
pffftd_transform(ptr, input.baseAddress, output.baseAddress, workAddress, dir)
|
||||
}
|
||||
|
||||
|
||||
public func forward(input: borrowing Buffer<Double>, output: borrowing Buffer<Double>, work: borrowing Buffer<Double>?) {
|
||||
transform(input, output, work, PFFFT_FORWARD)
|
||||
}
|
||||
|
||||
public func inverse(input: borrowing Buffer<Double>, output: borrowing Buffer<Double>, work: borrowing Buffer<Double>?) {
|
||||
transform(input, output, work, PFFFT_BACKWARD)
|
||||
}
|
||||
|
||||
public func forwardUnordered(input: borrowing Buffer<Double>, output: borrowing Buffer<Double>, work: borrowing Buffer<Double>?) {
|
||||
transformUnordered(input, output, work, PFFFT_FORWARD)
|
||||
}
|
||||
|
||||
public func inverseUnordered(input: borrowing Buffer<Double>, output: borrowing Buffer<Double>, work: borrowing Buffer<Double>?) {
|
||||
transformUnordered(input, output, work, PFFFT_BACKWARD)
|
||||
}
|
||||
|
||||
|
||||
public func reorderSpectrum(input: borrowing Buffer<Double>, output: borrowing Buffer<Double>) {
|
||||
checkFftBufferCounts(n: n, type: type, input: input, output: output, work: nil)
|
||||
pffftd_zreorder(ptr, input.baseAddress, output.baseAddress, PFFFT_FORWARD)
|
||||
}
|
||||
|
||||
public func convolveAccumulate(dftA: borrowing Buffer<Double>, dftB: borrowing Buffer<Double>, dftAB: borrowing Buffer<Double>, 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 convolve(dftA: borrowing Buffer<Double>, dftB: borrowing Buffer<Double>, dftAB: borrowing Buffer<Double>, 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)
|
||||
}
|
||||
|
||||
public static func minFftSize(for type: FFTType) -> Int {
|
||||
Int(pffftd_min_fft_size(pffft_transform_t(type)))
|
||||
}
|
||||
|
||||
public static func isValidSize(_ n: Int, for type: FFTType) -> Bool {
|
||||
pffftd_is_valid_size(Int32(n), pffft_transform_t(type)) != 0
|
||||
}
|
||||
|
||||
public static func nearestValidSize(_ n: Int, for type: FFTType, higher: Bool) -> Int {
|
||||
Int(pffftd_nearest_transform_size(Int32(n), pffft_transform_t(type), higher ? 1 : 0))
|
||||
}
|
||||
|
||||
deinit {
|
||||
pffftd_destroy_setup(ptr)
|
||||
}
|
||||
}
|
||||
|
||||
extension Double: FFTElement {
|
||||
public typealias FFTImpl = FFTDouble
|
||||
}
|
||||
@@ -1,75 +0,0 @@
|
||||
internal import PFFFTLib
|
||||
|
||||
public final class FFTFloat: PFFFTProtocol {
|
||||
let ptr: OpaquePointer
|
||||
let n: Int
|
||||
let type: FFTType
|
||||
|
||||
static let sharedCache = SetupCache<Element>()
|
||||
|
||||
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
|
||||
self.n = n
|
||||
self.type = type
|
||||
}
|
||||
|
||||
public func forward(input: borrowing Buffer<Float>, output: borrowing Buffer<Float>, work: borrowing Buffer<Float>?, sign: FFTSign) {
|
||||
checkFftBufferCounts(n: n, type: type, input: input, output: output, work: work)
|
||||
|
||||
let workAddress: UnsafeMutablePointer<Float>! = switch work {
|
||||
case let .some(b): b.baseAddress
|
||||
case .none: nil
|
||||
}
|
||||
pffft_transform_ordered(ptr, input.baseAddress, output.baseAddress, workAddress, pffft_direction_t(sign))
|
||||
}
|
||||
|
||||
public func fftUnordered(input: borrowing Buffer<Float>, output: borrowing Buffer<Float>, work: borrowing Buffer<Float>?, sign: FFTSign) {
|
||||
checkFftBufferCounts(n: n, type: type, input: input, output: output, work: work)
|
||||
|
||||
let workAddress: UnsafeMutablePointer<Float>! = switch work {
|
||||
case let .some(b): b.baseAddress
|
||||
case .none: nil
|
||||
}
|
||||
pffft_transform(ptr, input.baseAddress, output.baseAddress, workAddress, pffft_direction_t(sign))
|
||||
}
|
||||
|
||||
public func reorderSpectrum(input: borrowing Buffer<Float>, output: borrowing Buffer<Float>) {
|
||||
checkFftBufferCounts(n: n, type: type, input: input, output: output, work: nil)
|
||||
pffft_zreorder(ptr, input.baseAddress, output.baseAddress, PFFFT_FORWARD)
|
||||
}
|
||||
|
||||
public func convolveAccumulate(dftA: borrowing Buffer<Float>, dftB: borrowing Buffer<Float>, dftAB: borrowing Buffer<Float>, 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 convolve(dftA: borrowing Buffer<Float>, dftB: borrowing Buffer<Float>, dftAB: borrowing Buffer<Float>, 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)
|
||||
}
|
||||
|
||||
public static func minFftSize(for type: FFTType) -> Int {
|
||||
Int(pffft_min_fft_size(pffft_transform_t(type)))
|
||||
}
|
||||
|
||||
public static func isValidSize(_ n: Int, for type: FFTType) -> Bool {
|
||||
pffft_is_valid_size(Int32(n), pffft_transform_t(type)) != 0
|
||||
}
|
||||
|
||||
public static func nearestValidSize(_ n: Int, for type: FFTType, higher: Bool) -> Int {
|
||||
Int(pffft_nearest_transform_size(Int32(n), pffft_transform_t(type), higher ? 1 : 0))
|
||||
}
|
||||
|
||||
deinit {
|
||||
pffft_destroy_setup(ptr)
|
||||
}
|
||||
}
|
||||
|
||||
extension Float: FFTElement {
|
||||
public typealias FFTImpl = FFTFloat
|
||||
}
|
||||
@@ -1,27 +1,12 @@
|
||||
internal import PFFFTLib
|
||||
import ComplexModule
|
||||
import RealModule
|
||||
|
||||
@frozen
|
||||
public enum FFTType {
|
||||
case real
|
||||
case complex
|
||||
}
|
||||
|
||||
@frozen
|
||||
public enum FFTSign: Int {
|
||||
case forward = -1
|
||||
case backward = 1
|
||||
}
|
||||
|
||||
public enum FFTError: Error {
|
||||
case invalidSize
|
||||
}
|
||||
|
||||
public protocol FFTElement {
|
||||
associatedtype FFTImpl: PFFFTProtocol
|
||||
}
|
||||
|
||||
public protocol PFFFTProtocol<Element> {
|
||||
associatedtype Element
|
||||
public protocol PFFFTProtocol<T> {
|
||||
associatedtype T: FFTElemental
|
||||
associatedtype ComplexType = T.ComplexType
|
||||
associatedtype ScalarType = T.ScalarType
|
||||
|
||||
/// Get an FFT interface for the given size and FFT Type.
|
||||
///
|
||||
@@ -31,14 +16,14 @@ public protocol PFFFTProtocol<Element> {
|
||||
/// - 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
|
||||
static func setup(for n: Int) throws -> Self
|
||||
|
||||
/// Initialize the FFT implementation with the given size and type.
|
||||
/// - Parameters:
|
||||
/// - n: The size of the FFT.
|
||||
/// - type: The type of FFT.
|
||||
/// - Throws: `FFTError.invalidSize` if the size is invalid.
|
||||
init(n: Int, type: FFTType) throws
|
||||
init(n: Int) throws
|
||||
|
||||
/// Perform a forward FFT on the input buffer.
|
||||
///
|
||||
@@ -60,9 +45,9 @@ public protocol PFFFTProtocol<Element> {
|
||||
/// - 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 forward(input: borrowing Buffer<Element>, output: borrowing Buffer<Element>, work: borrowing Buffer<Element>?)
|
||||
func forward(input: borrowing Buffer<T>, output: borrowing Buffer<ComplexType>)
|
||||
|
||||
func inverse(input: borrowing Buffer<Element>, output: borrowing Buffer<Element>, work: borrowing Buffer<Element>?)
|
||||
func inverse(input: borrowing Buffer<ComplexType>, output: borrowing Buffer<T>)
|
||||
|
||||
/// Perform a forward FFT on the input buffer, with implementation defined order.
|
||||
///
|
||||
@@ -73,22 +58,11 @@ public protocol PFFFTProtocol<Element> {
|
||||
/// - 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 forwardUnordered(input: borrowing Buffer<Element>, output: borrowing Buffer<Element>, work: borrowing Buffer<Element>?)
|
||||
func forwardToInternalLayout(input: borrowing Buffer<T>, output: borrowing Buffer<ScalarType>)
|
||||
|
||||
func inverseUnordered(input: borrowing Buffer<Element>, output: borrowing Buffer<Element>, work: borrowing Buffer<Element>?)
|
||||
func inverseFromInternalLayout(input: borrowing Buffer<ScalarType>, output: borrowing Buffer<T>)
|
||||
|
||||
func reorder(input: borrowing Buffer<Element>, output: borrowing Buffer<Element>, sign: FFTSign)
|
||||
|
||||
/// Perform a convolution of two complex signals in the frequency domain.
|
||||
///
|
||||
/// Multiplies frequency domain components of `dftA` and `dftB` and accumulates the result in `dftAB`.
|
||||
/// The operation performed is `dftAB += (dftA * dftB) * scaling`.
|
||||
/// - Parameters:
|
||||
/// - dftA: The first input buffer of frequency domain data.
|
||||
/// - 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 convolveAccumulate(dftA: borrowing Buffer<Element>, dftB: borrowing Buffer<Element>, dftAB: borrowing Buffer<Element>, scaling: Element)
|
||||
func reorder(spectrum: borrowing Buffer<ScalarType>, output: borrowing Buffer<ComplexType>)
|
||||
|
||||
/// Perform a convolution of two complex signals in the frequency domain.
|
||||
///
|
||||
@@ -99,12 +73,24 @@ public protocol PFFFTProtocol<Element> {
|
||||
/// - 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 convolve(dftA: borrowing Buffer<Element>, dftB: borrowing Buffer<Element>, dftAB: borrowing Buffer<Element>, scaling: Element)
|
||||
func convolve(dftA: borrowing Buffer<ScalarType>, dftB: borrowing Buffer<ScalarType>, dftAB: borrowing Buffer<ScalarType>, scaling: ScalarType)
|
||||
|
||||
/// Perform a convolution of two complex signals in the frequency domain.
|
||||
///
|
||||
/// Multiplies frequency domain components of `dftA` and `dftB` and accumulates the result in `dftAB`.
|
||||
/// The operation performed is `dftAB += (dftA * dftB) * scaling`.
|
||||
/// - Parameters:
|
||||
/// - dftA: The first input buffer of frequency domain data.
|
||||
/// - 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 convolveAccumulate(dftA: borrowing Buffer<ScalarType>, dftB: borrowing Buffer<ScalarType>, dftAB: borrowing Buffer<ScalarType>, scaling: ScalarType)
|
||||
|
||||
|
||||
/// Returns the minimum FFT size for the given type.
|
||||
///
|
||||
/// - Parameter type: The type of FFT.
|
||||
static func minFftSize(for type: FFTType) -> Int
|
||||
static func minFftSize() -> Int
|
||||
|
||||
/// Returns whether the given size is valid for the given type.
|
||||
///
|
||||
@@ -112,10 +98,10 @@ public protocol PFFFTProtocol<Element> {
|
||||
/// - Parameters:
|
||||
/// - n: The size to check.
|
||||
/// - type: The type of FFT.
|
||||
static func isValidSize(_ n: Int, for type: FFTType) -> Bool
|
||||
static func isValidSize(_ n: Int) -> Bool
|
||||
|
||||
/// Returns the nearest valid size for the given type.
|
||||
static func nearestValidSize(_ n: Int, for type: FFTType, higher: Bool) -> Int
|
||||
static func nearestValidSize(_ n: Int, higher: Bool) -> Int
|
||||
}
|
||||
|
||||
public var simdArch: String {
|
||||
@@ -142,37 +128,3 @@ extension pffft_direction_t {
|
||||
}
|
||||
}
|
||||
|
||||
@inline(__always)
|
||||
func checkFftBufferCounts<Element>(n: Int, type: FFTType, input: borrowing Buffer<Element>, output: borrowing Buffer<Element>, work: borrowing Buffer<Element>?) {
|
||||
let minSize = type == .real ? n : 2 * n
|
||||
|
||||
let work: UnsafeMutableBufferPointer<Element>! = switch work {
|
||||
case let .some(b): b.buffer
|
||||
case .none: nil
|
||||
}
|
||||
|
||||
guard input.count >= minSize else {
|
||||
fatalError("input buffer too small")
|
||||
}
|
||||
guard output.count >= minSize else {
|
||||
fatalError("output buffer too small")
|
||||
}
|
||||
guard work == nil || work.count >= minSize else {
|
||||
fatalError("work buffer too small")
|
||||
}
|
||||
}
|
||||
|
||||
@inline(__always)
|
||||
func checkConvolveBufferCounts<Element>(n: Int, type: FFTType, a: borrowing Buffer<Element>, b: borrowing Buffer<Element>, ab: borrowing Buffer<Element>) {
|
||||
let minSize = type == .real ? n : 2 * n
|
||||
|
||||
guard a.count >= minSize else {
|
||||
fatalError("a buffer too small")
|
||||
}
|
||||
guard b.count >= minSize else {
|
||||
fatalError("b buffer too small")
|
||||
}
|
||||
guard ab.count >= minSize else {
|
||||
fatalError("ab buffer too small")
|
||||
}
|
||||
}
|
||||
|
||||
@@ -1,28 +1,24 @@
|
||||
import Foundation
|
||||
|
||||
struct CacheKey : Hashable {
|
||||
let n: Int
|
||||
let type: FFTType
|
||||
}
|
||||
|
||||
public class SetupCache<Element: FFTElement> : @unchecked Sendable {
|
||||
private var cache: [CacheKey: Element.FFTImpl?] = [:]
|
||||
public class SetupCache<T: FFTElemental> : @unchecked Sendable {
|
||||
typealias FFTType = FFT<T>
|
||||
private var cache: [Int: FFT<T>] = [:]
|
||||
|
||||
private let queue = DispatchQueue(label: String(describing: SetupCache.self), attributes: .concurrent)
|
||||
|
||||
public init() {}
|
||||
|
||||
public func get(for n: Int, type: FFTType) throws -> Element.FFTImpl {
|
||||
var setup: Element.FFTImpl??
|
||||
public func get(for n: Int) throws -> FFT<T> {
|
||||
var setup: FFTType??
|
||||
queue.sync {
|
||||
setup = cache[CacheKey(n: n, type: type)]
|
||||
setup = cache[n]
|
||||
}
|
||||
if setup == nil {
|
||||
queue.sync(flags: .barrier) {
|
||||
setup = cache[CacheKey(n: n, type: type)]
|
||||
setup = cache[n]
|
||||
if setup == nil {
|
||||
let entry = try? Element.FFTImpl(n: n, type: type)
|
||||
cache[CacheKey(n: n, type: type)] = entry
|
||||
let entry = try? FFTType(n: n)
|
||||
cache[n] = entry
|
||||
setup = entry
|
||||
}
|
||||
}
|
||||
|
||||
@@ -1,7 +1,9 @@
|
||||
import Algorithms
|
||||
import ComplexModule
|
||||
import Foundation
|
||||
import PFFFTLib
|
||||
import PFFFT
|
||||
import PFFFTLib
|
||||
import RealModule
|
||||
|
||||
struct Parameters {
|
||||
/// The number of beats to simulate.
|
||||
@@ -63,24 +65,21 @@ func stdev(_ data: [Double]) -> Double {
|
||||
return sqrt(data.lazy.map { ($0 - mean) * ($0 - mean) }.reduce(0.0, +) / (n - 1))
|
||||
}
|
||||
|
||||
struct RRProcess : ~Copyable {
|
||||
struct RRProcess: ~Copyable {
|
||||
let nrr: Int
|
||||
|
||||
let input: Buffer<Double>
|
||||
let output: Buffer<Double>
|
||||
let work: Buffer<Double>
|
||||
let fft: FFTDouble
|
||||
let spectrum: Buffer<Complex<Double>>
|
||||
let signal: Buffer<Double>
|
||||
let fft: FFT<Double>
|
||||
|
||||
init(nrr: Int) {
|
||||
self.nrr = nrr
|
||||
self.input = Buffer<Double>(capacity: nrr + 2)
|
||||
self.output = Buffer<Double>(capacity: nrr)
|
||||
self.work = Buffer<Double>(capacity: nrr)
|
||||
self.fft = try! FFTDouble.setup(for: nrr, type: .real)
|
||||
fft = try! FFT<Double>(n: nrr)
|
||||
spectrum = fft.makeSpectrumBuffer()
|
||||
signal = fft.makeSignalBuffer()
|
||||
}
|
||||
|
||||
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
|
||||
@@ -95,7 +94,7 @@ struct RRProcess : ~Copyable {
|
||||
let sf = Double(params.sfInternal)
|
||||
let df = sf / Double(nrr)
|
||||
|
||||
input.withUnsafeMutableBufferPointer { swc in
|
||||
spectrum.withUnsafeMutableBufferPointer { swc in
|
||||
for i in 0 ..< nrr / 2 + 1 {
|
||||
let w = df * Double(i) * 2.0 * .pi
|
||||
let dw1 = w - w1
|
||||
@@ -106,18 +105,18 @@ struct RRProcess : ~Copyable {
|
||||
let sw = (sf / 2.0) * sqrt(hw)
|
||||
let ph = 2.0 * .pi * Double.random(in: 0.0 ..< 1.0)
|
||||
|
||||
swc[i * 2] = sw * cos(ph)
|
||||
swc[i * 2 + 1] = sw * sin(ph)
|
||||
swc[i].real = sw * cos(ph)
|
||||
swc[i].imaginary = sw * sin(ph)
|
||||
}
|
||||
|
||||
// pack Nyquist frequency real to imaginary of DC
|
||||
swc[1] = swc[nrr]
|
||||
swc[0].imaginary = swc[nrr / 2].real
|
||||
}
|
||||
|
||||
self.fft.forward(input: input, output: output, work: nil, sign: .backward)
|
||||
fft.inverse(spectrum: spectrum, signal: signal)
|
||||
|
||||
var rr = output.withUnsafeMutableBufferPointer { outptr in
|
||||
return outptr.map { $0 * 1.0 / Double(nrr) }
|
||||
var rr = signal.withUnsafeMutableBufferPointer { outptr in
|
||||
outptr.map { $0 * 1.0 / Double(nrr) }
|
||||
}
|
||||
|
||||
let xstd = stdev(rr)
|
||||
@@ -127,11 +126,9 @@ struct RRProcess : ~Copyable {
|
||||
rr[i] = rr[i] * ratio + rrmean
|
||||
}
|
||||
return rr
|
||||
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
// func compute(params: Parameters) {
|
||||
// // adjust extrema parameters for mean heart rate
|
||||
// let hrFact = sqrt(params.hrMean / 60)
|
||||
|
||||
@@ -6,7 +6,8 @@ import Testing
|
||||
let p = Parameters()
|
||||
|
||||
|
||||
let rr = rrprocess(params: p, nrr: 131072)
|
||||
let rrp = RRProcess(nrr: 131072)
|
||||
let rr = rrp.generate(params: p)
|
||||
//fft(data: &a, isign: -1)
|
||||
//print("out: \(rr)")
|
||||
|
||||
|
||||
Reference in New Issue
Block a user