filter.h

/* ======== IGeoS ===== Distribition full ===== IGeoS =========
**                                                                         
** filter.h     : Created by 'AMSS-developer' on Fri Jan 12 19:49:03 2018
**                                                                         
** Product      : IGeoS - Integrated Geoscience Software
**                                                                         
** Description  : System for seismic, well log, and potential-field data analysis
**                                                                         
** =================== Limited License: ===================================
** 
** This software is provided free under the following terms and conditions:
** 
** 1. Permission to use, copy, and modify this software 
** for non-commercial purposes without fee is hereby granted, provided  
** that this copyright notice, the warranty disclaimer, and this
** permission notice appear in all copies.
** 
** 2. Distribution of this software or any part of it "bundled" in with
** any product is considered to be a 'commercial purpose'.
** 
** 3. Any new or adapted code developed to operate as a part of this
** software shall be contributed to the authors and distributed under
** the same license.
** 
** ================== Warranty Disclaimer: ================================
** 
** This software is provided "as is", with no support and without
** obligation on the part of the author to assist in its use, correction,
** modification, or enhancement. No guarantees or warranties,
** either express or implied, and regarding the accuracy, safety, or
** fitness for any particular purpose are provided by any contributor
** to this software package.
** 
** ======== IGeoS ===== Distribition full ===== IGeoS ========= */

#ifndef MOD_FILTER
#define MOD_FILTER

#include "../store/store.h"     /* typedef of STORE */
#include "../wllog/wllog.h"
#include "sia_fft.h"

/*--------- frequency bands -------------*/

#define POLYGON         1       /* all low-pass, band-pass, high-pass,
                                                        and notch filters */
#define PREDEFINED      2       /* filter retrieved from elsewhere */
#define BAND_SLOPE      3       /* band slope filter */
#define GAUSSIAN        4       /* gaussian-shaped */
#define SP_BALANCE      5       /* spectral balancing */
#define BUTTERWORTH     6       // Butterworth
#define F_WAVELET       7       // standard wavelet
#define RESPONSE        8       // response correction

/*--- other lists ----------------*/

#define RAMP            9

/*---- types of polygon filters -------*/

#define GENERAL         200
#define BAND_PASS       201
#define BAND_CUT        202
#define LOW_PASS        203
#define HIGH_PASS       204
#define NOTCH           205

struct Filter_Shape
{
  Filter_Shape() {}

  Filter_Shape(const Filter_Shape &r)
        : frequency(r.frequency),amplitude(r.amplitude) {}

  Filter_Shape &operator=(const Filter_Shape &r) 
        { frequency=r.frequency; amplitude=r.amplitude; return *this; }

  HEADER_PARAM  frequency, amplitude;
};

struct FILTER_SHAPES : public ARRAY<Filter_Shape> 
{
  FILTER_SHAPES() : ARRAY<Filter_Shape>() {}
  
  void read_preset_band ( int n, float *amps ); 

  void show() const;    
};

class FILTER;

struct FILTER_FFT : public SIA_FFT
{
  FILTER_FFT() : SIA_FFT(),_cfilter(NULL),alloc_freq(0) {}
  FILTER_FFT( TRACE *t, double pad = 0.25 );

  ~FILTER_FFT() { free(); }

  void alloc( int num_samples, double pad=0.25 ) 
        { num_points(num_samples,pad); 
          realloc(num_freq()); 
        }       


  void filter ( DATA_SAMPLE *data, int num_data_samples, boolean real_filter ) 
                                        { filter(data,num_data_samples,_cfilter,real_filter); }

  void filter ( DATA_SAMPLE *data, int num_data_samples, COMPLEX *filter, boolean real_filter );

  void time_shift ( double dt );        

  COMPLEX *cfilter() const { return _cfilter; } 

  COMPLEX cfilter( double freq ) const;         

  void clear_cfilter() { memset(_cfilter,0,num_freq()*sizeof(COMPLEX)); }

  void inverse_waterlevel ( double wlevel, double cutoff_freq_times_dt );

private:

  COMPLEX *_cfilter;    
  int alloc_freq;       

  void free() { if ( _cfilter ) delete[] _cfilter; _cfilter = NULL; }
  void realloc( int nfreq );
};

struct Filter : public SIA_List, WLLOG_PROC
{
    Filter ( Filter *pred=NULL, int _type=0, boolean _zero_phase=FALSE );
    ~Filter();

    int type;   

    /*----- flags to generate filter for each trace ---------------*/

    boolean     individual_gate,        
                individual_filter;      

    HEADER_PARAM        time_start, time_end,   
                        f_base,                 
                        f_scalar;               

    double              df;                     

    boolean read_time_gate ();  

    virtual void show() const;  

    virtual void plotmtv_filter( FILE *f ) const;       

    virtual void proc ( FILTER *FIL, TRACE *t );                                
    virtual void proc ( FILTER *FIL, DATA_SAMPLE *d, int num_samples );         

    virtual boolean build ( FILTER *FIL, TRACE *t );    

    virtual boolean build ( FILTER *F, byte *header, 
                        int samples_per_trace, double sample_interval )
                        { return FAIL; }

    /*------- work arrays and parameters --------------------------*/

    FILTER_FFT  fft;            

    boolean     zero_phase;     
    int         start,          
                nt,             
                start_freq, end_freq;   

    COMPLEX     *filter;        
};

#include "filter_polygon.h"
#include "filter_bslope.h"
#include "filter_spbal.h"
#include "filter_bworth.h"
#include "filter_wavelet.h"
#include "filter_response.h"
#include "filter_defined.h"

struct FILTER_TRANGE : public SIA_List
{
  FILTER_TRANGE(FILTER_TRANGE *pred=NULL)
        : SIA_List(pred),filter(NULL)
        { ramp_start.set(REAL,-1.0);ramp_end.set(REAL,-1.0); }

  FILTER_TRANGE( double f_locut, double f_lopass, 
                        double f_hipass, double f_hicut, 
                        boolean bandpass=TRUE,
                        FILTER_TRANGE *pred=NULL );

  FILTER_TRANGE( double f_locut, double f_lopass, boolean high=TRUE,
                                FILTER_TRANGE *pred=NULL );
                
  FILTER_TRANGE( const ARRAY<POINT> &shape, FILTER_TRANGE *pred=NULL );
                
  FILTER_TRANGE( Filter *f );   

  ~FILTER_TRANGE();

  Filter *filter;       

  HEADER_PARAM  ramp_start, ramp_end;

  Filter *add_filter( Filter *f ) 
                { if (filter) return (Filter*)filter->add_last(f); return filter = f; }

  /*---- read parameter lists ---------*/

  void  read_ramp (),
        read_band ( int filter_type ),
        read_preset_band ( Filter_Shape *s, int n, float *amps ),
        read_sp_balance ( int filter_type ),
        read_polygon ( int filter_type, int rep ),
        read_predefined ( int filter_type ),
        read_butterworth ( int filter_type ),
        read_wavelet ( int filter_type ),
        read_bslope ( int filter_type ),
        read_response ( int filter_type );

  void proc ( FILTER *FIL, TRACE *t, int num_filter  ); 


  void proc ( FILTER *FIL, DATA_SAMPLE *d, int num_samples );           

  void show() const;    

  void plotmtv_filter( FILE *f ) const;         

private:

  boolean remove_filter() 
        { 
          if ( filter )
            if ( Filter *f = ((Filter*)(filter->last())) )
            {
              f->_remove();

              if ( filter == f )
                filter = NULL;

              delete f;
              return OK;
            }

          return FAIL;
        }
};

struct FILTER : public SIA_MODULE
{
  FILTER();     

  FILTER( double f_locut, double f_lopass, double f_hipass, double f_hicut, 
                                                boolean bandpass=TRUE );

  FILTER( double f_locut, double f_lopass, boolean highpass=TRUE );
        
  FILTER( const ARRAY<POINT> &shape);

  FILTER( FILTER_TRANGE *f );   

  ~FILTER();

  boolean time_domain;          

  boolean panels;               

  AHEADER panel_count;                  
  HEADER_PARAM panel_count_start;       

  FILTER_TRANGE *trange;        

  int   max_allocated;          


  void add_filter( FILTER_TRANGE *f ) { if ( trange ) trange->add_last(f); else trange=f; }

  void build( GATHER *g );      


  void build( TRACE *t );       


  void allocate ( int num_samples );

  boolean apply( TRACE *t, int count=-1 );

  boolean apply( DATA_SAMPLE *srray, int num_samples);  

  void show( const char *title = NULL ) const;  

  void plotmtv_filter( const char *file_name, const char *heading=NULL, 
                        const char *line_color="1", boolean append=FALSE ) const;       

  CHARSTR module_name( boolean active );        
  int edit();                                   
  boolean process();                            
  
  double sample_interval;       

  /*------- work areas ------------------------------------*/

  DATA_SAMPLE   *rt;                            

  COMPLEX       *ct() const { return (COMPLEX*)rt; }    

  int           *work_n;        
  double        *work;

private:

  boolean remove_trange() 
        { 
          if ( trange )
            if ( FILTER_TRANGE *f = ((FILTER_TRANGE*)(trange->last())) )
            {
              f->_remove();

              if ( trange == f )
                trange = NULL;

              delete f;
              return OK;
            }

          return FAIL;
        }
};

struct FILTER_EDIT : public UI_X_PROPERTY_EDIT
{
  FILTER_EDIT( int _mode, double f_lo, double f_hi, 
                const char *comment, const char *status );      

  SELECTION_EDIT mode;  

  REAL_EDIT     f_locut, f_hicut;       

  FILTER filter;        

  ARRAY<UI_X_PROPERTY_EDIT*> sub_items();

//  boolean set_value( const char *new_value );
//      ///< Set value of the sub-item #'ind' using the 'new_value' string.
//      ///< This method should be called after a value has been edited
//      ///< by the user and provided in 'new_value'.

  void curr_value();    
};

/*=============================================================================
        Applies predefined real or complex filter to a time
        sequence.
        nfft    is the number of samples in 'data'
        work    is complex work area for fft (may coincide with 'data')
==============================================================================*/

void FILTER_apply ( DATA_SAMPLE *data, COMPLEX *filter, 
                        int nfft, boolean real_filter, COMPLEX *work );

#endif

---