intbutterworthhighpassfilter.java

一个用java写的地震分析软件(无源码)-used to write a seismic analysis software (without source)

package org.trinet.util;

public class IntButterworthHighPassFilter implements IntFilterIF {

    protected static boolean reverse = false;
    public static final int NZEROS = 4;
    public static final int NPOLES = 4;

    private double [] xv = new double [NZEROS+1];
    private double [] yv = new double [NPOLES+1];

    private double gain;
    private double coef1;
    private double coef2;
    private double coef3;
    private double coef4;

    protected static final IntButterworthHighPassFilter H_SEEDCHAN_FILTER = new IntButterworthHighPassFilter("H");
    protected static final IntButterworthHighPassFilter E_SEEDCHAN_FILTER = new IntButterworthHighPassFilter("E");

    public IntButterworthHighPassFilter () { }

    public IntButterworthHighPassFilter (boolean enableReverseFiltering) {
         this();
         reverse = enableReverseFiltering;
    }

    private IntButterworthHighPassFilter (String seedchan) {
        String substr = seedchan.substring(0,1).toUpperCase();
        if (substr.equals("H")) {
// 4 pole butterworth hp 80 sps corner at 1 Hz
            gain  =  1.108101438;
            coef1 = -0.8144059977;
            coef2 =  3.4247473473;
            coef3 = -5.4051668617;
            coef4 =  3.7947911031;
// 4 pole butterworth hp 80 sps corner at .5 Hz
/* 
            gain  = 1.052651765;
            coef1 = -0.9024653874;
            coef2 = 3.7024671485;
            coef3 = -5.6973900039 ;
            coef4 = 3.8973859825;
*/
        }
        else if (substr.equals("E")) {
// 4 pole butterworth hp 100 sps corner at 1 Hz
            gain  =  1.085574782;
            coef1 = -0.8485559993;
            coef2 =  3.5335352195;
            coef3 = -5.5208191366;
            coef4 =  3.8358255406;
        }
    }

    public FilterIF getFilter(String seedchan) {
        String substr = seedchan.substring(0,1).toUpperCase();
        if (substr.equals("H")) return H_SEEDCHAN_FILTER; 
        else if (substr.equals("E"))  return E_SEEDCHAN_FILTER; 
        else return null;
    }
    /** This is bogus but needed for compatiblity.
     * @deprecated */

    public FilterIF getFilter(int sps) {
        if (sps == 80) return H_SEEDCHAN_FILTER; 
        else if (sps == 100)  return E_SEEDCHAN_FILTER; 
        else return null;
    }

    public Object filter(Object in) {
       return filter((int []) in);
    }

    public int [] filter(int [] in) {
       int [] values = demean(in);
       values = org.trinet.util.AmplitudeCosineTaper.taper(values, 0.05);
       values = doFilter(values);
       return (reverse) ? doReverseFilter(values) : values;
    }

    public int [] demean(int [] values) {
        
        if (values == null) return values;
        int size = values.length;
        if (size == 0) return values;

        double sum = 0.;
        for (int idx = 0; idx < size; idx++) {
            sum += values[idx];
        }
        int mean = (int) Math.round(sum/size);
        for (int idx = 0; idx < size; idx++) {
            values[idx] = values[idx] - mean;
        }
        return values;
        
    }

    private void zeroBuffers() {
        java.util.Arrays.fill(xv, 0.);
        java.util.Arrays.fill(yv, 0.);
    }

    private int [] doFilter(int [] in) {
        if (in == null || gain == 0.) return in;
        zeroBuffers();
        int size = in.length;
        for (int idx=0; idx < size; idx++) {
            in[idx] = filterSample(in[idx]);
        }
        return in;
    }

    private int [] doReverseFilter(int [] in) {
        if (in == null || gain == 0.) return in;
        zeroBuffers();
        int size = in.length;
        for (int idx=size-1; idx > -1; idx--) {
            in[idx] = filterSample(in[idx]);
        }
        return in;
    }

    private int filterSample(int value) {
        double inputValue = value;
        xv[0] = xv[1];
        xv[1] = xv[2];
        xv[2] = xv[3];
        xv[3] = xv[4]; 
        xv[4] = inputValue / gain;
        yv[0] = yv[1];
        yv[1] = yv[2];
        yv[2] = yv[3];
        yv[3] = yv[4]; 
        yv[4] = (xv[0] + xv[4]) - 4 * (xv[1] + xv[3]) + 6 * xv[2]
                 + ( coef1 * yv[0]) + (  coef2 * yv[1])
                 + ( coef3 * yv[2]) + (  coef4 * yv[3]);
        return (int) Math.rint(yv[4]);
    }

    public static void setReverseFiltering(boolean value) {
        reverse = value;
    }

}