subroutine xbpf( datas, nsamps, fl, fh, ts )                   
c
c  xbpf (Apply Butterworth band pass filter) 
c  modified by a. tento 03/03 from
C
C     IIR FILTER DESIGN AND IMPLEMENTATION    
C  AUTHOR:  Dave Harris
C  LAST MODIFIED:  September 12, 1990
C 
C  ARGUMENTS:
C  ---------- 
C       
C    datas          REAL ARRAY CONTAINING SEQUENCE TO BE FILTERED
C                   ORIGINAL DATA DESTROYED, REPLACED BY FILTERED DATA
C
C    nsamps         NUMBER OF SAMPLES IN DATA
C
C
C    fl             LOW FREQUENCY CUTOFF OF FILTER (HERTZ)
C
C    fh             HIGH FREQUENCY CUTOFF OF FILTER (HERTZ)
C
C    ts             SAMPLING INTERVAL (SECONDS)
C
C
C INTERNAL PARAMETERS
C -------------------
C
C    IORD           ORDER (#POLES) OF ANALOG PROTOTYPE
C                   NOT TO EXCEED 10 IN THIS CONFIGURATION.  4 - 5
C                   SHOULD BE AMPLE. 
C    PASSES         INTEGER VARIABLE CONTAINING THE NUMBER OF PASSES
C                   1 -- FORWARD FILTERING ONLY
C                   2 -- FORWARD AND REVERSE (I.E. ZERO PHASE) FILTERING
C
        implicit none
        real datas, fl, fh, ts
        integer nsamps
        dimension datas(*)
c
        integer nsects, passes, iord
        real  flw, fhw, sn(30), sd(30), dcvalue, warp
        complex p(10), z(10)
        character*3 stype(10)
        logical zp
c
c
        parameter ( passes = 2,  iord = 2 )
C
C  Filter designed
C
        call buroots( p, stype, dcvalue, nsects, iord ) 
        flw = warp( fl*ts/2., 2. )
        fhw = warp( fh*ts/2., 2. )
        call lptbp( p, z, stype, dcvalue, nsects, flw, fhw, sn, sd )
        call bilin2( sn, sd, nsects )
C
C  Filter data
C
        if (   passes .eq. 1 ) then
          zp = .false.
        else
          zp = .true.
        end if
        call apply( datas, nsamps, zp, sn, sd, nsects )
C
      return
C
      end