diff --git a/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/index.js b/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/index.js
new file mode 100644
index 000000000000..d4e7f52daa0a
--- /dev/null
+++ b/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/index.js
@@ -0,0 +1,49 @@
+/**
+* @license Apache-2.0
+*
+* Copyright (c) 2026 The Stdlib Authors.
+*
+* Licensed under the Apache License, Version 2.0 (the "License");
+* you may not use this file except in compliance with the License.
+* You may obtain a copy of the License at
+*
+*    http://www.apache.org/licenses/LICENSE-2.0
+*
+* Unless required by applicable law or agreed to in writing, software
+* distributed under the License is distributed on an "AS IS" BASIS,
+* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+* See the License for the specific language governing permissions and
+* limitations under the License.
+*/
+
+'use strict';
+
+/**
+* Compute the forward real-valued fast Fourier transform (FFT).
+*
+* @module @stdlib/fft/base/fftpack/rfftf
+*
+* @example
+* var Float64Array = require( '@stdlib/array/float64' );
+* var rffti = require( '@stdlib/fft/base/fftpack/rffti' );
+* var rfftf = require( '@stdlib/fft/base/fftpack/rfftf' );
+*
+* var N = 4;
+*
+* var workspace = new Float64Array( ( 2*N ) + 34 );
+* rffti( N, workspace, 1, 0 );
+*
+* var r = new Float64Array( [ 1.0, 2.0, 3.0, 4.0 ] );
+*
+* rfftf( N, r, 1, 0, workspace, 1, 0 );
+* // r => <Float64Array>[ 10.0, -2.0, 2.0, -2.0 ]
+*/
+
+// MODULES //
+
+var main = require( './main.js' );
+
+
+// EXPORTS //
+
+module.exports = main;
diff --git a/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/main.js b/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/main.js
new file mode 100644
index 000000000000..0e54338933a1
--- /dev/null
+++ b/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/main.js
@@ -0,0 +1,124 @@
+/**
+* @license Apache-2.0
+*
+* Copyright (c) 2026 The Stdlib Authors.
+*
+* Licensed under the Apache License, Version 2.0 (the "License");
+* you may not use this file except in compliance with the License.
+* You may obtain a copy of the License at
+*
+*    http://www.apache.org/licenses/LICENSE-2.0
+*
+* Unless required by applicable law or agreed to in writing, software
+* distributed under the License is distributed on an "AS IS" BASIS,
+* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+* See the License for the specific language governing permissions and
+* limitations under the License.
+*
+*
+* ## Notice
+*
+* The original C code and copyright notice are from the [PFFFT library]{@link https://github.com/marton78/pffft/blob/master/fftpack.c}. The implementation follows the original, but has been modified for JavaScript.
+*
+* ```text
+* Copyright (c) 2004 the University Corporation for Atmospheric
+* Research ("UCAR"). All rights reserved. Developed by NCAR's
+* Computational and Information Systems Laboratory, UCAR,
+* www.cisl.ucar.edu.
+*
+* Redistribution and use of the Software in source and binary forms,
+* with or without modification, is permitted provided that the
+* following conditions are met:
+*
+*     - Neither the names of NCAR's Computational and Information Systems
+*       Laboratory, the University Corporation for Atmospheric Research,
+*       nor the names of its sponsors or contributors may be used to
+*       endorse or promote products derived from this Software without
+*       specific prior written permission.
+*
+*     - Redistributions of source code must retain the above copyright
+*       notices, this list of conditions, and the disclaimer below.
+*
+*     - Redistributions in binary form must reproduce the above copyright
+*       notice, this list of conditions, and the disclaimer below in the
+*       documentation and/or other materials provided with the
+*       distribution.
+*
+* THIS SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
+* EXPRESS OR IMPLIED, INCLUDING, BUT NOT LIMITED TO THE WARRANTIES OF
+* MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
+* NONINFRINGEMENT. IN NO EVENT SHALL THE CONTRIBUTORS OR COPYRIGHT
+* HOLDERS BE LIABLE FOR ANY CLAIM, INDIRECT, INCIDENTAL, SPECIAL,
+* EXEMPLARY, OR CONSEQUENTIAL DAMAGES OR OTHER LIABILITY, WHETHER IN AN
+* ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
+* CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS WITH THE
+* SOFTWARE.
+* ```
+*/
+
+'use strict';
+
+// MODULES //
+
+var isNonNegativeInteger = require( '@stdlib/assert/is-nonnegative-integer' );
+var isFloat64Array = require( '@stdlib/assert/is-float64array' );
+var isInteger = require( '@stdlib/assert/is-integer' );
+var rfftf1 = require( './rfftf1.js' );
+
+
+// MAIN //
+
+/**
+* Computes the forward real-valued fast Fourier transform (FFT).
+*
+* @param {NonNegativeInteger} N - length of the sequence to transform
+* @param {Float64Array} r - input array containing the sequence to transform
+* @param {integer} strideR - stride length for `r`
+* @param {NonNegativeInteger} offsetR - starting index for `r`
+* @param {Float64Array} workspace - workspace array containing pre-computed values
+* @param {integer} strideW - stride length for `workspace`
+* @param {NonNegativeInteger} offsetW - starting index for `workspace`
+* @returns {void}
+*
+* @example
+* var Float64Array = require( '@stdlib/array/float64' );
+* var rffti = require( '@stdlib/fft/base/fftpack/rffti' );
+*
+* var N = 4;
+*
+* var workspace = new Float64Array( ( 2*N ) + 34 );
+* rffti( N, workspace, 1, 0 );
+*
+* var r = new Float64Array( [ 1.0, 2.0, 3.0, 4.0 ] );
+*
+* rfftf( N, r, 1, 0, workspace, 1, 0 );
+* // r => <Float64Array>[ 10.0, -2.0, 2.0, -2.0 ]
+*/
+function rfftf( N, r, strideR, offsetR, workspace, strideW, offsetW ) {
+	var offsetT;
+	var offsetF;
+
+	if ( !isNonNegativeInteger( N ) || !isFloat64Array( r ) ||
+		!isInteger( strideR ) || !isNonNegativeInteger( offsetR ) ||
+		!isFloat64Array( workspace ) || !isInteger( strideW ) ||
+		!isNonNegativeInteger( offsetW ) ||
+		workspace.length < offsetW + ( ( ( 2*N ) + 34 ) * strideW ) ) {
+		return;
+	}
+
+	// When a sub-sequence is a single data point, the FFT is the identity, so no transformation necessary...
+	if ( N === 1 ) {
+		return;
+	}
+
+	// Resolve the starting indices for storing twiddle factors and factorization results:
+	offsetT = offsetW + (N * strideW); // index offset for twiddle factors
+	offsetF = offsetT + (N * strideW); // index offset for factors describing the sub-transforms
+
+	rfftf1( N, r, offsetR, workspace, offsetW, workspace, offsetT, workspace, offsetF ); // eslint-disable-line max-len
+}
+
+
+// EXPORTS //
+
+module.exports = rfftf;
diff --git a/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/radf2.js b/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/radf2.js
new file mode 100644
index 000000000000..0d0af15bbf79
--- /dev/null
+++ b/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/radf2.js
@@ -0,0 +1,366 @@
+/**
+* @license Apache-2.0
+*
+* Copyright (c) 2026 The Stdlib Authors.
+*
+* Licensed under the Apache License, Version 2.0 (the "License");
+* you may not use this file except in compliance with the License.
+* You may obtain a copy of the License at
+*
+*    http://www.apache.org/licenses/LICENSE-2.0
+*
+* Unless required by applicable law or agreed to in writing, software
+* distributed under the License is distributed on an "AS IS" BASIS,
+* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+* See the License for the specific language governing permissions and
+* limitations under the License.
+*
+*
+* ## Notice
+*
+* The original C code and copyright notice are from the [PFFFT library]{@link https://github.com/marton78/pffft/blob/0b4ee12c4ba45a4a8e567550c16d96d1679f50ce/src/fftpack.c}. The implementation follows the original, but has been modified for JavaScript.
+*
+* ```text
+* Copyright (c) 2004 the University Corporation for Atmospheric
+* Research ("UCAR"). All rights reserved. Developed by NCAR's
+* Computational and Information Systems Laboratory, UCAR,
+* www.cisl.ucar.edu.
+*
+* Redistribution and use of the Software in source and binary forms,
+* with or without modification, is permitted provided that the
+* following conditions are met:
+*
+*     - Neither the names of NCAR's Computational and Information Systems
+*       Laboratory, the University Corporation for Atmospheric Research,
+*       nor the names of its sponsors or contributors may be used to
+*       endorse or promote products derived from this Software without
+*       specific prior written permission.
+*
+*     - Redistributions of source code must retain the above copyright
+*       notices, this list of conditions, and the disclaimer below.
+*
+*     - Redistributions in binary form must reproduce the above copyright
+*       notice, this list of conditions, and the disclaimer below in the
+*       documentation and/or other materials provided with the
+*       distribution.
+*
+* THIS SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
+* EXPRESS OR IMPLIED, INCLUDING, BUT NOT LIMITED TO THE WARRANTIES OF
+* MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
+* NONINFRINGEMENT. IN NO EVENT SHALL THE CONTRIBUTORS OR COPYRIGHT
+* HOLDERS BE LIABLE FOR ANY CLAIM, INDIRECT, INCIDENTAL, SPECIAL,
+* EXEMPLARY, OR CONSEQUENTIAL DAMAGES OR OTHER LIABILITY, WHETHER IN AN
+* ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
+* CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS WITH THE
+* SOFTWARE.
+* ```
+*/
+
+/* eslint-disable max-len */
+
+'use strict';
+
+// FUNCTIONS //
+
+/**
+* Resolves an index into the input array.
+*
+* ## Notes
+*
+* In a forward real FFT, the previous stage writes its results as two "rows" per sub-sequence.
+*
+* Thus, when reading from an input array stored in linear memory, we can reinterpret the array as a three-dimensional logical view containing `L` independent sub-sequences having two "rows" (`even + odd` and `even - odd`, respectively) and where each "row" is arranged as `M*L` contiguous elements corresponding to interleaved real and imaginary components.
+*
+* Accordingly, the following is a logical view of an input array (zero-based indexing) which contains `L = 3` transforms and in which each sub-sequence has length `M = 4`:
+*
+* ```text
+*                  │           k = 0                    k = 1                    k = 2
+*                  │ ──────────────────────────────────────────────────────────────────────────→ k
+* j = 0 (even+odd) │  cc(0,0,0) ... cc(3,0,0)  cc(0,1,0) ... cc(3,1,0)  cc(0,2,0) ... cc(3,2,0)
+*                  │
+* j = 1 (even-odd) │  cc(0,0,1) ... cc(3,0,1)  cc(0,1,1) ... cc(3,1,1)  cc(0,2,1) ... cc(3,2,1)
+*                  └───────────────────────────────────────────────────────────────────────────→ i
+*                         ↑             ↑          ↑             ↑          ↑             ↑
+*                       i = 0          M-1         0            M-1         0            M-1
+* ```
+*
+* In the above,
+*
+* -   `i` is the fastest varying index, which walks within one short sub-sequence corresponding to either the `even + odd` or `even - odd` row.
+* -   `j` selects between the `even + odd` and `even - odd` row.
+* -   `k` specifies the index of one of the `L` independent transforms we are processing.
+*
+* In linear memory, the three-dimensional logical view is arranged as follows:
+*
+* ```text
+* | cc(0,0,0)...cc(3,0,0) ... cc(0,2,0)...cc(3,2,0) | cc(0,0,1)...cc(3,0,1) ... cc(0,2,1)...cc(3,2,1) |
+*       ↑           ↑             ↑           ↑           ↑           ↑             ↑           ↑
+*       0          M-1           LM         LM-1       (L+1)M     (L+1)M-1       (2L-1)M      2LM-1
+* ```
+*
+* @private
+* @param {NonNegativeInteger} i - index of an element within a sub-sequence
+* @param {NonNegativeInteger} k - index of the sub-sequence being transformed
+* @param {NonNegativeInteger} j - input row
+* @param {NonNegativeInteger} L - number of sub-sequences
+* @param {NonNegativeInteger} M - sub-sequence length
+* @param {integer} stride - stride length of the input array
+* @param {NonNegativeInteger} offset - index specifying the first indexed element in the output array
+* @returns {NonNegativeInteger} computed index
+*
+* @example
+* var stride = 1;
+* var offset = 0;
+*
+* var M = 4; // sub-sequence length
+* var L = 3; // number of sub-sequences
+*
+* var idx = iptr( 0, 0, 0, L, M, stride, offset );
+* // returns 0
+*
+* idx = iptr( 1, 0, 0, L, M, stride, offset );
+* // returns 1
+*
+* idx = iptr( M-1, 0, 0, L, M, stride, offset );
+* // returns 3
+*
+* idx = iptr( 0, 1, 0, L, M, stride, offset );
+* // returns 4
+*
+* // ...
+*
+* idx = iptr( M-1, L-1, 1, L, M, stride, offset );
+* // returns 23
+*/
+function iptr( i, k, j, L, M, stride, offset ) {
+	var n = i + ( ( k+(j*L) ) * M );
+	return ( n*stride ) + offset;
+}
+
+/**
+* Resolves an index into the output array.
+*
+* ## Notes
+*
+* When writing to an output array stored in linear memory, we can reinterpret the array as a three-dimensional logical view containing `L` independent sub-sequences having two "columns" corresponding to the even and odd half-spectrum of a folded complex vector (with real and imaginary parts of each half-spectrum interleaved along each sub-sequence) and where each "column" has `M` elements.
+*
+* Accordingly, the following is a logical view of an output array (zero-based indexing) which contains `L = 3` transforms and in which each "column" sub-sequence has length `M = 4`:
+*
+* ```text
+*                 j = 0 ("even" column)              j = 1 ("odd" column)
+* k = 0 ─┬───────────────────────────────────────┬───────────────────────────────────────┐
+*        │ out(0,0,0)  out(1,0,0) ... out(3,0,0) │ out(0,1,0)  out(1,1,0) ... out(3,1,0) │
+*        └───────────────────────────────────────┴───────────────────────────────────────┤
+* k = 1 ─┬───────────────────────────────────────┬───────────────────────────────────────┤
+*        │ out(0,0,1)  out(1,0,1) ... out(3,0,1) │ out(0,1,1)  out(1,1,1) ... out(3,1,1) │
+*        └───────────────────────────────────────┴───────────────────────────────────────┤
+* k = 2 ─┬───────────────────────────────────────┬───────────────────────────────────────┤
+*        │ out(0,0,2)  out(1,0,2) ... out(3,0,2) │ out(0,1,2)  out(1,1,2) ... out(3,1,2) │
+*        └───────────────────────────────────────┴───────────────────────────────────────┘
+*              ↑           ↑              ↑            ↑           ↑              ↑
+*            i = 0         1             M-1           0           1             M-1
+* ```
+*
+* In the above,
+*
+* -   `i` is the fastest varying index, which walks within one short "column" sub-sequence.
+* -   `j` selects between the even (0) and odd (1) column.
+* -   `k` specifies the index of one of the `L` independent transforms we are processing.
+*
+* In linear memory, the three-dimensional logical view is arranged as follows:
+*
+* ```text
+* | out(0,0,0)...out(3,0,0) | out(0,1,0)...out(3,1,0) | out(0,0,1)...out(3,0,1) | ... | out(0,1,2)...out(3,1,2) |
+*       ↑            ↑            ↑            ↑            ↑            ↑                  ↑            ↑
+*       0           M-1           M          2M-1          2M          3M-1              (2L-1)M       2LM-1
+* ```
+*
+* As may be observed, when resolving an index in the output array, the `j` and `k` dimensions are swapped relative index resolution in the input array. This stems from `radf2` being only one stage in a multi-stage driver which alternates between using `cc` and `out` as workspace buffers. After each stage, the next stage reads what the previous stage wrote.
+*
+* Each stage expects a transpose, and, in order to avoid explicit transposition between the stages, we swap the last two logical dimensions while still maintaining cache locality within the inner loop logical dimension, as indexed by `i`.
+*
+* @private
+* @param {NonNegativeInteger} i - index of an element within a sub-sequence
+* @param {NonNegativeInteger} j - index specifying which of the two complex halves we are in (either `0` or `1`)
+* @param {NonNegativeInteger} k - index of the sub-sequence being transformed
+* @param {NonNegativeInteger} M - sub-sequence length
+* @param {integer} stride - stride length of the output array
+* @param {NonNegativeInteger} offset - index specifying the first indexed element in the output array
+* @returns {NonNegativeInteger} computed index
+*
+* @example
+* var stride = 1;
+* var offset = 0;
+*
+* var M = 4; // sub-sequence length
+* var L = 3; // number of sub-sequences
+*
+* var idx = optr( 0, 0, 0, M, stride, offset );
+* // returns 0
+*
+* idx = optr( 1, 0, 0, M, stride, offset );
+* // returns 1
+*
+* idx = optr( M-1, 0, 0, M, stride, offset );
+* // returns 3
+*
+* idx = optr( 0, 1, 0, M, stride, offset );
+* // returns 4
+*
+* // ...
+*
+* idx = optr( M-1, 1, L-1, M, stride, offset );
+* // returns 23
+*/
+function optr( i, j, k, M, stride, offset ) {
+	var n = i + ( ( j+(k*2) ) * M );
+	return ( n*stride ) + offset;
+}
+
+
+// MAIN //
+
+/**
+* Performs one radix-2 stage within a forward Fourier transform for a real-valued sequence.
+*
+* @private
+* @param {NonNegativeInteger} M - number of elements in each sub-sequence to be transformed
+* @param {NonNegativeInteger} L - number of sub-sequences to be transformed
+* @param {Float64Array} cc - input array containing the sub-sequences to be transformed
+* @param {integer} sc - stride length for `cc`
+* @param {NonNegativeInteger} oc - index offset for `cc`
+* @param {Float64Array} out - output array containing transformed sub-sequences
+* @param {integer} so - stride length for `out`
+* @param {NonNegativeInteger} oo - index offset for `out`
+* @param {Float64Array} twiddles - array containing twiddle factors
+* @param {integer} st - stride length for `twiddles`
+* @param {NonNegativeInteger} ot - index offset for `twiddles`
+* @returns {void}
+*/
+function radf2( M, L, cc, sc, oc, out, so, oo, twiddles, st, ot ) { // eslint-disable-line max-params
+	var tr2;
+	var ti2;
+	var ip1;
+	var ip2;
+	var it1;
+	var it2;
+	var io;
+	var im;
+	var i;
+	var k;
+
+	/*
+	* First, perform the core butterfly for each sub-sequence being transformed.
+	*
+	* In the following loop, we handle harmonic `n = 0` for every transform. As described for `iptr`, the input array is interpreted as a three-dimensional array, containing two "rows" per sub-sequence.
+	*
+	*     row0 = even + odd    (j = 0)
+	*     row1 = even - odd    (j = 1)
+	*
+	* For a radix-2 forward FFT of a real-valued signal `x` and `n = 0`,
+	*
+	*     x[0] = row0 + row1 = 2⋅even
+	*     x[M] = row0 - row1 = 2⋅odd
+	*
+	* Because `W_2^0 = 1`, no twiddle multiplication is necessary.
+	*/
+	for ( k = 0; k < L; k++ ) {
+		ip1 = iptr( 0, k, 0, L, M, sc, oc ); // real part row 0
+		ip2 = iptr( 0, k, 1, L, M, sc, oc ); // real part row 1
+
+		io = optr( 0, 0, k, M, so, oo );
+		out[ io ] = cc[ ip1 ] + cc[ ip2 ]; // even + odd
+
+		io = optr( M-1, 1, k, M, so, oo );
+		out[ io ] = cc[ ip1 ] - cc[ ip2 ]; // even - odd
+	}
+	// When the number of elements in a sub-sequence is less than `2`, there is nothing more to do, as the above butterfly produced the full result...
+	if ( M < 2 ) {
+		return;
+	}
+	/*
+	* Next, apply the general case where we need to loop through the non-trivial harmonics.
+	*
+	* For each harmonic `n = 1, ..., M/2-1`, we need to
+	*
+	* -   recover even/odd spectra from the two input rows.
+	* -   multiply the odd part by the twiddle factor `W_n = cos(θ) - j sin(θ)`.
+	* -   form the following
+	*
+	*         x[2n]   = E_n + W_n⋅O_n    => column 0
+	*         x[2n+1] = E_n - W_n⋅O_n    => column 1
+	*
+	* The mirror index `im = M - i` selects the conjugate-symmetric partner, thus allowing the routine to read each symmetry pair only once.
+	*/
+	if ( M >= 3 ) {
+		// Loop over each sub-sequence to be transformed...
+		for ( k = 0; k < L; k++ ) {
+			// Loop over the elements in each sub-sequence...
+			for ( i = 2; i < M; i += 2 ) {
+				im = M - i; // "mirror" index
+
+				// Resolve twiddle factor indices:
+				it1 = ( (i-2)*st ) + ot; // cos(θ)
+				it2 = ( (i-1)*st ) + ot; // sin(θ)
+
+				// Load the `even-odd` row...
+				ip1 = iptr( i-1, k, 1, L, M, sc, oc );   // c = Re(row1)
+				ip2 = iptr( i, k, 1, L, M, sc, oc ); // d = Im(row1)
+
+				// eslint-disable-next-line stdlib/capitalized-comments, stdlib/empty-line-before-comment
+				// tmp = W_n ⋅ (c + j⋅d)
+				tr2 = ( twiddles[ it1 ] * cc[ ip1 ] ) + ( twiddles[ it2 ] * cc[ ip2 ] ); // Re(tmp)
+				ti2 = ( twiddles[ it1 ] * cc[ ip2 ] ) - ( twiddles[ it2 ] * cc[ ip1 ] ); // Im(tmp)
+
+				// Load the `even+odd` row...
+				ip1 = iptr( i, k, 0, L, M, sc, oc ); // b = Im(row0)
+				io = optr( i, 0, k, M, so, oo );     // Im(x[2n])
+				out[ io ] = cc[ ip1 ] + ti2;
+
+				io = optr( im, 1, k, M, so, oo );      // Im(x[2n+1])
+				out[ io ] = ti2 - cc[ ip1 ];
+
+				ip1 = iptr( i-1, k, 0, L, M, sc, oc );   // a = Re(row0)
+				io = optr( i-1, 0, k, M, so, oo );       // Re(x[2n])
+				out[ io ] = cc[ ip1 ] + tr2;
+
+				io = optr( im-1, 1, k, M, so, oo );    // Re(x[2n+1])
+				out[ io ] = cc[ ip1 ] - tr2;
+			}
+		}
+		// When `M` is odd, there is no Nyquist pair to process, so we do not need to perform any further transformations...
+		if ( M%2 === 1 ) {
+			return;
+		}
+	}
+	/*
+	* Lastly, handle the Nyquist frequency where `i = M-1` (i.e., the last element of each sub-sequence).
+	*
+	* When `M` is even, the Nyquist index is `i = M/2`. In this stage, we've stored that element at the end of each sub-sequence (i.e., `i = M-1` in the packed layout).
+	*
+	* At this point, the cosine term is -1 and the sine term is 0, so the twiddle multiplication collapses to simple addition/subtraction.
+	*
+	* In this case,
+	*
+	*     W_n⋅O_n = ±Re(O_n)
+	*
+	* where
+	*
+	*     row0 =  Re(E_n)
+	*     row1 = -Re(O_n)
+	*/
+	for ( k = 0; k < L; k++ ) {
+		ip1 = iptr( M-1, k, 1, L, M, sc, oc ); // -Re(O_n)
+		io = optr( 0, 1, k, M, so, oo );
+		out[ io ] = -cc[ ip1 ];                // -(-Re(O_n)) = Re(O_n)
+
+		ip1 = iptr( M-1, k, 0, L, M, sc, oc ); // Re(E_n)
+		io = optr( M-1, 0, k, M, so, oo );
+		out[ io ] = cc[ ip1 ];
+	}
+}
+
+
+// EXPORTS //
+
+module.exports = radf2;
diff --git a/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/radf3.js b/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/radf3.js
new file mode 100644
index 000000000000..2fced4b2c683
--- /dev/null
+++ b/lib/node_modules/@stdlib/fft/base/fftpack/rfftf/lib/radf3.js
@@ -0,0 +1,397 @@
+/**
+* @license Apache-2.0
+*
+* Copyright (c) 2026 The Stdlib Authors.
+*
+* Licensed under the Apache License, Version 2.0 (the "License");
+* you may not use this file except in compliance with the License.
+* You may obtain a copy of the License at
+*
+*    http://www.apache.org/licenses/LICENSE-2.0
+*
+* Unless required by applicable law or agreed to in writing, software
+* distributed under the License is distributed on an "AS IS" BASIS,
+* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+* See the License for the specific language governing permissions and
+* limitations under the License.
+*
+*
+* ## Notice
+*
+* The original C code and copyright notice are from the [PFFFT library]{@link https://github.com/marton78/pffft/blob/0b4ee12c4ba45a4a8e567550c16d96d1679f50ce/src/fftpack.c}. The implementation follows the original, but has been modified for JavaScript.
+*
+* ```text
+* Copyright (c) 2004 the University Corporation for Atmospheric
+* Research ("UCAR"). All rights reserved. Developed by NCAR's
+* Computational and Information Systems Laboratory, UCAR,
+* www.cisl.ucar.edu.
+*
+* Redistribution and use of the Software in source and binary forms,
+* with or without modification, is permitted provided that the
+* following conditions are met:
+*
+*     - Neither the names of NCAR's Computational and Information Systems
+*       Laboratory, the University Corporation for Atmospheric Research,
+*       nor the names of its sponsors or contributors may be used to
+*       endorse or promote products derived from this Software without
+*       specific prior written permission.
+*
+*     - Redistributions of source code must retain the above copyright
+*       notices, this list of conditions, and the disclaimer below.
+*
+*     - Redistributions in binary form must reproduce the above copyright
+*       notice, this list of conditions, and the disclaimer below in the
+*       documentation and/or other materials provided with the
+*       distribution.
+*
+* THIS SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
+* EXPRESS OR IMPLIED, INCLUDING, BUT NOT LIMITED TO THE WARRANTIES OF
+* MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
+* NONINFRINGEMENT. IN NO EVENT SHALL THE CONTRIBUTORS OR COPYRIGHT
+* HOLDERS BE LIABLE FOR ANY CLAIM, INDIRECT, INCIDENTAL, SPECIAL,
+* EXEMPLARY, OR CONSEQUENTIAL DAMAGES OR OTHER LIABILITY, WHETHER IN AN
+* ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
+* CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS WITH THE
+* SOFTWARE.
+* ```
+*/
+
+/* eslint-disable max-len */
+
+'use strict';
+
+// VARIABLES //
+
+var TAUR = -0.5; // cos( 2π/3 )
+var TAUI = 0.866025403784439; // sin( 2π/3 )
+
+
+// FUNCTIONS //
+
+/**
+* Resolves an index into the input array.
+*
+* ## Notes
+*
+* In a forward real FFT, the previous stage writes its results as three "rows" per sub-sequence.
+*
+* Thus, when reading from an input array stored in linear memory, we can reinterpret the array as a three-dimensional logical view containing `L` independent sub-sequences having three "rows" (components 0, 1, and 2) and where each "row" is arranged as `M*L` contiguous elements corresponding to interleaved real and imaginary components.
+*
+* Accordingly, the following is a logical view of an input array (zero-based indexing) which contains `L = 3` transforms and in which each sub-sequence has length `M = 4`:
+*
+* ```text
+*       │           k = 0                    k = 1                    k = 2
+*       │ ──────────────────────────────────────────────────────────────────────────→ k
+* j = 0 │  cc(0,0,0) ... cc(3,0,0)  cc(0,1,0) ... cc(3,1,0)  cc(0,2,0) ... cc(3,2,0)
+*       │
+* j = 1 │  cc(0,0,1) ... cc(3,0,1)  cc(0,1,1) ... cc(3,1,1)  cc(0,2,1) ... cc(3,2,1)
+*       │
+* j = 2 │  cc(0,0,2) ... cc(3,0,2)  cc(0,1,2) ... cc(3,1,2)  cc(0,2,2) ... cc(3,2,2)
+*       └───────────────────────────────────────────────────────────────────────────→ i
+*              ↑             ↑          ↑             ↑          ↑             ↑
+*            i = 0          M-1         0            M-1         0            M-1
+* ```
+*
+* In the above,
+*
+* -   `i` is the fastest varying index, which walks within one short sub-sequence corresponding to one of the three component rows.
+* -   `j` selects between the three component rows (0, 1, or 2).
+* -   `k` specifies the index of one of the `L` independent transforms we are processing.
+*
+* In linear memory, the three-dimensional logical view is arranged as follows:
+*
+* ```text
+* | cc(0,0,0)...cc(3,0,0) ... cc(0,2,0)...cc(3,2,0) | cc(0,0,1)...cc(3,0,1) ... cc(0,2,1)...cc(3,2,1) | cc(0,0,2)...cc(3,0,2) ... cc(0,2,2)...cc(3,2,2) |
+*       ↑           ↑             ↑           ↑           ↑           ↑             ↑           ↑           ↑           ↑             ↑           ↑
+*       0          M-1           LM         LM-1       (L+1)M     (L+1)M-1       (2L-1)M      2LM-1        2LM      (2L+1)M-1       (3L-1)M      3LM-1
+* ```
+*
+* @private
+* @param {NonNegativeInteger} i - index of an element within a sub-sequence
+* @param {NonNegativeInteger} k - index of the sub-sequence being transformed
+* @param {NonNegativeInteger} j - input row
+* @param {NonNegativeInteger} L - number of sub-sequences
+* @param {NonNegativeInteger} M - sub-sequence length
+* @param {integer} stride - stride length of the input array
+* @param {NonNegativeInteger} offset - index specifying the first indexed element in the output array
+* @returns {NonNegativeInteger} computed index
+*
+* @example
+* var stride = 1;
+* var offset = 0;
+*
+* var M = 4; // sub-sequence length
+* var L = 3; // number of sub-sequences
+*
+* var idx = iptr( 0, 0, 0, L, M, stride, offset );
+* // returns 0
+*
+* idx = iptr( 1, 0, 0, L, M, stride, offset );
+* // returns 1
+*
+* idx = iptr( M-1, 0, 0, L, M, stride, offset );
+* // returns 3
+*
+* idx = iptr( 0, 1, 0, L, M, stride, offset );
+* // returns 4
+*
+* // ...
+*
+* idx = iptr( M-1, L-1, 1, L, M, stride, offset );
+* // returns 23
+*/
+function iptr( i, k, j, L, M, stride, offset ) {
+	var n = i + ( ( k+(j*L) ) * M );
+	return ( n*stride ) + offset;
+}
+
+/**
+* Resolves an index into the output array.
+*
+* ## Notes
+*
+* When writing to an output array stored in linear memory, we can reinterpret the array as a three-dimensional logical view containing `L` independent sub-sequences having three "columns" corresponding to the three components of a radix-3 stage (with real and imaginary parts of each component interleaved along each sub-sequence) and where each "column" has `M` elements.
+*
+* Accordingly, the following is a logical view of an output array (zero-based indexing) which contains `L = 3` transforms and in which each "column" sub-sequence has length `M = 4`:
+*
+* ```text
+*                 j = 0 (component 0)                     j = 1 (component 1)                    j = 2 (component 2)
+* k = 0 ─┬───────────────────────────────────────┬───────────────────────────────────────┬───────────────────────────────────────┐
+*        │ out(0,0,0)  out(1,0,0) ... out(3,0,0) │ out(0,1,0)  out(1,1,0) ... out(3,1,0) │ out(0,2,0)  out(1,2,0) ... out(3,2,0) │
+*        └───────────────────────────────────────┴───────────────────────────────────────┴───────────────────────────────────────┤
+* k = 1 ─┬───────────────────────────────────────┬───────────────────────────────────────┬───────────────────────────────────────┤
+*        │ out(0,0,1)  out(1,0,1) ... out(3,0,1) │ out(0,1,1)  out(1,1,1) ... out(3,1,1) │ out(0,2,1)  out(1,2,1) ... out(3,2,1) │
+*        └───────────────────────────────────────┴───────────────────────────────────────┴───────────────────────────────────────┤
+* k = 2 ─┬───────────────────────────────────────┬───────────────────────────────────────┬───────────────────────────────────────┤
+*        │ out(0,0,2)  out(1,0,2) ... out(3,0,2) │ out(0,1,2)  out(1,1,2) ... out(3,1,2) │ out(0,2,2)  out(1,2,2) ... out(3,2,2) │
+*        └───────────────────────────────────────┴───────────────────────────────────────┴───────────────────────────────────────┘
+*              ↑           ↑              ↑            ↑           ↑              ↑            ↑           ↑              ↑
+*            i = 0         1             M-1           0           1             M-1           0           1             M-1
+* ```
+*
+* In the above,
+*
+* -   `i` is the fastest varying index, which walks within one short "column" sub-sequence.
+* -   `j` selects between one of the three components (0, 1, or 2).
+* -   `k` specifies the index of one of the `L` independent transforms we are processing.
+*
+* In linear memory, the three-dimensional logical view is arranged as follows:
+*
+* ```text
+* | out(0,0,0)...out(3,0,0) | out(0,1,0)...out(3,1,0) | out(0,2,0)...out(3,2,0) | out(0,0,1)...out(3,0,1) | ... | out(0,2,2)...out(3,2,2) |
+*       ↑            ↑            ↑            ↑            ↑            ↑            ↑            ↑                  ↑            ↑
+*       0           M-1           M          2M-1          2M          3M-1          3M         4M-1              (3L-1)M       3LM-1
+* ```
+*
+* As may be observed, when resolving an index in the output array, the `j` and `k` dimensions are swapped relative index resolution in the input array. This stems from `radf3` being only one stage in a multi-stage driver which alternates between using `cc` and `out` as workspace buffers. After each stage, the next stage reads what the previous stage wrote.
+*
+* Each stage expects a transpose, and, in order to avoid explicit transposition between the stages, we swap the last two logical dimensions while still maintaining cache locality within the inner loop logical dimension, as indexed by `i`.
+*
+* @private
+* @param {NonNegativeInteger} i - index of an element within a sub-sequence
+* @param {NonNegativeInteger} j - index specifying which of the three complex components we are in (0, 1, or 2)
+* @param {NonNegativeInteger} k - index of the sub-sequence being transformed
+* @param {NonNegativeInteger} M - sub-sequence length
+* @param {integer} stride - stride length of the output array
+* @param {NonNegativeInteger} offset - index specifying the first indexed element in the output array
+* @returns {NonNegativeInteger} computed index
+*
+* @example
+* var stride = 1;
+* var offset = 0;
+*
+* var M = 4; // sub-sequence length
+* var L = 3; // number of sub-sequences
+*
+* var idx = optr( 0, 0, 0, M, stride, offset );
+* // returns 0
+*
+* idx = optr( 1, 0, 0, M, stride, offset );
+* // returns 1
+*
+* idx = optr( M-1, 0, 0, M, stride, offset );
+* // returns 3
+*
+* idx = optr( 0, 1, 0, M, stride, offset );
+* // returns 4
+*
+* // ...
+*
+* idx = optr( M-1, 2, L-1, M, stride, offset );
+* // returns 35
+*/
+function optr( i, j, k, M, stride, offset ) {
+	var n = i + ( ( j+(k*3) ) * M );
+	return ( n*stride ) + offset;
+}
+
+
+// MAIN //
+
+/**
+* Performs one radix-3 stage within a forward Fourier transform for a real-valued sequence.
+*
+* @private
+* @param {NonNegativeInteger} M - number of elements in each sub-sequence to be transformed
+* @param {NonNegativeInteger} L - number of sub-sequences to be transformed
+* @param {Float64Array} cc - input array containing the sub-sequences to be transformed
+* @param {integer} sc - stride length for `cc`
+* @param {NonNegativeInteger} oc - index offset for `cc`
+* @param {Float64Array} out - output array containing transformed sub-sequences
+* @param {integer} so - stride length for `out`
+* @param {NonNegativeInteger} oo - index offset for `out`
+* @param {Float64Array} twiddles1 - first array containing twiddle factors
+* @param {integer} st1 - stride length for `twiddles1`
+* @param {NonNegativeInteger} ot1 - index offset for `twiddles1`
+* @param {Float64Array} twiddles2 - second array containing twiddle factors
+* @param {integer} st2 - stride length for `twiddles2`
+* @param {NonNegativeInteger} ot2 - index offset for `twiddles2`
+* @returns {void}
+*/
+function radf3( M, L, cc, sc, oc, out, so, oo, twiddles1, st1, ot1, twiddles2, st2, ot2 ) { // eslint-disable-line max-params
+	var it11;
+	var it12;
+	var it21;
+	var it22;
+	var tr2;
+	var ti2;
+	var tr3;
+	var ti3;
+	var cr2;
+	var ci2;
+	var dr2;
+	var di2;
+	var dr3;
+	var di3;
+	var ip1;
+	var ip2;
+	var ip3;
+	var io;
+	var im;
+	var i;
+	var k;
+
+	/*
+	* First, perform the core butterfly for each sub-sequence being transformed.
+	*
+	* In the following loop, we handle harmonic `n = 0` for every transform. As described for `iptr`, the input array is interpreted as a three-dimensional array, containing three "rows" per sub-sequence.
+	*
+	*     row0 = input component 0    (j = 0)
+	*     row1 = input component 1    (j = 1)
+	*     row2 = input component 2    (j = 2)
+	*
+	* For a radix-3 forward FFT of a real-valued signal `x` and `n = 0`,
+	*
+	*     x[0] = row0 + row1 + row2
+	*     x[2M] = taui * ( row2-row1 )
+	*     x[M] = row0 + ( taur * ( row1+row2 ) )
+	*
+	* where `taur = cos(2π/3)` and `taui = sin(2π/3)`.
+	*
+	* Because `W_3^0 = 1`, no twiddle multiplication is necessary beyond these constant factors (`taur` and `taui`).
+	*/
+	for ( k = 0; k < L; k++ ) {
+		ip1 = iptr( 0, k, 0, L, M, sc, oc ); // real part row 0
+		ip2 = iptr( 0, k, 1, L, M, sc, oc ); // real part row 1
+		ip3 = iptr( 0, k, 2, L, M, sc, oc ); // real part row 2
+
+		cr2 = cc[ ip2 ] + cc[ ip3 ];
+
+		io = optr( 0, 0, k, M, so, oo );
+		out[ io ] = cc[ ip1 ] + cr2; // row0 + row1 + row2
+
+		io = optr( 0, 2, k, M, so, oo );
+		out[ io ] = TAUI * ( cc[ ip3 ] - cc[ ip2 ] ); // taui * ( row2-row1 )
+
+		io = optr( M-1, 1, k, M, so, oo );
+		out[ io ] = cc[ ip1 ] + ( TAUR * cr2 ); // row0 + ( taur * ( row1+row2 ) )
+	}
+	// When the number of elements in a sub-sequence is less than `2`, there is nothing more to do, as the above butterfly produced the full result...
+	if ( M < 2 ) {
+		return;
+	}
+	/*
+	* Next, apply the general case where we need to loop through the non-trivial harmonics.
+	*
+	* For each harmonic `n = 1, ..., M/2-1`, we need to
+	*
+	* -   recover three spectra from the three input rows.
+	* -   apply radix-3 twiddle factors to rows 1 and 2, then combine with row 0 to form three output columns.
+	* -   form the following
+	*
+	*         x[3n]   = X₀ + (W₁⋅X₁ + W₂⋅X₂)                                         => column 0
+	*         x[3n+1] = X₀ + taur⋅(W₁⋅X₁ + W₂⋅X₂) + taui⋅i⋅(W₁⋅X₁ - W₂⋅X₂)           => column 2
+	*         x[3n+2] = X₀ + taur⋅(W₁⋅X₁ + W₂⋅X₂) - taui⋅i⋅(W₁⋅X₁ - W₂⋅X₂)           => column 1
+	*
+	* The mirror index `im = M - i` selects the conjugate-symmetric partner, thus allowing the routine to read each symmetry pair only once.
+	*/
+	// Loop over each sub-sequence to be transformed...
+	for ( k = 0; k < L; k++ ) {
+		// Loop over the elements in each sub-sequence...
+		for ( i = 2; i < M; i += 2 ) {
+			im = M - i; // "mirror" index
+
+			// Resolve twiddle factor indices for component 1:
+			it11 = ( (i-2)*st1 ) + ot1; // cos(θ) for component 1
+			it12 = ( (i-1)*st1 ) + ot1; // sin(θ) for component 1
+
+			// Resolve twiddle factor indices for component 2:
+			it21 = ( (i-2)*st2 ) + ot2; // cos(θ) for component 2
+			it22 = ( (i-1)*st2 ) + ot2; // sin(θ) for component 2
+
+			// Load component 1 data...
+			ip1 = iptr( i-1, k, 1, L, M, sc, oc ); // real part component 1
+			ip2 = iptr( i, k, 1, L, M, sc, oc );   // imag part component 1
+
+			// Apply twiddle factor for component 1
+			dr2 = ( twiddles1[ it11 ] * cc[ ip1 ] ) + ( twiddles1[ it12 ] * cc[ ip2 ] );
+			di2 = ( twiddles1[ it11 ] * cc[ ip2 ] ) - ( twiddles1[ it12 ] * cc[ ip1 ] );
+
+			// Load component 2 data...
+			ip1 = iptr( i-1, k, 2, L, M, sc, oc ); // real part component 2
+			ip2 = iptr( i, k, 2, L, M, sc, oc );   // imag part component 2
+
+			// Apply twiddle factor for component 2
+			dr3 = ( twiddles2[ it21 ] * cc[ ip1 ] ) + ( twiddles2[ it22 ] * cc[ ip2 ] );
+			di3 = ( twiddles2[ it21 ] * cc[ ip2 ] ) - ( twiddles2[ it22 ] * cc[ ip1 ] );
+
+			// Combine the twiddled components
+			cr2 = dr2 + dr3; // sum of real parts
+			ci2 = di2 + di3; // sum of imag parts
+
+			// Load component 0 data...
+			ip1 = iptr( i-1, k, 0, L, M, sc, oc ); // real part component 0
+			io = optr( i-1, 0, k, M, so, oo );
+			out[ io ] = cc[ ip1 ] + cr2;
+
+			ip2 = iptr( i, k, 0, L, M, sc, oc );   // imag part component 0
+			io = optr( i, 0, k, M, so, oo );
+			out[ io ] = cc[ ip2 ] + ci2;
+
+			// Intermediate results for components 1 and 2
+			tr2 = cc[ ip1 ] + ( TAUR * cr2 );
+			ti2 = cc[ ip2 ] + ( TAUR * ci2 );
+			tr3 = TAUI * ( di2 - di3 );
+			ti3 = TAUI * ( dr3 - dr2 );
+
+			// Output component 1
+			io = optr( i-1, 2, k, M, so, oo );
+			out[ io ] = tr2 + tr3;
+
+			io = optr( im-1, 1, k, M, so, oo );
+			out[ io ] = tr2 - tr3;
+
+			// Output component 2
+			io = optr( i, 2, k, M, so, oo );
+			out[ io ] = ti2 + ti3;
+
+			io = optr( im, 1, k, M, so, oo );
+			out[ io ] = ti3 - ti2;
+		}
+	}
+}
+
+
+// EXPORTS //
+
+module.exports = radf3;