haskins logo

SineWave Synthesis

Source Code

C****
C
C Note:
C
C What follows is slightly modified FORTRAN source code
C that is being used as part of the SWS weblet.
C This source code is intended only as an illustration.
C Included is the SWS sinewave sound generation code subroutines.
C Modifications have been made for purposes of illustration.
C (c) 1980-1996, Haskins Laboratories, New Haven, CT.
C
C****
C
C SUBROUTINE SWS_GENV1
C
C P. RUBIN 2/26/80, 2/5/81
C 6/16/81 FOR VAX
C 2/13/82
C 10/27/82 CORRECT HEADER WRITE
C 5/26/83 HANDLE PREEMPHASIS,CONTIG
C 6/19/84 CHANGES FOR AMP_TYPE C C PURPOSE:
C
C GENERATE PCM DATA FROM ARRAYS OF TIME SLICE (TIMSL),
C SINE WAVE FREQUENCY (SWF) AND SINE WAVE AMPLITUDE (SWA)
C INFORMATION. C C VOICE 1 : SINE WAVE
C
C ARGUMENTS PASSED IN COMMON:
C
C OUTCHN LOGICAL OUTPUT CHANNEL FLAG, WHERE:
C .TRUE. = CHANNEL ASSIGNED C .FALSE.= CHANNEL NOT ASSIGNED
C
C ERR INTEGER ERROR FLAG, WHERE:
C 0 = NO ERROR
C <0 = ERROR
C
C NOTES:
C
C FREQUENCY IS CONSTRAINED TO BE BETWEEN 1 HZ AT THE MINIMUM
C AND THE NYQUIST LIMIT ( (SAMPLING FREQ./2) - 1 ) AT THE MAX.
C
C================================================
PARAMETER MXRECL = 64 ! MAX. RECORD LENGTH
PARAMETER HEADER = 4 ! # HEADER RECORDS
PARAMETER MIDSCL = 2048 ! PCM MIDSCALE
PARAMETER PI = 3.14159 ! PI
PARAMETER TWOPI = 6.28319 ! 2. * PI
C===============================================
C
C*****
C
C INCLUDE 'SWS.INC/NOLIST' ! MAIN SWS COMMON
C
C*****
C*****
C S W S . I N C C C P. RUBIN 9/16/80
C 7/21/81 FOR VAX C 6/15/84
C 7/18/84 CHANGES FOR MAC AND INC FEATURES
C 8/8/84 C 9/19/85 V4.1 WITH MANY MORE SLICES
C 3/4/88 V4.2 supports DATA TRANSLATION output
C 10/15/91 changes for additional output systems
C 6/4/92 changes for AUDIO flag C C SINEWAVE SYNTHESIZER COMMON
C*****
PARAMETER MAXSW = 50 ! MAX # SINE WAVES
PARAMETER MAXSL = 2000 ! MAX # TIME SLICES
PARAMETER MAXDISP= 999 ! MAX. # TIME SLICES TO DISPLAY
PARAMETER CTRLZ = -99 ! CONTROL-Z CODE
INTEGER MAXMDEF ! MAX. # OF MACRO DEFINITIONS
PARAMETER ( MAXMDEF = 100 )
COMMON /SWSM/ LIST,TTO,TTI,LP,CFLUN,SPFIL,DTFIL,PCFIL,
1 SWF(0:MAXSL,MAXSW),SWA(0:MAXSL,MAXSW),
2 TIMSL(0:MAXSL),NSLIS,NSW,SR,DEF_SR,
3 NRECS,NHEDER,TEXT1,LTEX1,TEXT2,LTEX2,
4 NCODE, CODEC, CODEV, PROMPT, DEFV, ICHAN,
5 SWDAT, SPOPEN, DISPLA, OUTCHN, PCMDAT,
6 LNSP, FAXIS, TAXIS, NCKPAG, TIK_FLAG, ERR,
7 IVOICE, VOIARG(12), EMPH, LPCMEXT, SC_CHAN,
8 AMP_MAX, NAMESP, PCMEXT, AMP_TYPE, AUDIO,
9 INI_LUN, ARGS, FPARS, NUMMDEF, MDEF, LMDEF,
1 PCM_DT INTEGER LIST ! LISTING DEVICE LUN

INTEGER TTO ! TERMINAL OUT LUN
INTEGER TTI ! TERMINAL IN LUN
INTEGER LP ! LINEPRINTER LUN
INTEGER CFLUN ! COMMAND FILE LUN
INTEGER SPFIL ! SPEECH FILE LUN
INTEGER DTFIL ! SW DATA FILE LUN
INTEGER PCFIL ! PCM DATA FILE LUN REAL SWF ! SINE WAVE FREQUENCIES
REAL SWA ! SINE WAVE AMPLITUDES
REAL TIMSL ! TIME SLICES
INTEGER NSLIS ! # SLICES
INTEGER NSW ! # SINE WAVES
REAL SR ! SAMPLING RATE
REAL DEF_SR ! SAMPLING RATE DEFAULT
INTEGER NRECS ! # RECS IN PCM FILE
INTEGER NHEDER ! # HEADER BLOCKS IN PCM
CHARACTER*100 TEXT1 ! TEXT 1
INTEGER LTEX1 ! LENGTH OF TEXT 1
CHARACTER*100 TEXT2 ! TEXT 2
INTEGER LTEX2 ! LENGTH OF TEXT 2
INTEGER NCODE ! # INPUT CODES
BYTE CODEC(2) ! CHARACTER CODES
REAL CODEV(2) ! CODE VALUES
REAL PROMPT ! PROMPT VALUE
REAL DEFV ! DEFAULT VALUE
INTEGER ICHAN ! CHANNEL NUMBER
LOGICAL SWDAT ! SINEWAVE DATA FLAG
LOGICAL SPOPEN ! SPEECH FILE OPENED ?
LOGICAL DISPLA ! DISPLAY ALLOWED ?
LOGICAL OUTCHN ! OUTPUT CHANNEL ASSIGNED?
LOGICAL PCMDAT ! PCM DATA EXIST ?
INTEGER LNSP ! LENGTH OF NAMESP
REAL FAXIS ! FREQ. AXIS SCALE DEF.
REAL TAXIS ! TIME AXIS SCALE DEF.
INTEGER NCKPAG ! LINE # FOR CKPAGE
LOGICAL TIK_FLAG ! TIME TICKS FOR DISPLAY
INTEGER ERR ! ERROR FLAG
INTEGER IVOICE ! VOICE QUALITY
REAL VOIARG ! VOICE QUALITY ARGUMENTS
LOGICAL*1 EMPH ! PRE-EMPHASIS FLAG
INTEGER LPCMEXT ! LENGTH OF PCM EXTENSION
INTEGER*2 SC_CHAN ! SCA0 DEVICE CHANNEL
REAL AMP_MAX ! MAX. AMPLITUDE VALUE
CHARACTER*100 NAMESP ! NAME OF SWS PCM FILE
CHARACTER*4 PCMEXT ! PCM EXTENSION
CHARACTER*1 AMP_TYPE ! AMPLITUDE: LINEAR OR DB
CHARACTER*4 AUDIO ! AUDIO SYSTEM !
NONE, DT, GRAD or SIOP
INTEGER INI_LUN ! INITIALIZATION FILE LUN
LOGICAL*1 ARGS ! .TRUE. IF AT LEAST 1 ARGUMENT
REAL FPARS(12) ! F.P. ARGUMENTS
INTEGER NUMMDEF ! # OF MACROS CURRENTLY DEFINED
CHARACTER MDEF(MAXMDEF)*100 ! MACRO DEFINITIONS
INTEGER LMDEF(MAXMDEF) ! LENGTH OF EACH MACRO DEF.
LOGICAL PCM_DT ! Data Translation PCM flag
C.....
C
C End of SWS.INC
C
C============================================
EXTERNAL CONTIG
REAL FREQW(MAXSW) ! working frequencies
REAL AMPLW(MAXSW) ! working amplitudes
REAL RADIAN(MAXSW) ! angle of sine waves in radians
REAL FNYQL ! Nyquist freq. limit
INTEGER LEN ! total # samples
INTEGER NSMPIN ! total # samples
INTEGER FOR_STAT INTEGER*2 IB(MXRECL) ! record buffer
LOGICAL RAMPA(MAXSW) ! amplitude ramp (true or false)
LOGICAL RAMPF(MAXSW) ! frequency ramp (true or false)
EQUIVALENCE (IB(2),LEN) ! put len in IB(2) and (3)
C:::::
C
C Begin routine
C
C:::::
ERR = 0 ! error flag
FNYQL = ( SR/2. ) - 1. ! Nyquist frequency
IF ( .NOT. SWDAT ) THEN
WRITE (TTO,15)
15 FORMAT ('0 *** SINE WAVE DATA DOES NOT EXIST !! ***',/)
ERR = -1
GO TO 9000
END IF
C.....
C Open new disk file to write pcm data to.
C.....
LEN = TIMSL(NSLIS) * ( SR/1000. ) + .5 ! length
NRECS = ( LEN + (MXRECL-1) ) / MXRECL ! # PCM records
NRECS = NRECS + 1 ! extra leeway
NDB = ( NRECS+3 ) / 4 + 1 ! # disk blocks
OPEN (UNIT=PCFIL,NAME=NAMESP,TYPE='NEW',ACCESS='DIRECT',
1 USEROPEN=CONTIG,
2 RECORDSIZE=32,INITIALSIZE=NDB,ERR=50)
SPOPEN = .TRUE.
GO TO 100
C.....
C FILE OPEN ERROR
C.....
50 WRITE (TTO,60)
60 FORMAT ('0 ***
WRITE PCM OPEN DISK FILE ERROR !! ***')
WRITE (TTO,62) PCFIL,NAMESP(1:LNSP)
62 FORMAT (' *** FILE ON LUN ',I2,' IS: ',A,' ***',/)
ERR = -1
CLOSE (UNIT=PCFIL)
GO TO 9000
C.....
C Generate the sine waves from time slice to time slice
C.....
100 IBPTR = 0 ! Buffer pointer
NRECS = 0 ! # of PCM records
NSMPIN = 0
IF ( NSLIS .LT. 2 ) THEN ! Too few slices to work with
WRITE (TTO,110)
110 FORMAT ('0 *** TOO FEW TIME SLICES !! ***',/)
ERR = -1
GO TO 9000
END IF
DO I = 1, NSW ! NSW is the # of sinewaves
RADIAN(I) = 0. ! Initialize radians
END DO
C.....
C If the first time slice does not start at 0. msec.,
C output midscale data.
C.....
IF ( TIMSL(1) .EQ. 0. ) GO TO 250
NSMSL = TIMSL(1) * ( SR/1000. )
IF ( NSMSL .LT. 1 ) GO TO 250
DO J = 0, NSMSL-1
IBPTR = IBPTR + 1 IB(IBPTR) = MIDSCL
IF ( IBPTR .GE. MXRECL ) THEN
NRECS = NRECS + 1
WRITE (PCFIL' NRECS+HEADER) IB
IBPTR = 0
END IF
END DO
NSMPIN = NSMPIN + NSMSL
250 DO 1000 I = 1,NSLIS-1
IF ( TIMSL(I+1) .LE. TIMSL(I) ) GO TO 1000
NSMSL = ( TIMSL(I+1) - TIMSL(I) ) * (SR/1000.)
IF ( NSMSL .LT. 1 ) THEN
WRITE (TTO,270) TIMSL(I),TIMSL(I+1)
270 FORMAT ('0 *** DURATION FROM ',F10.3,' MSEC TO ',F10.3,
1 ' MSEC ')
WRITE (TTO,272)
272 FORMAT (' IS TOO SMALL FOR THIS SAMPLING RATE ',
1 14X,'***',/)
ERR = -1
GO TO 1100
END IF
DO J = 1, NSW
FREQW(J) = SWF(I,J) ! working frequency
AWORK = SWA(I,J) ! working amplitude
C.....
C Amplitude data is specified, as linear from 0 to 1.0.
C A flag, called AMP_TYPE, can be set so that the
C amplitude values are specified in DB.
C If this flag is set, the routine SWS_DB is called.
C.....
IF ( AMP_TYPE .EQ. 'D' ) CALL SWS_DB ('L',AWORK,AWORK)
AMPLW(J) = AWORK ! working amplitude
RAMPA(J) = .FALSE.
RAMPF(J) = .FALSE.
IF ( NSMSL .GE. 2 ) THEN
IF ( SWF(I,J) .NE. SWF(I+1,J) ) RAMPF(J) = .TRUE.
IF ( SWA(I,J) .NE. SWA(I+1,J) ) RAMPA(J) = .TRUE.
END IF
END DO
DO J = 0, NSMSL-1
XSUM = 0.
DO K = 1, NSW
C..... PERCNT is the percent advance into the current slice
IF ( RAMPA(K) .OR. RAMPF(K) ) PERCNT = FLOAT (J) /
1 FLOAT (NSMSL-1)
C..... If changing, compute new frequency and amplitude.
IF ( RAMPA(K) ) THEN
AWORK = SWA(I,K)
AWORK2 = SWA(I+1,K)
IF ( AMP_TYPE .EQ. 'D' ) THEN
CALL SWS_DB ('L',AWORK,AWORK)
CALL SWS_DB ('L',AWORK2,AWORK2)
END IF
AMPLW(K) = AWORK + PERCNT * ( AWORK2 - AWORK )
END IF
IF ( RAMPF(K) ) FREQW(K) = SWF(I,K) + PERCNT *
1 ( SWF(I+1,K)-SWF(I,K) )
C.....
C Window frequency: MIN = 1 HZ; MAX = NYQUIST limit
C.....
IF ( FREQW(K) .LT. 1. ) FREQW(K) = 1. IF ( FREQW(K) .GE. FNYQL ) FREQW(K) = FNYQL
C.....
C Generate
C..... RADIAN(K) = RADIAN(K) + TWOPI * FREQW(K) / SR
IF ( RADIAN(K) .GE. TWOPI ) RADIAN(K) =
1 AMOD(RADIAN(K),TWOPI)
X = AMPLW(K) * SIN( RADIAN(K) )
XSUM = XSUM + X
END DO
IBPTR = IBPTR + 1
IOUT = IFIX ( XSUM/NSW * 2047. ) + MIDSCL
IB(IBPTR) = IOUT
IF ( IBPTR .GE. MXRECL ) THEN
NRECS = NRECS + 1
WRITE (PCFIL' NRECS+HEADER)
IB IBPTR = 0
END IF
END DO
NSMPIN = NSMPIN + NSMSL
1000 CONTINUE
1100 IF ( IBPTR .NE. 0 ) THEN
DO I = IBPTR+1, MXRECL
IB(I) = MIDSCL
END DO
NRECS = NRECS + 1
WRITE (PCFIL' NRECS+HEADER) IB
END IF
C.....
C Rewrite header block.
C.....
DO I = 1, MXRECL
IB(I) = 0
END DO
DO IREC = 2, 4
WRITE (PCFIL' IREC) IB
END DO
LEN = NSMPIN
IB(1) = 1
IB(4) = IFIX(SR)
IF ( .NOT. EMPH ) IB(5) = 1
IB(62) = IB(4)
IB(63) = -32000
NHEDER = HEADER
WRITE (PCFIL' 1) IB
9000 IF ( SPOPEN ) CLOSE (UNIT=PCFIL)
SPOPEN = .FALSE.
RETURN
END

SUBROUTINE SWS_DB (TYPE,FIN,FOUT)
C
C P. RUBIN 6/19/84
C
C
C PURPOSE:
C
C Convert between linear (0-1.) and C db (0-100db) values.
C
C
C ARGUMENTS:
C
C TYPE CHARACTER*1 type of conversion, where:
C 'L' = convert to linear
C 'D' = convert to DB
C FIN REAL value to be converted C FOUT REAL converted value
C C============================================
PARAMETER SCALE = 100000. ! scale for 0 to 1
PARAMETER FACDB = 50. ! constant factor
REAL FIN ! value to be converted
REAL FOUT ! converted value C
HARACTER*1 TYPE ! type of conversion
C:::::
C
C BEGIN ROUTINE
C
C:::::
IF ( TYPE .EQ. 'D' ) THEN
FIN = FIN * SCALE ! scale value up
IF ( FIN .EQ. 0. ) THEN
FOUT = 0. ! convert to db
ELSE
FOUT = 20. * ALOG10( FIN ) F
OUT = FOUT - FACDB ! subtract constant
IF ( FOUT .LT. 0. ) FOUT = 0.
END IF
ELSE
FTMP = FIN + FACDB ! add constant factor
FOUT = 10.** (FTMP/20.) ! convert to linear
FOUT = FOUT / SCALE
END IF
RETURN
END