aboutsummaryrefslogtreecommitdiffstats
path: root/src/main/java/com/softsynth/math
diff options
context:
space:
mode:
Diffstat (limited to 'src/main/java/com/softsynth/math')
-rw-r--r--src/main/java/com/softsynth/math/AudioMath.java82
-rw-r--r--src/main/java/com/softsynth/math/ChebyshevPolynomial.java45
-rw-r--r--src/main/java/com/softsynth/math/FourierMath.java254
-rw-r--r--src/main/java/com/softsynth/math/JustRatio.java47
-rw-r--r--src/main/java/com/softsynth/math/Polynomial.java259
-rw-r--r--src/main/java/com/softsynth/math/PolynomialTableData.java64
-rw-r--r--src/main/java/com/softsynth/math/PrimeFactors.java244
7 files changed, 995 insertions, 0 deletions
diff --git a/src/main/java/com/softsynth/math/AudioMath.java b/src/main/java/com/softsynth/math/AudioMath.java
new file mode 100644
index 0000000..06eb45b
--- /dev/null
+++ b/src/main/java/com/softsynth/math/AudioMath.java
@@ -0,0 +1,82 @@
+/*
+ * Copyright 1998 Phil Burk, Mobileer Inc
+ *
+ * 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.
+ */
+
+package com.softsynth.math;
+
+/**
+ * Miscellaneous math functions useful in Audio
+ *
+ * @author (C) 1998 Phil Burk
+ */
+public class AudioMath {
+ // use scalar to convert natural log to log_base_10
+ private final static double a2dScalar = 20.0 / Math.log(10.0);
+ public static final int CONCERT_A_PITCH = 69;
+ public static final double CONCERT_A_FREQUENCY = 440.0;
+ private static double mConcertAFrequency = CONCERT_A_FREQUENCY;
+
+ /**
+ * Convert amplitude to decibels. 1.0 is zero dB. 0.5 is -6.02 dB.
+ */
+ public static double amplitudeToDecibels(double amplitude) {
+ return Math.log(amplitude) * a2dScalar;
+ }
+
+ /**
+ * Convert decibels to amplitude. Zero dB is 1.0 and -6.02 dB is 0.5.
+ */
+ public static double decibelsToAmplitude(double decibels) {
+ return Math.pow(10.0, decibels / 20.0);
+ }
+
+ /**
+ * Calculate MIDI pitch based on frequency in Hertz. Middle C is 60.0.
+ */
+ public static double frequencyToPitch(double frequency) {
+ return CONCERT_A_PITCH + 12 * Math.log(frequency / mConcertAFrequency) / Math.log(2.0);
+ }
+
+ /**
+ * Calculate frequency in Hertz based on MIDI pitch. Middle C is 60.0. You can use fractional
+ * pitches so 60.5 would give you a pitch half way between C and C#.
+ */
+ public static double pitchToFrequency(double pitch) {
+ return mConcertAFrequency * Math.pow(2.0, ((pitch - CONCERT_A_PITCH) * (1.0 / 12.0)));
+ }
+
+ /**
+ * This can be used to globally adjust the tuning in JSyn from Concert A at 440.0 Hz to
+ * a slightly different frequency. Some orchestras use a higher frequency, eg. 441.0.
+ * This value will be used by pitchToFrequency() and frequencyToPitch().
+ *
+ * @param concertAFrequency
+ */
+ public static void setConcertAFrequency(double concertAFrequency) {
+ mConcertAFrequency = concertAFrequency;
+ }
+
+ public static double getConcertAFrequency() {
+ return mConcertAFrequency;
+ }
+
+ /** Convert a delta value in semitones to a frequency multiplier.
+ * @param semitones
+ * @return scaler For example 2.0 for an input of 12.0 semitones.
+ */
+ public static double semitonesToFrequencyScaler(double semitones) {
+ return Math.pow(2.0, semitones / 12.0);
+ }
+}
diff --git a/src/main/java/com/softsynth/math/ChebyshevPolynomial.java b/src/main/java/com/softsynth/math/ChebyshevPolynomial.java
new file mode 100644
index 0000000..bc0e854
--- /dev/null
+++ b/src/main/java/com/softsynth/math/ChebyshevPolynomial.java
@@ -0,0 +1,45 @@
+/*
+ * Copyright 1997 Phil Burk, Mobileer Inc
+ *
+ * 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.
+ */
+
+package com.softsynth.math;
+
+/**
+ * ChebyshevPolynomial<br>
+ * Used to generate data for waveshaping table oscillators.
+ *
+ * @author Nick Didkovsky (C) 1997 Phil Burk and Nick Didkovsky
+ */
+
+public class ChebyshevPolynomial {
+ static final Polynomial twoX = new Polynomial(2, 0);
+ static final Polynomial one = new Polynomial(1);
+ static final Polynomial oneX = new Polynomial(1, 0);
+
+ /**
+ * Calculates Chebyshev polynomial of specified integer order. Recursively generated using
+ * relation Tk+1(x) = 2xTk(x) - Tk-1(x)
+ *
+ * @return Chebyshev polynomial of specified order
+ */
+ public static Polynomial T(int order) {
+ if (order == 0)
+ return one;
+ else if (order == 1)
+ return oneX;
+ else
+ return Polynomial.minus(Polynomial.mult(T(order - 1), (twoX)), T(order - 2));
+ }
+}
diff --git a/src/main/java/com/softsynth/math/FourierMath.java b/src/main/java/com/softsynth/math/FourierMath.java
new file mode 100644
index 0000000..d133d7f
--- /dev/null
+++ b/src/main/java/com/softsynth/math/FourierMath.java
@@ -0,0 +1,254 @@
+/*
+ * Copyright 2009 Phil Burk, Mobileer Inc
+ *
+ * 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.
+ */
+
+package com.softsynth.math;
+
+//Simple Fast Fourier Transform.
+public class FourierMath {
+ static private final int MAX_SIZE_LOG_2 = 16;
+ static BitReverseTable[] reverseTables = new BitReverseTable[MAX_SIZE_LOG_2];
+ static DoubleSineTable[] sineTables = new DoubleSineTable[MAX_SIZE_LOG_2];
+ static FloatSineTable[] floatSineTables = new FloatSineTable[MAX_SIZE_LOG_2];
+
+ private static class DoubleSineTable {
+ double[] sineValues;
+
+ DoubleSineTable(int numBits) {
+ int len = 1 << numBits;
+ sineValues = new double[1 << numBits];
+ for (int i = 0; i < len; i++) {
+ sineValues[i] = Math.sin((i * Math.PI * 2.0) / len);
+ }
+ }
+ }
+
+ private static double[] getDoubleSineTable(int n) {
+ DoubleSineTable sineTable = sineTables[n];
+ if (sineTable == null) {
+ sineTable = new DoubleSineTable(n);
+ sineTables[n] = sineTable;
+ }
+ return sineTable.sineValues;
+ }
+
+ private static class FloatSineTable {
+ float[] sineValues;
+
+ FloatSineTable(int numBits) {
+ int len = 1 << numBits;
+ sineValues = new float[1 << numBits];
+ for (int i = 0; i < len; i++) {
+ sineValues[i] = (float) Math.sin((i * Math.PI * 2.0) / len);
+ }
+ }
+ }
+
+ private static float[] getFloatSineTable(int n) {
+ FloatSineTable sineTable = floatSineTables[n];
+ if (sineTable == null) {
+ sineTable = new FloatSineTable(n);
+ floatSineTables[n] = sineTable;
+ }
+ return sineTable.sineValues;
+ }
+
+ private static class BitReverseTable {
+ int[] reversedBits;
+
+ BitReverseTable(int numBits) {
+ reversedBits = new int[1 << numBits];
+ for (int i = 0; i < reversedBits.length; i++) {
+ reversedBits[i] = reverseBits(i, numBits);
+ }
+ }
+
+ static int reverseBits(int index, int numBits) {
+ int i, rev;
+
+ for (i = rev = 0; i < numBits; i++) {
+ rev = (rev << 1) | (index & 1);
+ index >>= 1;
+ }
+
+ return rev;
+ }
+ }
+
+ private static int[] getReverseTable(int n) {
+ BitReverseTable reverseTable = reverseTables[n];
+ if (reverseTable == null) {
+ reverseTable = new BitReverseTable(n);
+ reverseTables[n] = reverseTable;
+ }
+ return reverseTable.reversedBits;
+ }
+
+ /**
+ * Calculate the amplitude of the sine wave associated with each bin of a complex FFT result.
+ *
+ * @param ar
+ * @param ai
+ * @param magnitudes
+ */
+ public static void calculateMagnitudes(double ar[], double ai[], double[] magnitudes) {
+ for (int i = 0; i < magnitudes.length; ++i) {
+ magnitudes[i] = Math.sqrt((ar[i] * ar[i]) + (ai[i] * ai[i]));
+ }
+ }
+
+ /**
+ * Calculate the amplitude of the sine wave associated with each bin of a complex FFT result.
+ *
+ * @param ar
+ * @param ai
+ * @param magnitudes
+ */
+ public static void calculateMagnitudes(float ar[], float ai[], float[] magnitudes) {
+ for (int i = 0; i < magnitudes.length; ++i) {
+ magnitudes[i] = (float) Math.sqrt((ar[i] * ar[i]) + (ai[i] * ai[i]));
+ }
+ }
+
+ public static void transform(int sign, int n, double ar[], double ai[]) {
+ double scale = (sign > 0) ? (2.0 / n) : (0.5);
+
+ int numBits = FourierMath.numBits(n);
+ int[] reverseTable = getReverseTable(numBits);
+ double[] sineTable = getDoubleSineTable(numBits);
+ int mask = n - 1;
+ int cosineOffset = n / 4; // phase offset between cos and sin
+
+ int i, j;
+ for (i = 0; i < n; i++) {
+ j = reverseTable[i];
+ if (j >= i) {
+ double tempr = ar[j] * scale;
+ double tempi = ai[j] * scale;
+ ar[j] = ar[i] * scale;
+ ai[j] = ai[i] * scale;
+ ar[i] = tempr;
+ ai[i] = tempi;
+ }
+ }
+
+ int mmax, stride;
+ int numerator = sign * n;
+ for (mmax = 1, stride = 2 * mmax; mmax < n; mmax = stride, stride = 2 * mmax) {
+ int phase = 0;
+ int phaseIncrement = numerator / (2 * mmax);
+ for (int m = 0; m < mmax; ++m) {
+ double wr = sineTable[(phase + cosineOffset) & mask]; // cosine
+ double wi = sineTable[phase];
+
+ for (i = m; i < n; i += stride) {
+ j = i + mmax;
+ double tr = (wr * ar[j]) - (wi * ai[j]);
+ double ti = (wr * ai[j]) + (wi * ar[j]);
+ ar[j] = ar[i] - tr;
+ ai[j] = ai[i] - ti;
+ ar[i] += tr;
+ ai[i] += ti;
+ }
+
+ phase = (phase + phaseIncrement) & mask;
+ }
+ mmax = stride;
+ }
+ }
+
+ public static void transform(int sign, int n, float ar[], float ai[]) {
+ float scale = (sign > 0) ? (2.0f / n) : (0.5f);
+
+ int numBits = FourierMath.numBits(n);
+ int[] reverseTable = getReverseTable(numBits);
+ float[] sineTable = getFloatSineTable(numBits);
+ int mask = n - 1;
+ int cosineOffset = n / 4; // phase offset between cos and sin
+
+ int i, j;
+ for (i = 0; i < n; i++) {
+ j = reverseTable[i];
+ if (j >= i) {
+ float tempr = ar[j] * scale;
+ float tempi = ai[j] * scale;
+ ar[j] = ar[i] * scale;
+ ai[j] = ai[i] * scale;
+ ar[i] = tempr;
+ ai[i] = tempi;
+ }
+ }
+
+ int mmax, stride;
+ int numerator = sign * n;
+ for (mmax = 1, stride = 2 * mmax; mmax < n; mmax = stride, stride = 2 * mmax) {
+ int phase = 0;
+ int phaseIncrement = numerator / (2 * mmax);
+ for (int m = 0; m < mmax; ++m) {
+ float wr = sineTable[(phase + cosineOffset) & mask]; // cosine
+ float wi = sineTable[phase];
+
+ for (i = m; i < n; i += stride) {
+ j = i + mmax;
+ float tr = (wr * ar[j]) - (wi * ai[j]);
+ float ti = (wr * ai[j]) + (wi * ar[j]);
+ ar[j] = ar[i] - tr;
+ ai[j] = ai[i] - ti;
+ ar[i] += tr;
+ ai[i] += ti;
+ }
+
+ phase = (phase + phaseIncrement) & mask;
+ }
+ mmax = stride;
+ }
+ }
+
+ /**
+ * Calculate log2(n)
+ *
+ * @param powerOf2 must be a power of two, for example 512 or 1024
+ * @return for example, 9 for an input value of 512
+ */
+ public static int numBits(int powerOf2) {
+ int i;
+ assert ((powerOf2 & (powerOf2 - 1)) == 0); // is it a power of 2?
+ for (i = -1; powerOf2 > 0; powerOf2 = powerOf2 >> 1, i++)
+ ;
+ return i;
+ }
+
+ /**
+ * Calculate an FFT in place, modifying the input arrays.
+ *
+ * @param n
+ * @param ar
+ * @param ai
+ */
+ public static void fft(int n, double ar[], double ai[]) {
+ transform(1, n, ar, ai); // TODO -1 or 1
+ }
+
+ /**
+ * Calculate an inverse FFT in place, modifying the input arrays.
+ *
+ * @param n
+ * @param ar
+ * @param ai
+ */
+ public static void ifft(int n, double ar[], double ai[]) {
+ transform(-1, n, ar, ai); // TODO -1 or 1
+ }
+}
diff --git a/src/main/java/com/softsynth/math/JustRatio.java b/src/main/java/com/softsynth/math/JustRatio.java
new file mode 100644
index 0000000..f4070b4
--- /dev/null
+++ b/src/main/java/com/softsynth/math/JustRatio.java
@@ -0,0 +1,47 @@
+/*
+ * Copyright 2009 Phil Burk, Mobileer Inc
+ *
+ * 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.
+ */
+
+package com.softsynth.math;
+
+public class JustRatio {
+ public long numerator;
+ public long denominator;
+
+ public JustRatio(long numerator, long denominator) {
+ this.numerator = numerator;
+ this.denominator = denominator;
+ }
+
+ public JustRatio(int numerator, int denominator) {
+ this.numerator = numerator;
+ this.denominator = denominator;
+ }
+
+ public double getValue() {
+ return (double) numerator / denominator;
+ }
+
+ public void invert() {
+ long temp = denominator;
+ denominator = numerator;
+ numerator = temp;
+ }
+
+ @Override
+ public String toString() {
+ return numerator + "/" + denominator;
+ }
+}
diff --git a/src/main/java/com/softsynth/math/Polynomial.java b/src/main/java/com/softsynth/math/Polynomial.java
new file mode 100644
index 0000000..6c6f96a
--- /dev/null
+++ b/src/main/java/com/softsynth/math/Polynomial.java
@@ -0,0 +1,259 @@
+/*
+ * Copyright 1997 Phil Burk, Mobileer Inc
+ *
+ * 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.
+ */
+
+package com.softsynth.math;
+
+import org.slf4j.Logger;
+import org.slf4j.LoggerFactory;
+
+import java.util.Vector;
+
+/**
+ * Polynomial<br>
+ * Implement polynomial using Vector as coefficient holder. Element index is power of X, value at a
+ * given index is coefficient.<br>
+ * <br>
+ *
+ * @author Nick Didkovsky, (C) 1997 Phil Burk and Nick Didkovsky
+ */
+
+public class Polynomial {
+
+ private static final Logger LOGGER = LoggerFactory.getLogger(Polynomial.class);
+
+ private final Vector terms;
+
+ // TODO: Does this need to exist?
+ static class DoubleHolder {
+ double value;
+
+ public DoubleHolder(double val) {
+ value = val;
+ }
+
+ public double get() {
+ return value;
+ }
+
+ public void set(double val) {
+ value = val;
+ }
+ }
+
+ /** create a polynomial with no terms */
+ public Polynomial() {
+ terms = new Vector();
+ }
+
+ /** create a polynomial with one term of specified constant */
+ public Polynomial(double c0) {
+ this();
+ appendTerm(c0);
+ }
+
+ /** create a polynomial with two terms with specified coefficients */
+ public Polynomial(double c1, double c0) {
+ this(c0);
+ appendTerm(c1);
+ }
+
+ /** create a polynomial with specified coefficients */
+ public Polynomial(double c2, double c1, double c0) {
+ this(c1, c0);
+ appendTerm(c2);
+ }
+
+ /** create a polynomial with specified coefficients */
+ public Polynomial(double c3, double c2, double c1, double c0) {
+ this(c2, c1, c0);
+ appendTerm(c3);
+ }
+
+ /** create a polynomial with specified coefficients */
+ public Polynomial(double c4, double c3, double c2, double c1, double c0) {
+ this(c3, c2, c1, c0);
+ appendTerm(c4);
+ }
+
+ /**
+ * Append a term with specified coefficient. Power will be next available order (ie if the
+ * polynomial is of order 2, appendTerm will supply the coefficient for x^3
+ */
+ public void appendTerm(double coefficient) {
+ terms.addElement(new DoubleHolder(coefficient));
+ }
+
+ /** Set the coefficient of given term */
+ public void setTerm(double coefficient, int power) {
+ // If setting a term greater than the current order of the polynomial, pad with zero terms
+ int size = terms.size();
+ if (power >= size) {
+ for (int i = 0; i < (power - size + 1); i++) {
+ appendTerm(0);
+ }
+ }
+ ((DoubleHolder) terms.elementAt(power)).set(coefficient);
+ }
+
+ /**
+ * Add the coefficient of given term to the specified coefficient. ex. addTerm(3, 1) add 3x to a
+ * polynomial, addTerm(4, 3) adds 4x^3
+ */
+ public void addTerm(double coefficient, int power) {
+ setTerm(coefficient + get(power), power);
+ }
+
+ /** @return coefficient of nth term (first term=0) */
+ public double get(int power) {
+ if (power >= terms.size())
+ return 0.0;
+ else
+ return ((DoubleHolder) terms.elementAt(power)).get();
+ }
+
+ /** @return number of terms in this polynomial */
+ public int size() {
+ return terms.size();
+ }
+
+ /**
+ * Add two polynomials together
+ *
+ * @return new Polynomial that is the sum of p1 and p2
+ */
+ public static Polynomial plus(Polynomial p1, Polynomial p2) {
+ Polynomial sum = new Polynomial();
+ for (int i = 0; i < Math.max(p1.size(), p2.size()); i++) {
+ sum.appendTerm(p1.get(i) + p2.get(i));
+ }
+ return sum;
+ }
+
+ /**
+ * Subtract polynomial from another. (First arg - Second arg)
+ *
+ * @return new Polynomial p1 - p2
+ */
+ public static Polynomial minus(Polynomial p1, Polynomial p2) {
+ Polynomial sum = new Polynomial();
+ for (int i = 0; i < Math.max(p1.size(), p2.size()); i++) {
+ sum.appendTerm(p1.get(i) - p2.get(i));
+ }
+ return sum;
+ }
+
+ /**
+ * Multiply two Polynomials
+ *
+ * @return new Polynomial that is the product p1 * p2
+ */
+
+ public static Polynomial mult(Polynomial p1, Polynomial p2) {
+ Polynomial product = new Polynomial();
+ for (int i = 0; i < p1.size(); i++) {
+ for (int j = 0; j < p2.size(); j++) {
+ product.addTerm(p1.get(i) * p2.get(j), i + j);
+ }
+ }
+ return product;
+ }
+
+ /**
+ * Multiply a Polynomial by a scaler
+ *
+ * @return new Polynomial that is the product p1 * p2
+ */
+
+ public static Polynomial mult(double scaler, Polynomial p1) {
+ Polynomial product = new Polynomial();
+ for (int i = 0; i < p1.size(); i++) {
+ product.appendTerm(p1.get(i) * scaler);
+ }
+ return product;
+ }
+
+ /** Evaluate this polynomial for x */
+ public double evaluate(double x) {
+ double result = 0.0;
+ for (int i = 0; i < terms.size(); i++) {
+ result += get(i) * Math.pow(x, i);
+ }
+ return result;
+ }
+
+ @Override
+ public String toString() {
+ String s = "";
+ if (size() == 0)
+ s = "empty polynomial";
+ boolean somethingPrinted = false;
+ for (int i = size() - 1; i >= 0; i--) {
+ if (get(i) != 0.0) {
+ if (somethingPrinted)
+ s += " + ";
+ String coeff = "";
+ // if (get(i) == (int)(get(i)))
+ // coeff = (int)(get(i)) + "";
+ if ((get(i) != 1.0) || (i == 0))
+ coeff += get(i);
+ if (i == 0)
+ s += coeff;
+ else {
+ String power = "";
+ if (i != 1)
+ power = "^" + i;
+ s += coeff + "x" + power;
+ }
+ somethingPrinted = true;
+ }
+ }
+ return s;
+ }
+
+ public static void main(String[] args) {
+ Polynomial p1 = new Polynomial();
+ LOGGER.debug("p1=" + p1);
+ Polynomial p2 = new Polynomial(3);
+ LOGGER.debug("p2=" + p2);
+ Polynomial p3 = new Polynomial(2, 3);
+ LOGGER.debug("p3=" + p3);
+ Polynomial p4 = new Polynomial(1, 2, 3);
+ LOGGER.debug("p4=" + p4);
+ LOGGER.debug("p4*5=" + Polynomial.mult(5.0, p4));
+
+ LOGGER.debug("{}", p4.evaluate(10));
+
+ LOGGER.debug("{}", Polynomial.plus(p4, p1));
+ LOGGER.debug("{}", Polynomial.minus(p4, p3));
+ p4.setTerm(12.2, 5);
+ LOGGER.debug("{}", p4);
+ p4.addTerm(0.8, 5);
+ LOGGER.debug("{}", p4);
+ p4.addTerm(0.8, 7);
+ LOGGER.debug("{}", p4);
+ LOGGER.debug("{}", Polynomial.mult(p3, p2));
+ LOGGER.debug("{}", Polynomial.mult(p3, p3));
+ LOGGER.debug("{}", Polynomial.mult(p2, p2));
+
+ Polynomial t2 = new Polynomial(2, 0, -1); // 2x^2-1, Chebyshev Polynomial of order 2
+ Polynomial t3 = new Polynomial(4, 0, -3, 0); // 4x^3-3x, Chebyshev Polynomial of order 3
+ // Calculate Chebyshev Polynomial of order 4 from relation Tk+1(x) = 2xTk(x) - Tk-1(x)
+ Polynomial t4 = Polynomial.minus(Polynomial.mult(t3, (new Polynomial(2, 0))), t2);
+ LOGGER.debug(t2 + "\n" + t3 + "\n" + t4);
+ // com.softsynth.jmsl.util
+
+ }
+}
diff --git a/src/main/java/com/softsynth/math/PolynomialTableData.java b/src/main/java/com/softsynth/math/PolynomialTableData.java
new file mode 100644
index 0000000..8110151
--- /dev/null
+++ b/src/main/java/com/softsynth/math/PolynomialTableData.java
@@ -0,0 +1,64 @@
+/*
+ * Copyright 1997 Phil Burk, Mobileer Inc
+ *
+ * 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.
+ */
+
+package com.softsynth.math;
+
+/**
+ * PolynomialTableData<br>
+ * Provides an array of double[] containing data generated by a polynomial.<br>
+ * This is typically used with ChebyshevPolynomial. Input to Polynomial is -1..1, output is -1..1.
+ *
+ * @author Nick Didkovsky (C) 1997 Phil Burk and Nick Didkovsky
+ * @see ChebyshevPolynomial
+ * @see Polynomial
+ */
+
+public class PolynomialTableData {
+
+ double[] data;
+ Polynomial polynomial;
+
+ /**
+ * Constructor which fills double[numFrames] with Polynomial data -1..1<br>
+ * Note that any Polynomial can plug in here, just make sure output is -1..1 when input ranges
+ * from -1..1
+ */
+ public PolynomialTableData(Polynomial polynomial, int numFrames) {
+ data = new double[numFrames];
+ this.polynomial = polynomial;
+ buildData();
+ }
+
+ public double[] getData() {
+ return data;
+ }
+
+ void buildData() {
+ double xInterval = 2.0 / (data.length - 1); // FIXED, added "- 1"
+ double x;
+ for (int i = 0; i < data.length; i++) {
+ x = i * xInterval - 1.0;
+ data[i] = polynomial.evaluate(x);
+ // LOGGER.debug("x = " + x + ", p(x) = " + data[i] );
+ }
+
+ }
+
+ public static void main(String[] args) {
+ PolynomialTableData chebData = new PolynomialTableData(ChebyshevPolynomial.T(2), 8);
+ }
+
+}
diff --git a/src/main/java/com/softsynth/math/PrimeFactors.java b/src/main/java/com/softsynth/math/PrimeFactors.java
new file mode 100644
index 0000000..06c0d55
--- /dev/null
+++ b/src/main/java/com/softsynth/math/PrimeFactors.java
@@ -0,0 +1,244 @@
+/*
+ * Copyright 2011 Phil Burk, Mobileer Inc
+ *
+ * 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.
+ */
+
+package com.softsynth.math;
+
+import java.util.ArrayList;
+
+/**
+ * Tool for factoring primes and prime ratios. This class contains a static array of primes
+ * generated using the Sieve of Eratosthenes.
+ *
+ * @author Phil Burk (C) 2011 Mobileer Inc
+ */
+public class PrimeFactors {
+ private static final int SIEVE_SIZE = 1000;
+ private static int[] primes;
+ private final int[] factors;
+
+ static {
+ // Use Sieve of Eratosthenes to fill Prime table
+ boolean[] sieve = new boolean[SIEVE_SIZE];
+ ArrayList<Integer> primeList = new ArrayList<Integer>();
+ int i = 2;
+ while (i < (SIEVE_SIZE / 2)) {
+ if (!sieve[i]) {
+ primeList.add(i);
+ int multiple = 2 * i;
+ while (multiple < SIEVE_SIZE) {
+ sieve[multiple] = true;
+ multiple += i;
+ }
+ }
+ i += 1;
+ }
+ primes = primeListToArray(primeList);
+ }
+
+ private static int[] primeListToArray(ArrayList<Integer> primeList) {
+ int[] primes = new int[primeList.size()];
+ for (int i = 0; i < primes.length; i++) {
+ primes[i] = primeList.get(i);
+ }
+ return primes;
+ }
+
+ public PrimeFactors(int[] factors) {
+ this.factors = factors;
+ }
+
+ public PrimeFactors(int numerator, int denominator) {
+ int[] topFactors = factor(numerator);
+ int[] bottomFactors = factor(denominator);
+ factors = subtract(topFactors, bottomFactors);
+ }
+
+ public PrimeFactors subtract(PrimeFactors pf) {
+ return new PrimeFactors(subtract(factors, pf.factors));
+ }
+
+ public PrimeFactors add(PrimeFactors pf) {
+ return new PrimeFactors(add(factors, pf.factors));
+ }
+
+ public static int[] subtract(int[] factorsA, int[] factorsB) {
+ int max;
+ int min;
+ if (factorsA.length > factorsB.length) {
+ max = factorsA.length;
+ min = factorsB.length;
+ } else {
+
+ min = factorsA.length;
+ max = factorsB.length;
+ }
+ ArrayList<Integer> primeList = new ArrayList<Integer>();
+ int i;
+ for (i = 0; i < min; i++) {
+ primeList.add(factorsA[i] - factorsB[i]);
+ }
+ if (factorsA.length > factorsB.length) {
+ for (; i < max; i++) {
+ primeList.add(factorsA[i]);
+ }
+ } else {
+ for (; i < max; i++) {
+ primeList.add(0 - factorsB[i]);
+ }
+ }
+ trimPrimeList(primeList);
+ return primeListToArray(primeList);
+ }
+
+ public static int[] add(int[] factorsA, int[] factorsB) {
+ int max;
+ int min;
+ if (factorsA.length > factorsB.length) {
+ max = factorsA.length;
+ min = factorsB.length;
+ } else {
+ min = factorsA.length;
+ max = factorsB.length;
+ }
+ ArrayList<Integer> primeList = new ArrayList<Integer>();
+ int i;
+ for (i = 0; i < min; i++) {
+ primeList.add(factorsA[i] + factorsB[i]);
+ }
+ if (factorsA.length > factorsB.length) {
+ for (; i < max; i++) {
+ primeList.add(factorsA[i]);
+ }
+ } else if (factorsB.length > factorsA.length) {
+ for (; i < max; i++) {
+ primeList.add(factorsB[i]);
+ }
+ }
+ trimPrimeList(primeList);
+ return primeListToArray(primeList);
+ }
+
+ private static void trimPrimeList(ArrayList<Integer> primeList) {
+ int i;
+ // trim zero factors off end.
+ for (i = primeList.size() - 1; i >= 0; i--) {
+ if (primeList.get(i) == 0) {
+ primeList.remove(i);
+ } else {
+ break;
+ }
+ }
+ }
+
+ public static int[] factor(int n) {
+ ArrayList<Integer> primeList = new ArrayList<Integer>();
+ int i = 0;
+ int p = primes[i];
+ int exponent = 0;
+ while (n > 1) {
+ // does the prime number divide evenly into n?
+ int d = n / p;
+ int m = d * p;
+ if (m == n) {
+ n = d;
+ exponent += 1;
+ } else {
+ primeList.add(exponent);
+ exponent = 0;
+ i += 1;
+ p = primes[i];
+ }
+ }
+ if (exponent > 0) {
+ primeList.add(exponent);
+ }
+ return primeListToArray(primeList);
+ }
+
+ /**
+ * Get prime from table.
+ *
+ *
+ * @param n Warning: Do not exceed getPrimeCount()-1.
+ * @return Nth prime number, the 0th prime is 2
+ */
+ public static int getPrime(int n) {
+ return primes[n];
+ }
+
+ /**
+ * @return the number of primes stored in the table
+ */
+ public static int getPrimeCount() {
+ return primes.length;
+ }
+
+ public JustRatio getJustRatio() {
+ long n = 1;
+ long d = 1;
+ for (int i = 0; i < factors.length; i++) {
+ int exponent = factors[i];
+ int p = primes[i];
+ if (exponent > 0) {
+ for (int k = 0; k < exponent; k++) {
+ n = n * p;
+ }
+ } else if (exponent < 0) {
+ exponent = 0 - exponent;
+ for (int k = 0; k < exponent; k++) {
+ d = d * p;
+ }
+ }
+ }
+ return new JustRatio(n, d);
+ }
+
+ public int[] getFactors() {
+ return factors.clone();
+ }
+
+ @Override
+ public String toString() {
+ StringBuffer buffer = new StringBuffer();
+ printFactors(buffer, 1);
+ buffer.append("/");
+ printFactors(buffer, -1);
+ return buffer.toString();
+ }
+
+ private void printFactors(StringBuffer buffer, int sign) {
+ boolean gotSome = false;
+ for (int i = 0; i < factors.length; i++) {
+ int pf = factors[i] * sign;
+ if (pf > 0) {
+ if (gotSome)
+ buffer.append('*');
+ int prime = primes[i];
+ if (pf == 1) {
+ buffer.append("" + prime);
+ } else if (pf == 2) {
+ buffer.append(prime + "*" + prime);
+ } else if (pf > 2) {
+ buffer.append("(" + prime + "^" + pf + ")");
+ }
+ gotSome = true;
+ }
+ }
+ if (!gotSome) {
+ buffer.append("1");
+ }
+ }
+}
generated by cgit v1.2.3 (git 2.39.1) at 2024-11-28 13:58:22 +0000