refactor: revert to old implementation

---
type: pre_commit_static_analysis_report
description: Results of running static analysis checks when committing changes.
report:
  - task: lint_filenames
    status: passed
  - task: lint_editorconfig
    status: passed
  - task: lint_markdown
    status: na
  - task: lint_package_json
    status: na
  - task: lint_repl_help
    status: na
  - task: lint_javascript_src
    status: passed
  - task: lint_javascript_cli
    status: na
  - task: lint_javascript_examples
    status: na
  - task: lint_javascript_tests
    status: na
  - task: lint_javascript_benchmarks
    status: na
  - task: lint_python
    status: na
  - task: lint_r
    status: na
  - task: lint_c_src
    status: na
  - task: lint_c_examples
    status: na
  - task: lint_c_benchmarks
    status: na
  - task: lint_c_tests_fixtures
    status: na
  - task: lint_shell
    status: na
  - task: lint_typescript_declarations
    status: passed
  - task: lint_typescript_tests
    status: na
  - task: lint_license_headers
    status: passed
---

Author: gunjjoshi

lib/node_modules/@stdlib/fft/base/fftpack/rffti/lib/rffti1.js

@@ -60,15 +60,15 @@
 
 // MODULES //
 
+var sincos = require( '@stdlib/math/base/special/sincos' ).assign;
 var TWO_PI = require( '@stdlib/constants/float64/two-pi' );
 var floor = require( '@stdlib/math/base/special/floor' );
-var sin = require( '@stdlib/math/base/special/sin' );
-var cos = require( '@stdlib/math/base/special/cos' );
 var decompose = require( '@stdlib/fft/base/fftpack/decompose' );
 
 
 // VARIABLES //
 
+// Define a list of initial trial factors for FFT factorization:
 var TRIAL_FACTORS = [ 4, 2, 3, 5 ];
 
 
@@ -78,104 +78,102 @@ var TRIAL_FACTORS = [ 4, 2, 3, 5 ];
 * Initializes the working arrays for applying a Fourier transform to a real-valued sequence.
 *
 * @private
-* @param {NonNegativeInteger} n - length of the sequence
-* @param {Collection} wa - array of twiddle factors
-* @param {integer} strideW - stride length for `wa`
-* @param {NonNegativeInteger} offsetW - starting index for `wa`
-* @param {Collection} ifac - array containing a radix factorization
-* @param {integer} strideF - stride length for `ifac`
-* @param {NonNegativeInteger} offsetF - starting index for `ifac`
+* @param {NonNegativeInteger} N - length of the sequence
+* @param {Collection} twiddles - array of twiddle factors
+* @param {integer} strideT - stride length for `twiddles`
+* @param {NonNegativeInteger} offsetT - starting index for `twiddles`
+* @param {Collection} factors - array containing a radix factorization
+* @param {integer} strideF - stride length for `factors`
+* @param {NonNegativeInteger} offsetF - starting index for `factors`
 * @returns {void}
 */
-function rffti1( n, wa, strideW, offsetW, ifac, strideF, offsetF ) {
+function rffti1( N, twiddles, strideT, offsetT, factors, strideF, offsetF ) {
+	var factor;
 	var argld;
 	var argh;
-	var nfm1;
-	var arg;
-	var ipm;
-	var ido;
+	var fidx;
 	var nf;
 	var fi;
-	var ip;
 	var l2;
 	var l1;
 	var ld;
-	var ii;
-	var is;
-	var k1;
-	var i;
+	var st;
+	var it;
+	var im;
+	var M;
+	var m;
+	var k;
 	var j;
 
 	// Decompose the sequence length into its radix factors:
-	nf = decompose( n, 4, TRIAL_FACTORS, 1, 0, ifac, strideF, offsetF );
+	nf = decompose( N, 4, TRIAL_FACTORS, 1, 0, factors, strideF, offsetF );
 
+	// If the number of radix factors is `1`, the only twiddle factor we need is `W_N^0 = 1`, which the main transform kernels already hard-code, so nothing to pre-compute...
+	if ( nf-1 === 0 ) {
+		return;
+	}
 	// Compute the master angular step (i.e., the basic angle that generates all twiddles):
-	argh = TWO_PI / n;
-
-	// Initialize the index into the twiddle factor array:
-	is = 0;
+	argh = TWO_PI / N;
 
-	// Compute the number of factors minus one:
-	nfm1 = nf - 1;
+	// Define the location of the first sine term we want to compute:
+	im = 1;
 
 	// Initialize a running product of factors already processed:
 	l1 = 1;
 
-	// If the number of radix factors is `1`, the only twiddle factor we need is `W_N^0 = 1`, which the main transform kernels already hard-code, so nothing to pre-compute...
-	if ( nfm1 === 0 ) {
-		return;
-	}
-	// Generate twiddle factors for each radix factor...
-	for ( k1 = 0; k1 < nfm1; k1++ ) {
-		// Resolve the next radix factor (C code: ifac[k1+2] with 1-based, JS: ifac[k1+2] with 0-based):
-		ip = ifac[ offsetF + ((k1+2)*strideF) ];
+	// Compute the index of the first radix factor:
+	fidx = offsetF + ( 2*strideF );
 
-		// Initialize a running offset used to step the angle:
-		ld = 0;
+	// Compute the stride length for each cosine/sine pair:
+	st = 2 * strideT;
+
+	// Generate twiddle factors...
+	for ( k = 0; k < nf-1; k++ ) {
+		// Resolve the next radix factor:
+		factor = factors[ fidx ];
 
 		// Compute the length of the transform after including the current radix:
-		l2 = l1 * ip;
+		l2 = factor * l1;
 
-		// Compute the number of data points in each sub-transform:
-		ido = floor( n / l2 );
+		// Compute the number of the "butterfly wings" (at this stage, the data is viewed as a `factor`-point transform of sub-vectors of length `l1`; `M` describes how many such vectors fit in the full array of length `N`):
+		M = floor( N / l2 );
 
-		// Compute the number of iterations for the inner loop:
-		ipm = ip - 1;
+		// Initialize a running offset used to step the angle:
+		ld = 0;
 
 		// Iterate over each butterfly column within the current radix...
-		for ( j = 0; j < ipm; j++ ) {
+		for ( j = 1; j < factor; j++ ) {
 			// Advance to the next column:
 			ld += l1;
 
-			// Initialize the index into the twiddle factor array:
-			i = is;
-
 			// Compute the angle step for this column:
-			argld = ld * argh;
+			argld = ld * argh; // 2π ⋅ j ⋅ li / N, which is the j-th column's base angle
 
-			// Initialize the harmonic counter:
-			fi = 0.0;
+			// Initialize a counter which counts from `1` to `(M/2)-1` (i.e., all non-trivial harmonics):
+			fi = 1.0;
 
-			// Iterate over each non-trivial harmonic in the column...
-			for ( ii = 2; ii < ido; ii += 2 ) {
-				// Advance the index by 2 (for cosine and sine):
-				i += 2;
+			// Compute the index of the sine term:
+			it = offsetT + ( (im+1)*strideT );
 
-				// Update the harmonic counter:
-				fi += 1.0;
+			// Iterate over each non-trivial harmonic in the column (note: we skip the `m=0` and `m=1` (i.e., the first cosine/sine pair, the first harmonic `W_N^0 = cos(0) + j sin(0) = 1 + j⋅0` is trivial and hard-coded by transform kernels)...
+			for ( m = 2; m < M; m += 2 ) {
+				// Compute `(cos(fi*argld), sin(fi*argld))`:
+				sincos( fi*argld, twiddles, -strideT, it ); // note: `sincos` returns the sine and cosine, in that order, so we need to use a negative stride when assigning the values to `twiddles` to ensure that cosine comes before sine in the twiddles table
 
-				// Compute the angle for this harmonic:
-				arg = fi * argld;
+				// Update for the next harmonic:
+				fi += 1.0;
 
-				// Store the cosine and sine values:
-				wa[ offsetW + ( ( i-1 ) * strideW ) ] = cos( arg );
-				wa[ offsetW + ( i*strideW ) ] = sin( arg );
+				// Resolve the index of the next sine term:
+				it += st;
 			}
-			// Advance the index by the number of data points in each sub-transform:
-			is += ido;
+			// Advance `im` by the size of the current block, skipping the cosine/sine pairs we have just written:
+			im += M;
 		}
 		// Update the running product of factors already processed:
 		l1 = l2;
+
+		// Resolve the index of the next radix factor:
+		fidx += strideF;
 	}
 }