home *** CD-ROM | disk | FTP | other *** search
/ Oakland CPM Archive / oakcpm.iso / sigm / vol152 / woof.for < prev   
Encoding:
Text File  |  1984-04-29  |  5.9 KB  |  185 lines
  1.       PROGRAM WOOF
  2. C Program to plot the frequency response of a loudspeaker enclosure
  3. C given the speaker's parameters. From IEEE Trans. Electroacoustics.
  4.       DIMENSION RESP(49),ZEL(49),JY(4),JZ(4),J(6),MATRIX(49,29),MK(3)
  5.      1,J1(6)
  6.       REAL M,M1,M2,M2OPT
  7.       DATA A1,V,FAR,S1,BL,R,Q0,M2/.051,.45,40.6,1777.,7.9,5.5,
  8.      15.09,.013/
  9.       DATA MK(1),MK(2),MK(3)/1H ,1H+,1H./
  10.       DATA JY/-10,-20,-30,-40/
  11.       DATA JZ/30,20,10,0/
  12.       DATA J/0,40,80,120,160,200/
  13.       DATA J1/0,10,20,30,40,50/
  14.     6 WRITE (3, 90)
  15.    90 FORMAT('1LOUDSPEAKER ENCLOSURE ANALYSIS PROGRAM.'/
  16.      1'0For information about this program contact'/
  17.      2' Trevor MARSHALL, SYSOP, T.O. Tech RBBS'/
  18.      3'           805-492-5472, (voice)805-492-3693')
  19.        WRITE (3, 567)
  20.   567 FORMAT(' '/' For HELP type 1,otherwise 0'/)
  21.       READ (3,568) IHELP
  22.   568 FORMAT(I1)
  23.       IF(IHELP.EQ.1)GO TO 676
  24.       WRITE (3, 452)
  25.   452 FORMAT(/
  26.      1'0WHAT IS THE EFFECTIVE PISTON AREA IN SQUARE METRES?')
  27.       READ (3,1)A1
  28.       WRITE (3, 91)
  29.    91 FORMAT(' ADIABATIC ENCLOSURE VOLUME IN CUBIC METRES?'/
  30.      1'(SET TO 99999 IF AN INFINITE BAFFLE IS IN USE)')
  31.       READ (3,1)V
  32.       WRITE (3, 92)
  33.    92 FORMAT(' WOOFER FREE AIR RESONANCE?')
  34.       READ (3,1)FAR
  35.       WRITE (3, 93)
  36.    93 FORMAT(' SUSPENSION STIFFNESS IN NEWTON/METRE?')
  37.       READ (3,1)S1
  38.       WRITE (3, 94)
  39.    94 FORMAT(' B L PRODUCT IN WEBERS/METER?')
  40.       READ (3,1)BL
  41.       WRITE (3, 95)
  42.    95 FORMAT(' VOICE COIL RESISTANCE?')
  43.       READ (3,1)R
  44.       WRITE (3, 96)
  45.    96 FORMAT(' SUSPENSION Q FACTOR?')
  46.       READ (3,1)Q0
  47.       WRITE (3, 97)
  48.    97 FORMAT(' VENT AIR MASS IN KILOGRAMS'/
  49.      1' (SET TO 9999 IF ENCLOSURE IS SEALED)')
  50.       READ (3,1)M2
  51.     1 FORMAT(F10.0)
  52.       IF(BL.EQ.0.)STOP
  53.   676 CONTINUE
  54. C      WRITE (3,  999)
  55.  999  FORMAT(1H1)
  56.       WRITE (3,  10)
  57.   10  FORMAT(21X,26HWOOFER PERFORMANCE PLOTTER  /20X,28(1H'),/'0THIS PRO
  58.      1GRAM CALCULATES IN TWO FREQUENCY RANGES,'/' 0 TO 200 HZ IN 4 HZ STE
  59.      2PS '/' AND 0 TO 50 HZ IN 1 HZ STEPS.')
  60.       WRITE (3, 20)A1,V,FAR,S1,BL,R,Q0,M2
  61.    20 FORMAT(4X,24HSYSTEM INPUT DATA ARE...  /
  62.      228H EFFECTIVE PISTON AREA       ,F6.3,17H SQUARE METRES   /
  63.      328H ADIABATIC ENCLOSURE VOLUME  ,F6.3,17H CUBIC METRES    /
  64.      428H WOOFER FREE AIR RESONANCE   ,F6.1,17H HERTZ           /
  65.      528H SUSPENSION STIFFNESS        ,F6.0,17H NEWTON/METRE    /
  66.      628H B L PRODUCT                 ,F6.1,17H WEBER/METRE     /
  67.      728H VOICE COIL RESISTANCE       ,F6.1,17H OHMS            /
  68.      828H SUSPENSION Q FACTOR         ,F6.2,17H                 /
  69.      928H REFLEX VENT AIR MASS        ,F6.3,' KILOGRAM')
  70.       IF(M2.GT.999.) WRITE (3, 2222)
  71. 2222  FORMAT(' THIS ENCLOSURE WILL BE TREATED AS FULLY SEALED (INFINITELY
  72.      1Y SMALL VENT).  ')
  73. CONCLUDES INPUT CHECK
  74. CALCULATE PARAMETERS
  75.       M1=S1/(6.283*FAR)**2
  76.       S2=143000.*A1**2/V
  77.       M2OPT=M1*S2/S1
  78.       REM=BL**2/R
  79.       FH=0.159*SQRT(S2/M2)
  80.       ASYMP=0.0552*REM*(A1/M1)**2
  81.       X1=SQRT(S1*M1)
  82.       R1=X1/Q0
  83.       Q=X1/(R1+REM)
  84.       S=S2/S1
  85.       M=M2/M1
  86.       A=Q**(-2.)-2.-2.*S-2.*S/M
  87.       B=1.+2.*S+S*S+4.*S/M+2.*S*S/M+S*S/M/M-2.*S/Q/Q/M
  88.       C=S*S/Q/Q/M/M-2.*S/M-2.*S*S/M-2.*S*S/M/M
  89.       D=S*S/M/M
  90.       WRITE (3,  30)M1,S2,REM,FH,M2OPT,ASYMP,Q,S,M,A,B,C,D
  91.   30  FORMAT(4X,43HFROM THESE DATA, CONSTANT PARAMETERS ARE... /
  92.      228H MASS OF CONE AND AIR LOAD    ,F6.3,17H KILOGRAM         /
  93.      328H ENCLOSURE AIR STIFFNESS     ,F6.0,17H NEWTON/METRE      /
  94.      428H ELECTRODYNAMIC DRAG         ,F6.2,17H NEWTON SEC/METRE  /
  95.      528H HELMHOLTZ FREQUENCY         ,F6.1,17H HERTZ             /
  96.      628H SUGGESTED VENT AIR MASS     ,F6.3,17H KILOGRAM          /
  97.      728H ASYMPTOTIC EFFICIENCY     ,F6.3,17H PERCENT             /
  98.      8' DIAGNOSTIC PARAMETERS Q,S,M,A,B,C,D,ARE'7(1X,F5.2))
  99. COMMENCE EVALUATING RESPONSE AND IMPEDANCE
  100.       DO 998 ISTEP=1,2
  101. C      IF(ISTEP.EQ.2)WRITE (3, 999)
  102.       DO100 I=1,49
  103.       IF(ISTEP.EQ.1)FR=4*I
  104.       IF(ISTEP.EQ.2)FR=I
  105.       W=6.283*FR
  106.       RRAD=.022*(FR*A1)**2
  107.       R2=RRAD
  108.     FR=W*M2
  109.     FI=FR-S2/W
  110.     CALL CDIV(R2,FR,R2,FI,FR1,FI1)
  111.     ZR=(R1+RRAD)
  112.     ZI=(W*M1-S1/W)
  113.     CALL CMULT(0.,1.,FR1,FI1,FR2,FI2)
  114.     ZR=ZR-FR2*S2/W
  115.     ZI=ZI-(FI2*S2/W)
  116.     ZR2=ZR+REM
  117.     RESP(I)=RRAD*REM*CABS(FR1,FI1)**2/CABS(ZR2,ZI)**2
  118.     CALL CDIV(REM,0.,ZR,ZI,ZR1,ZI1)
  119.     ZR1=ZR1+1.
  120.     ZEL(I)=R*CABS(ZR1,ZI1)
  121. 100    CONTINUE
  122. CLEAR AND FILL DISPLAY MATRIX
  123.       DO 200 IFREQ=1,49
  124.       DO150 LEVEL=1,29
  125.  150  MATRIX(IFREQ,LEVEL)=MK(1)
  126.       IY=40.5+10.*ALOG10(RESP(IFREQ))
  127.       IF((IY.LT.1).OR.(IY.GT.29))GO TO 162
  128.   162 MATRIX(IFREQ,IY)=MK(2)
  129.  170  IZ=0.5+ZEL(IFREQ)
  130.       IF((IZ.LT.1).OR.(IZ.GT.29))GO TO 200
  131.  172  MATRIX(IFREQ,IZ)=MK(3)
  132.  200  CONTINUE
  133. C   NOW WRITE  DISPLAY MATRIX
  134.       WRITE (3,  600)
  135.  600  FORMAT('1',8X,'ABSOLUTE RESPONSE(+)IN DB,AND IMPEDANCE(.)IN OHMS')
  136.  
  137.       WRITE (3, 601)
  138.  601  FORMAT(6(9X,1H+))
  139.       WRITE (3, 602)JY(1),JZ(1)
  140.  602  FORMAT(I9,51(1H+),I3)
  141.       DO 670 LOSS=1,29
  142.       IY=30-LOSS
  143.       IF(LOSS-10)615,616,615
  144.  615  IF(LOSS-20) 617,618,617
  145.  617  WRITE (3,  620)(MATRIX(IFREQ,IY),IFREQ=1,49)
  146.  620  FORMAT(9X,1H+,49A1,1H+)
  147.       GO TO 670
  148.  616  WRITE (3,  630)JY(2),(MATRIX(IFREQ,IY),IFREQ=1,49),JZ(2)
  149.  630  FORMAT(I9,1H+,49A1,1H+,I3)
  150.       GO TO 670
  151.  618  WRITE (3,  630)JY(3),(MATRIX(IFREQ,IY),IFREQ=1,49),JZ(3)
  152.   670 CONTINUE
  153.       WRITE (3, 602)JY(4),JZ(4)
  154.       WRITE (3,  601)
  155.       IF(ISTEP.EQ.2)WRITE (3,  680)J1
  156.       IF(ISTEP.EQ.1)WRITE (3,  680)J
  157.  680  FORMAT(6(7X,I3))
  158.       WRITE (3, 690)
  159.  690  FORMAT(25X,17HFREQUENCY,  HERTZ   ///)
  160.  998  CONTINUE
  161.       WRITE (3, 1000)
  162.  1000 FORMAT(' TYPE 0 TO FINISH, 1 FOR ANOTHER RUN')
  163.       READ (3,1001)IRUN
  164.  1001 FORMAT(I1)
  165.       IF(IRUN.EQ.1)GO TO 6
  166.       STOP
  167.     END
  168.     FUNCTION CABS(X,Y)
  169.     CABS=SQRT(X*X+Y*Y)
  170.     RETURN
  171.     END
  172.     SUBROUTINE CMULT(OPR1R,OPR1I,OPR2R,OPR2I,RLTR,RLTI)
  173.     RLTR=OPR1R*OPR2R-OPR1I*OPR2I
  174.     RLTI=OPR1R*OPR2I+OPR1I*OPR2R
  175.     RETURN
  176.       END
  177.     SUBROUTINE CDIV(OPR1R,OPR1I,OPR2R,OPR2I,RLTR,RLTI)
  178.     DIV9=OPR2R*OPR2R+OPR2I*OPR2I
  179.     RLTR=OPR1R*OPR2R-OPR1I*(0.-OPR2I)
  180.     RLTI=OPR1R*(0.-OPR2I)+OPR1I*OPR2R
  181.     RLTR=RLTR/DIV9
  182.     RLTI=RLTI/DIV9
  183.     RETURN
  184.     END
  185.