test: add test.rfftf.js

Author: gunjjoshi

lib/node_modules/@stdlib/fft/base/fftpack/test/test.rfftf.js

@@ -0,0 +1,166 @@
+/**
+* @license Apache-2.0
+*
+* Copyright (c) 2025 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';
+
+// MODULES //
+
+var tape = require( 'tape' );
+var Float64Array = require( '@stdlib/array/float64' );
+var abs = require( '@stdlib/math/base/special/abs' );
+var sin = require( '@stdlib/math/base/special/sin' );
+var cos = require( '@stdlib/math/base/special/cos' );
+var sqrt = require( '@stdlib/math/base/special/sqrt' );
+var floor = require( '@stdlib/math/base/special/floor' );
+var max = require( '@stdlib/math/base/special/max' );
+var PI = require( '@stdlib/constants/float64/pi' );
+var rffti = require( './../lib/rffti.js' );
+var rfftf = require( './../lib/rfftf.js' );
+var rfftb = require( './../lib/rfftb.js' );
+
+
+// FIXTURES //
+
+var nd = [ 120, 91, 54, 49, 32, 28, 24, 8, 4, 3, 2 ];
+
+
+// TESTS //
+
+tape( 'main export is a function', function test( t ) {
+	t.ok( true, __filename );
+	t.strictEqual( typeof rfftf, 'function', 'main export is a function' );
+	t.end();
+});
+
+tape( 'rfftf computes the forward transform with accuracy', function test( t ) {
+	var tolerance = 1.0e-7;
+	var maxError = 0.0;
+	var sumReal = 0.0;
+	var sumImag = 0.0;
+	var sqrt2 = sqrt( 2.0 );
+	var scale = 0.0;
+	var angle = 0.0;
+	var sum1 = 0.0;
+	var sum2 = 0.0;
+	var arg = 0.0;
+	var err = 0.0;
+	var ref = null;
+	var ns2 = 0;
+	var dt = 0.0;
+	var xh = null;
+	var ni = 0;
+	var w = null;
+	var x = null;
+	var y = null;
+	var n = 0;
+	var j = 0;
+	var i = 0;
+	var k = 0;
+
+	for ( ni = 0; ni < nd.length; ni++ ) {
+		n = nd[ ni ];
+		w = new Float64Array( (2*n) + 15 );
+		x = new Float64Array( n );
+		xh = new Float64Array( n );
+		y = new Float64Array( n );
+
+		for ( j = 0; j < n; j++ ) {
+			x[ j ] = sin( (j+1) * sqrt2 );
+			xh[ j ] = x[ j ];
+			y[ j ] = x[ j ];
+		}
+
+		rffti( n, w );
+
+		dt = (2.0 * PI) / n;
+		ns2 = floor( (n + 1) / 2 );
+		ref = new Float64Array( n );
+
+		sum1 = 0.0;
+		sum2 = 0.0;
+		for ( i = 0; i < n-1; i += 2 ) {
+			sum1 += x[ i ];
+			sum2 += x[ i+1 ];
+		}
+		if ( n % 2 === 1 ) {
+			sum1 += x[ n-1 ];
+		}
+		ref[ 0 ] = sum1 + sum2;
+		if ( n % 2 === 0 ) {
+			ref[ n-1 ] = sum1 - sum2;
+		}
+
+		if ( ns2 >= 2 ) {
+			for ( k = 2; k <= ns2; k++ ) {
+				sumReal = 0.0;
+				sumImag = 0.0;
+				arg = (k-1) * dt;
+
+				for ( i = 0; i < n; i++ ) {
+					angle = arg * i;
+					sumReal += x[ i ] * cos( angle );
+					sumImag += x[ i ] * sin( angle );
+				}
+				ref[ (2*k)-3 ] = sumReal;
+				ref[ (2*k)-2 ] = -sumImag;
+			}
+		}
+
+		rfftf( n, x, 1, w, 1 );
+
+		maxError = 0.0;
+		for ( i = 0; i < n; i++ ) {
+			err = abs( x[i] - ref[i] );
+			maxError = max( maxError, err );
+		}
+		t.ok( maxError < tolerance, 'n=' + n + ' max error: ' + maxError );
+
+		rfftb( n, x, 1, w, 1 );
+		scale = 1.0 / n;
+		maxError = 0.0;
+		for ( i = 0; i < n; i++ ) {
+			err = abs( (x[i] * scale) - xh[i] );
+			maxError = max( maxError, err );
+		}
+		t.ok( maxError < tolerance, 'n=' + n + ' inverse consistency error: ' + maxError );
+	}
+	t.end();
+});
+
+tape( 'rfftf handles edge cases', function test( t ) {
+	var w = null;
+	var x = null;
+	var n = 0;
+
+	n = 1;
+	w = new Float64Array( 2*n + 15 );
+	x = new Float64Array( [ 5.0 ] );
+	rffti( n, w );
+	rfftf( n, x, 1, w, 1 );
+	t.equal( x[0], 5.0, 'n=1 transform equals input' );
+
+	n = 2;
+	w = new Float64Array( 2*n + 15 );
+	x = new Float64Array( [ 1.0, 2.0 ] );
+	rffti( n, w );
+	rfftf( n, x, 1, w, 1 );
+	t.ok( abs( x[0] - 3.0 ) < 1e-10, 'n=2 DC component' );
+	t.ok( abs( x[1] - (-1.0) ) < 1e-10, 'n=2 Nyquist component' );
+
+	t.end();
+});