      PROGRAM EDFB
CCC   160630
CCC   IES=1 FOR SURROUNDING EXTERNAL SPHERE CENTERED AT ORIGIN;
CCC   HOMOGENEOUS MATERIAL WITHIN THE EXT. SPHERE IS READ IN 
CCC   AS THAT OF THE (NSHL+1)-TH LAYER OF THE 1-ST SPHERE:
CCC   RCF(1,NSHL)=1.0D0, RCF(1,NSHL+1)>1.0D0 AND
CCC   ROS(1) IS RADIUS OF SPHERE 1
CCC
CCC   READ DATA FOR BUILDING VECTOR XIV FROM WITHIN SUB INXI
CCC   
CCC   IDFC>0 WHEN ALL DIEL. FUNCT. OF SPHERES ARE CONSTANTS;
CCC   IDFC=0 WHEN DIEL. FUNCT. OF SPHERES DEPEND ON XI;
CCC   IDFC<0 WHEN DIEL. FUNCT. OF SPHERES ARE AT XIP VALUE ONLY AND
CCC   XI IS SCALE FACTOR FOR DIMENSIONS
CCC   INSN CHOOSES THE VARIABLE THE DIEL. FUNCT. DEPEND ON:
CCC   (INSN=1)=XI;
CCC   (INSN=2)=WN (THE WAVENUMBER, IN m**-1);
CCC   (INSN=3)=WL (THE WAVELENGTH, IN m);
CCC   (INSN=4)=PU (THE ANGULAR FREQUENCY, IN s**-1);
CCC   (INSN=5)=EV (THE PHOTON ENERGY, IN ev);
CCC   INSTPC>0 WHEN VARIABLE INCREASES WITH A CONSTANT STEP;
CCC   INSTPC=0 WHEN VARIABLE IS SAVED IN A VECTOR YOU READ BEFORE
CCC   SUPPORTS EXPERIMENTAL DIELECTRIC FUNCTIONS ONLY
CCC   NSPH=6
      IMPLICIT REAL*8(A-H,O-Z)
CCC   COMMON/C1/DC0(NSHL-NTL+1),DC0M(NSHL-NTL+1,NSPH,NXI),
CCC  1ROS(NSPH),RCF(NSPH,NSHL+1),IOG(NSPH),NSHL(NSPH)
      COMMON/C1/DC0(5),DC0M(5,6,200),
     1ROS(6),RCF(6,9),IOG(6),NSHL(6)
      COMPLEX*16 DC0,DC0M
CCC   COMMON/C3/XIV(NXI),WNS(NXI),WLS(NXI),PUS(NXI),EVS(NXI),
CCC  1VSS(NXI),VNS(5)
      COMMON/C3/XIV(200),WNS(200),WLS(200),PUS(200),EVS(200),
     1VSS(200),VNS(5)
      CHARACTER*3 VNS
 5010 FORMAT(16I5)
 6005 FORMAT(' SPHERE N.',I4)
 6009 FORMAT(' NONTRANSITION LAYER N.',I2,', SCALE = ',A3)
 6010 FORMAT(I5,1X,1PD12.4,1PD12.4)
      IR=5
      IW=6
      IT=7
CCC
CCC   SIZE(I)=VK*ROS(I) IS SIZE PARAMETER in vacuo
CCC
CCC
CCC   READING OF DIELECTRIC FUNCTIONS DRIVEN BY ICI DEFINED BELOW
CCC
      OPEN(IR,FILE='../../test_data/sphere/DEDFB',STATUS='OLD')
      READ(IR,*)NSPH,IES
      IF(IES.NE.0)IES=1
      READ(IR,*)EXDC,WP,XIP,IDFC,NXI,INSTPC,INSN
      OPEN(IW,FILE='OEDFB',STATUS='UNKNOWN')
      CALL INXI(IR,IW,WP,XIP,IDFC,NXI,INSTPC,INSN)
      READ(IR,5010)(IOG(I),I=1,NSPH)
      DO 113 I=1,NSPH
      IF(IOG(I).LT.I)GO TO 113
      READ(IR,*)NSHL(I),ROS(I)
      NSH=NSHL(I)
      IF(I.EQ.1)NSH=NSH+IES
      DO 112 NS=1,NSH
  112 READ(IR,*)RCF(I,NS)
  113 CONTINUE
      OPEN(IT,FILE='TEDF',FORM='UNFORMATTED',STATUS='UNKNOWN')
      WRITE(IT)NSPH
      WRITE(IT)(IOG(I),I=1,NSPH)
      WRITE(IT)EXDC,WP,XIP,IDFC,NXI
      WRITE(IT)(XIV(I),I=1,NXI)
      DO 115 I=1,NSPH
      IF(IOG(I).LT.I)GO TO 115
      WRITE(IT)NSHL(I),ROS(I)
      NSH=NSHL(I)
      IF(I.EQ.1)NSH=NSH+IES
      WRITE(IT)(RCF(I,NS),NS=1,NSH)
  115 CONTINUE
      DO 468 JXI=1,NXI
      IF((IDFC.NE.0).AND.(JXI.GT.1))GO TO 468
      DO 162 I=1,NSPH
      IF(IOG(I).LT.I)GO TO 162
CCC
      NSH=NSHL(I)
      ICI=(NSH+1)/2
      IF(I.EQ.1)ICI=ICI+IES
      DO 157 IC=1,ICI
      READ(IR,*)DC0(IC)
  157 DC0M(IC,I,JXI)=DC0(IC)
CCC
      WRITE(IT)(DC0(IC),IC=1,ICI)
  162 CONTINUE
  468 CONTINUE
      IF(IDFC.EQ.0)GO TO 474
      WRITE(IW,*)' DIELECTRIC CONSTANTS'
      DO 473 I=1,NSPH
      IF(IOG(I).NE.I)GO TO 473
      ICI=(NSHL(I)+1)/2
      IF(I.EQ.1)ICI=ICI+IES
      WRITE(IW,6005)I
      DO 472 IC=1,ICI
      WRITE(IW,6010)IC,DC0M(IC,I,1)
  472 CONTINUE
  473 CONTINUE
      GO TO 499
  474 WRITE(IW,*)' DIELECTRIC FUNCTIONS'
      DO 478 I=1,NSPH
      IF(IOG(I).NE.I)GO TO 478
      ICI=(NSHL(I)+1)/2
      IF(I.EQ.1)ICI=ICI+IES
      WRITE(IW,6005)I
      DO 477 IC=1,ICI
      WRITE(IW,6009)IC,VNS(INSN)
      DO 476 JXI=1,NXI
      WRITE(IW,6010)JXI,DC0M(IC,I,JXI)
  476 CONTINUE
  477 CONTINUE
  478 CONTINUE
  499 CLOSE(IR)
      CLOSE(IW)
      CLOSE(IT)
      STOP
      END
      SUBROUTINE INXI(IR,IW,WP,XIP,IDFC,NXI,INSTPC,INSN)
      IMPLICIT REAL*8(A-H,O-Z)
      COMMON/C3/XIV(200),WNS(200),WLS(200),PUS(200),EVS(200),
     1VSS(200),VNS(5)
      CHARACTER*3 VNS
      PIGT=DACOS(0.0D0)*4.0D0
      EVC=6.5821188D-16
      IF(IDFC.LT.0)GO TO 300
      IF(INSTPC.EQ.0)GO TO 200
CCC   VLST=V+(NXI-1)*VSTP
      GO TO(105,125,145,165,185),INSN
      RETURN
  105 READ(IR,*)XI,XISTP
      DO 110 JXI=1,NXI
      PU=XI*WP
      WN=PU/3.0D08
      VNS(INSN)='XIV'
      VSS(JXI)=XI
      XIV(JXI)=XI
      PUS(JXI)=PU
      EVS(JXI)=PU*EVC
      WNS(JXI)=WN
      WLS(JXI)=PIGT/WN
  110 XI=XI+XISTP
      WRITE(IW,6601)
 6601 FORMAT
     1(2X,'JXI',5X,'XIV',10X,'WNS',10X,'WLS',10X,'PUS',10X,'EVS')
      WRITE(IW,6600)
     1(JXI,XIV(JXI),WNS(JXI),WLS(JXI),PUS(JXI),EVS(JXI),JXI=1,NXI)
 6600 FORMAT((I5,5(1PD13.4)))
      RETURN
  125 READ(IR,*)WN,WNSTP
      DO 130 JXI=1,NXI
      XI=3.0D08*WN/WP
      PU=XI*WP
      VNS(INSN)='WNS'
      VSS(JXI)=WN
      WNS(JXI)=WN
      WLS(JXI)=PIGT/WN
      XIV(JXI)=XI
      PUS(JXI)=PU
      EVS(JXI)=PU*EVC
  130 WN=WN+WNSTP
      WRITE(IW,6602)
 6602 FORMAT
     1(2X,'JXI',5X,'WNS',10X,'WLS',10X,'PUS',10X,'EVS',10X,'XIV')
      WRITE(IW,6600)
     1(JXI,WNS(JXI),WLS(JXI),PUS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
      RETURN
  145 READ(IR,*)WL,WLSTP
      DO 150 JXI=1,NXI
      WN=PIGT/WL
      XI=3.0D08*WN/WP
      PU=XI*WP
      VNS(INSN)='WLS'
      VSS(JXI)=WL
      WLS(JXI)=WL
      WNS(JXI)=WN
      XIV(JXI)=XI
      PUS(JXI)=PU
      EVS(JXI)=PU*EVC
  150 WL=WL+WLSTP
      WRITE(IW,6603)
 6603 FORMAT
     1(2X,'JXI',5X,'WLS',10X,'WNS',10X,'PUS',10X,'EVS',10X,'XIV')
      WRITE(IW,6600)
     1(JXI,WLS(JXI),WNS(JXI),PUS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
      RETURN
  165 READ(IR,*)PU,PUSTP
      DO 170 JXI=1,NXI
      XI=PU/WP
      WN=PU/3.0D08
      VNS(INSN)='PUS'
      VSS(JXI)=PU
      PUS(JXI)=PU
      XIV(JXI)=XI
      WNS(JXI)=WN
      WLS(JXI)=PIGT/WN
      EVS(JXI)=PU*EVC
  170 PU=PU+PUSTP
      WRITE(IW,6604)
 6604 FORMAT
     1(2X,'JXI',5X,'PUS',10X,'WNS',10X,'WLS',10X,'EVS',10X,'XIV')
      WRITE(IW,6600)
     1(JXI,PUS(JXI),WNS(JXI),WLS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
      RETURN
  185 READ(IR,*)EV,EVSTP
      DO 190 JXI=1,NXI
      PU=EV/EVC
      XI=PU/WP
      WN=PU/3.0D08
      VNS(INSN)='EVS'
      VSS(JXI)=EV
      EVS(JXI)=EV
      PUS(JXI)=PU
      XIV(JXI)=XI
      WNS(JXI)=WN
      WLS(JXI)=PIGT/WN
  190 EV=EV+EVSTP
      WRITE(IW,6605)
 6605 FORMAT
     1(2X,'JXI',5X,'EVS',10X,'WNS',10X,'WLS',10X,'PUS',10X,'XIV')
      WRITE(IW,6600)
     1(JXI,EVS(JXI),WNS(JXI),WLS(JXI),PUS(JXI),XIV(JXI),JXI=1,NXI)
      RETURN
  200 GO TO(205,225,245,265,285),INSN
      RETURN
  205 DO 210 JXI=1,NXI
      READ(IR,*)XI
      PU=XI*WP
      WN=PU/3.0D08
      VNS(INSN)='XIV'
      VSS(JXI)=XI
      XIV(JXI)=XI
      PUS(JXI)=PU
      EVS(JXI)=PU*EVC
      WNS(JXI)=WN
      WLS(JXI)=PIGT/WN
  210 CONTINUE
      WRITE(IW,6601)
      WRITE(IW,6600)
     1(JXI,XIV(JXI),WNS(JXI),WLS(JXI),PUS(JXI),EVS(JXI),JXI=1,NXI)
      RETURN
  225 DO 230 JXI=1,NXI
      READ(IR,*)WN
      XI=3.0D08*WN/WP
      PU=XI*WP
      VNS(INSN)='WNS'
      VSS(JXI)=WN
      WNS(JXI)=WN
      WLS(JXI)=PIGT/WN
      XIV(JXI)=XI
      PUS(JXI)=PU
      EVS(JXI)=PU*EVC
  230 CONTINUE
      WRITE(IW,6602)
      WRITE(IW,6600)
     1(JXI,WNS(JXI),WLS(JXI),PUS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
      RETURN
  245 DO 250 JXI=1,NXI
      READ(IR,*)WL
      WN=PIGT/WL
      XI=3.0D08*WN/WP
      PU=XI*WP
      VNS(INSN)='WLS'
      VSS(JXI)=WL
      WLS(JXI)=WL
      WNS(JXI)=WN
      XIV(JXI)=XI
      PUS(JXI)=PU
      EVS(JXI)=PU*EVC
  250 CONTINUE
      WRITE(IW,6603)
      WRITE(IW,6600)
     1(JXI,WLS(JXI),WNS(JXI),PUS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
      RETURN
  265 DO 270 JXI=1,NXI
      READ(IR,*)PU
      XI=PU/WP
      WN=PU/3.0D08
      VNS(INSN)='PUS'
      VSS(JXI)=PU
      PUS(JXI)=PU
      XIV(JXI)=XI
      WNS(JXI)=WN
      WLS(JXI)=PIGT/WN
      EVS(JXI)=PU*EVC
  270 CONTINUE
      WRITE(IW,6604)
      WRITE(IW,6600)
     1(JXI,PUS(JXI),WNS(JXI),WLS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
      RETURN
  285 DO 290 JXI=1,NXI
      READ(IR,*)EV
      PU=EV/EVC
      XI=PU/WP
      WN=PU/3.0D08
      VNS(INSN)='EVS'
      VSS(JXI)=EV
      EVS(JXI)=EV
      PUS(JXI)=PU
      XIV(JXI)=XI
      WNS(JXI)=WN
      WLS(JXI)=PIGT/WN
  290 CONTINUE
      WRITE(IW,6605)
      WRITE(IW,6600)
     1(JXI,EVS(JXI),WNS(JXI),WLS(JXI),PUS(JXI),XIV(JXI),JXI=1,NXI)
      RETURN
  300 IF(INSTPC.GT.0)GO TO 315
      DO 310 JXI=1,NXI
      READ(IR,*)XI
      VNS(INSN)='XIV'
      VSS(JXI)=XI
      XIV(JXI)=XI
  310 CONTINUE
      GO TO 330
  315 READ(IR,*)XI,XISTP
      DO 320 JXI=1,NXI
      VNS(INSN)='XIV'
      VSS(JXI)=XI
      XIV(JXI)=XI
  320 XI=XI+XISTP
  330 PU=XIP*WP
      WN=PU/3.0D08
      WRITE(IW,6611)
 6611 FORMAT
     1(10X,'XIP',10X,'WN ',10X,'WL ',10X,'PU ',10X,'EV')
      WRITE(IW,6610)XIP,WN,PIGT/WN,PU,PU*EVC
 6610 FORMAT((5X,5(1PD13.4)))
      WRITE(IW,*)' SCALE FACTORS XI'
      WRITE(IW,6612)(JXI,XIV(JXI),JXI=1,NXI)
 6612 FORMAT(I5,1PD13.4)
      RETURN
      END
CCC      
