
C-----------------------------------------------------------------------------
      SUBROUTINE CBLOB3
C-----------------------------------------------------------------------------
C     SOUND SPEED PERTURBATION MODEL WITH MULTIPLE BLOBS EACH A
C     MULTIPLICATIVE PERTURBATION WITH EXPONENTIALLY DECAYING 
C     EFFECT IN ALL THREE DIRECTIONS. THE FORM OF THE         
C     PERTURBATION IS: 
C     C = C(1+DEL)    
C     WHERE     DEL = SUM(DELSI)   AND I=1,2,3 FOR EACH BLOB.
C     GIVE LATITUDE INSTEAD OF CO-LATITUDE
	implicit double precision (a-h,o-z)
C     COMMON DECK "CONST" INSERTED HERE    
      COMMON/PCONST/CREF,RGAS,GAMMA        
      COMMON/MCONST/PI,PIT2,PID2,DEGS,RAD,ALN10 
C     COMMON DECK "RKAM" INSERTED HERE  
      DOUBLE PRECISION KR,KTH,KPH
      COMMON//R,TH,PH,KR,KTH,KPH,RKVARS(14),TPULSE,CSTEP,DRDT(20)
C     COMMON DECK "B4" INSERTED HERE 
      INTEGER DCMX,DCNTBL,DCITBL,DCFRMTB,IDSDC(10)   
      COMMON/B4/DCMX,DCNTBL(10),DCITBL(10),DCFRMTB(10),DCGP(10)   
      EQUIVALENCE (DCGP,IDSDC)
C     COMMON DECK "CC" INSERTED HERE
      DOUBLE PRECISION MODC
      COMMON/CC/MODC(4),CS,PCST,PCSR,PCSTH,PCSPH
      COMMON/WW/ dum(10), MAXW,W(1000)
      EQUIVALENCE (EARTHR,W(1))
C     DELTA SOUND SPEED 175-199
      DOUBLE PRECISION C0(3),Z0(3),CBTH0(3),TH0(3),PH0(3)
      DOUBLE PRECISION WZ(3),WTH(3),WPH(3)
      DOUBLE PRECISION FWZ(3),FWTH(3),FWPH(3)
      EQUIVALENCE (PN,W(178)),(C0,W(179)), (Z0,W(182)),(CBTH0,W(185))
      EQUIVALENCE (PH0,W(188)),(WZ,W(191)),(WTH,W(194)),(WPH,W(197))
      INTEGER DFMX    ,DFNTBL(10),DFITBL(10),DFFRMTB(10)
      DATA DFMX/1/ 
      DATA DFNTBL/1,11,8*0/  
      DATA DFITBL/1,9*0/ 
      DATA DFFRMTB/1,9*0/
C 
C 
      ENTRY SETPSP
      DCMX=DFMX 
      CALL IMOVE(DCNTBL,DFNTBL,10)
      CALL IMOVE(DCITBL,DFITBL,10)
      CALL IMOVE(DCFRMTB,DFFRMTB,10)
      RETURN 
C  
      ENTRY IPSPEED
      if (w(175) .ne. 3.0) then
           PRINT*,' Model check number mismatch: '
           STOP   ' W(175) should be 3 if you are using model CBLOB3.'
      endif

      MODC(3)=6HCBLOB3
      MODC(4)= w(177)  ! ID number for a particular model's parameters
      NP=PN
      DO 10 I=1,NP
      FWZ(I)=0.d0 
      FWTH(I)=0.d0
      FWPH(I)=0.d0
      TH0(I)= PID2-CBTH0(I)
      IF(WZ(I).NE.0.d0) FWZ(I)=2.d0/WZ(I)/WZ(I)
      IF(WTH(I).NE.0.d0) FWTH(I)=2.d0/WTH(I)/WTH(I)
      IF(WPH(I).NE.0.d0) FWPH(I)=2.d0/WPH(I)/WPH(I)
10    CONTINUE 
      RETURN 
C
      ENTRY PSPEED
      DEL=1.d0 
C     INITIALIZE PARTIALS OF DEL
      PDELR=0.d0
      PDELTH=0.d0
      PDELPH=0.d0
      DO 20 I=1,NP
      DZ=R-EARTHR-Z0(I)
      DTH=TH-TH0(I)
      DPH=PH-PH0(I)
      DEXPO=0.d0   
      EXPO=-0.5d0*(DZ*DZ*FWZ(I)+DTH*DTH*FWTH(I)+DPH*DPH*FWPH(I)) 
      IF(EXPO .GT. -2.d2) DEXPO=C0(I)*dEXP(EXPO)
      DEL=DEL+DEXPO 
      PDELR=PDELR-DEXPO*FWZ(I)*DZ
      PDELTH=PDELTH-DEXPO*FWTH(I)*DTH
20    PDELPH=PDELPH-DEXPO*FWPH(I)*DPH
C     IN THIS MODEL WE ARE USING THE SQUARE OF DEL 
C     AS THE CS PERTURBATION
      FTR=CS*2.d0*DEL
      DEL=DEL*DEL
      CS=CS*DEL
      PCSR=PCSR*DEL+FTR*PDELR
      PCSTH=PCSTH*DEL+FTR*PDELTH
      PCSPH=PCSPH*DEL+FTR*PDELPH
      RETURN
      END  



