!
! This file is part of SACAMOS, State of the Art CAble MOdels for Spice. 
! It was developed by the University of Nottingham and the Netherlands Aerospace 
! Centre (NLR) for ESA under contract number 4000112765/14/NL/HK.
! 
! Copyright (C) 2016-2018 University of Nottingham
! 
! SACAMOS is free software: you can redistribute it and/or modify it under the 
! terms of the GNU General Public License as published by the Free Software 
! Foundation, either version 3 of the License, or (at your option) any later 
! version.
! 
! SACAMOS is distributed in the hope that it will be useful, but 
! WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 
! or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License 
! for more details.
! 
! A copy of the GNU General Public License version 3 can be found in the 
! file GNU_GPL_v3 in the root or at <http://www.gnu.org/licenses/>.
! 
! SACAMOS uses the EISPACK library (in /SRC/EISPACK). EISPACK is subject to 
! the GNU Lesser General Public License. A copy of the GNU Lesser General Public 
! License version can be found in the file GNU_LGPL in the root of EISPACK 
! (/SRC/EISPACK ) or at <http://www.gnu.org/licenses/>.
! 
! The University of Nottingham can be contacted at: ggiemr@nottingham.ac.uk
!
!
! File Contents:
! SUBROUTINE spacewire_set_parameters
! SUBROUTINE spacewire_set_internal_domain_information
! SUBROUTINE spacewire_plot
!
! NAME
!     spacewire_set_parameters
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     Set the overall parameters for a spacewire cable
!
! COMMENTS
!      
!
! HISTORY
!
!     started 5/9/2016 CJS based on spacewire.F90
!     16/11/2017 CJS Include network synthesis process to replace s-domain transfer functions
!
!
SUBROUTINE spacewire_set_parameters(cable)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  type(cable_specification_type),intent(INOUT)    :: cable

! local variables

! START

  cable%cable_type=cable_geometry_type_spacewire
  cable%tot_n_conductors=13
  cable%tot_n_domains=10
  cable%n_external_conductors=1
  cable%n_internal_conductors=12
  cable%n_internal_domains=9
  cable%n_parameters=13
  cable%n_dielectric_filters=3
  cable%n_transfer_impedance_models=2

END SUBROUTINE spacewire_set_parameters
!
! NAME
!     spacewire_set_internal_domain_information
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     Set the overall parameters for a spacewire cable
!
! COMMENTS
!     Set the dimension of the domain transformation matrices to include an external reference conductor for the cable 
!      
!
! HISTORY
!
!     started 12/4/201 CJS 
!     8/9/2016 CJS common mode/ differential mode loss correction
!     19/9/2016 CJS frequency dependent dielectric in Laplace solver
!     2/11/2016 CJS inhomogeneous dielectric in twisted pair model
!     8/5/2017         CJS: Include references to Theory_Manual
!
!
SUBROUTINE spacewire_set_internal_domain_information(cable)

USE type_specifications
USE constants
USE general_module
USE maths
USE PUL_parameter_module

IMPLICIT NONE


! variables passed to subroutine

  type(cable_specification_type),intent(INOUT)    :: cable

! local variables

  integer :: n_conductors
  integer :: dim

  integer  :: domain
  integer  :: inner_cable
  real(dp) :: d,ctheta

  real(dp) :: L11,L12,L13
  real(dp) :: C11,C12
  
  real(dp) :: LC,LD,CC,CD
  logical  :: dielectric_is_homogeneous
  
  integer  :: conductor_1,conductor_2,reference_conductor
  integer  :: diff,com,inner_shield,outer_shield
 
  type(PUL_type)    :: PUL
  
  real(dp) :: C_air
  type(Sfilter) :: jw
 
  integer :: ierr,i
  
  real(dp) :: epsr
 
! variables for cable parameter checks 
  logical :: cable_spec_error

  real(dp) :: rw
  real(dp) :: rd
  real(dp) :: s
  real(dp) :: rs
  real(dp) :: rd2
  real(dp) :: stpr
  real(dp) :: rs2
  real(dp) :: rd3
  real(dp) :: t1
  real(dp) :: t2
  real(dp) :: sigma_w
  real(dp) :: sigma_s1
  real(dp) :: sigma_s2

  type(Sfilter) :: epsr1,epsr2,epsr3,ZT1,ZT2
  type(Sfilter) :: YC,YD
  
  character(LEN=error_message_length) :: message 
  
! START

! Check the cable parameters

  rw=cable%parameters(1)
  rd=cable%parameters(2)
  s=cable%parameters(3)
  rs=cable%parameters(4)
  t1=cable%parameters(5)  
  rd2=cable%parameters(6)
  stpr=cable%parameters(7)
  rs2=cable%parameters(8)
  t2=cable%parameters(9)  
  rd3=cable%parameters(10)
  
  sigma_w=cable%parameters(11)
  sigma_s1=cable%parameters(12)
  sigma_s2=cable%parameters(13)
  
  epsr1=cable%dielectric_filter(1)
  epsr2=cable%dielectric_filter(2)
  epsr3=cable%dielectric_filter(3)
  ZT1=cable%transfer_impedance(1)
  ZT2=cable%transfer_impedance(2)
  
  write(*,*)'CHECKING FOR ERRORS'
  
  cable_spec_error=.FALSE.    ! assume no errors initially
  message=''
  CALL spacewire_check(rw,rd,s,rs,rd2,stpr,rs2,rd3,cable_spec_error,cable%cable_name,message)
  CALL conductivity_check(sigma_w,cable_spec_error,cable%cable_name,message)
  CALL conductivity_check(sigma_s1,cable_spec_error,cable%cable_name,message)
  CALL conductivity_check(sigma_s2,cable_spec_error,cable%cable_name,message)
  CALL dielectric_check(epsr1,cable_spec_error,cable%cable_name,message)
  CALL dielectric_check(epsr2,cable_spec_error,cable%cable_name,message)
  CALL dielectric_check(epsr3,cable_spec_error,cable%cable_name,message)
  CALL transfer_impedance_check(ZT1,cable_spec_error,cable%cable_name,message)
  CALL transfer_impedance_check(ZT2,cable_spec_error,cable%cable_name,message)
  CALL surface_impedance_check(ZT1,sigma_s1,rs,t1,cable_spec_error,cable%cable_name,message)
  CALL surface_impedance_check(ZT2,sigma_s2,rs2,t2,cable_spec_error,cable%cable_name,message)
    
  if (cable_spec_error) then       
    run_status='ERROR in cable_model_builder, error on parameters for cable:'//trim(cable%cable_name)//'. '//trim(message)
    CALL write_program_status()
    STOP 1
  end if
  
  write(*,*)'Use laplace=',use_laplace
  
! pre-calculate inductance matrix elements for two conductors in a cylindrical shield
! See C.R. Paul, 1st edition, equation 3.67a,b with cos(thetaij)=-1

  epsr=evaluate_Sfilter_high_frequency_limit(cable%dielectric_filter(1))
  
  jw=jwA_filter(1d0)
  
  domain=1

  if (use_laplace) then
  
! allocate memory for the PUL parameter solver interface

    if(verbose) write(*,*)'Domain:',domain
    if(verbose) write(*,*)'Allocating PUL data structure for shielded twisted pairs'
    n_conductors=3
    
    CALL allocate_and_reset_PUL_data(PUL,n_conductors)
    
    PUL%shape(1:n_conductors)=circle
    
    PUL%x(1)=-s/2d0
    PUL%y(1)=0.0
    PUL%r(1)=rw
    PUL%rd(1)=rd
    PUL%epsr(1)=epsr1
    
    PUL%x(2)=s/2d0
    PUL%y(2)=0.0
    PUL%r(2)=rw
    PUL%rd(2)=rd
    PUL%epsr(2)=epsr1
    
    PUL%epsr_background = 1d0  ! permittivity of homogeneous dielectric medium surrounding the insulated conductors (air)
  
! no ground plane
    PUL%ground_plane_present=.FALSE.
      
! add overshield i.e. the twinax shield
    PUL%overshield_present=.TRUE.
    PUL%overshield_x = 0d0          ! shield is centred at the origin in this calculation
    PUL%overshield_y = 0d0
    PUL%overshield_r = rs           ! twinax shield radius
    
    CALL PUL_LC_Laplace(PUL,cable%cable_name,cable%Y_fit_model_order,cable%Y_fit_freq_spec,domain) 

! there may be slight asymmmetry due to meshing so average diagonal and off diagonal elements  
! Theory_Manual_Eqn 3.21
    L11=(PUL%L%mat(1,1)+PUL%L%mat(1,1))/2d0
    L12=(PUL%L%mat(1,2)+PUL%L%mat(2,1))/2d0
    
    C11=(PUL%C%mat(1,1)+PUL%C%mat(1,1))/2d0
    C12=(PUL%C%mat(1,2)+PUL%C%mat(2,1))/2d0
    
    dielectric_is_homogeneous=.FALSE.
    
    CALL shielded_twisted_pair_cm_dm_parameter_calculation(L11,L12,C11,C12,epsr,LC,LD,CC,CD,dielectric_is_homogeneous)    
 ! Theory_Manual_Eqn 3.22
 
    YD=0.5d0*( PUL%Yfilter%sfilter_mat(1,1)+((-1d0)*PUL%Yfilter%sfilter_mat(1,2)) )
    YC=2.0d0*( PUL%Yfilter%sfilter_mat(1,1)+PUL%Yfilter%sfilter_mat(1,2) )          
    
  else

! See C.R. Paul, 1st edition, equation 3.67a,b with cos(thetaij)=-1 ! Theory_Manual_Eqn 2.27, 2.28

    L11=(mu0/(2d0*pi))*log( (rs**2-(s/2d0)**2)/(rs*rw) )
    L12=(mu0/(2d0*pi))*log( (s/(2d0*rs)) * (rs**2+(s/2d0)**2)/(2d0*(s/2d0)**2) )
    
    dielectric_is_homogeneous=.TRUE.
    
    CALL shielded_twisted_pair_cm_dm_parameter_calculation(L11,L12,C11,C12,epsr,LC,LD,CC,CD,dielectric_is_homogeneous)    

    YD=CD*jw
    YC=CC*jw
    
  end if  
  
  if (use_laplace) CALL deallocate_PUL_data(PUL)  ! deallocate the PUL data structure
  
  domain=0
  
  do inner_cable=1,4
    
! DOMAIN 1 of this sub-cable: Set the parameters for the internal differential mode domain  
    domain=domain+1
    
    cable%n_internal_conductors_in_domain(domain)=2
    
! The number of modes in the internal differential mode domain is 1

    dim=1
    cable%L_domain(domain)%dim=dim
    ALLOCATE(cable%L_domain(domain)%mat(dim,dim))
    cable%C_domain(domain)%dim=dim
    ALLOCATE(cable%C_domain(domain)%mat(dim,dim))
  
    cable%Z_domain(domain)%dim=dim
    ALLOCATE(cable%Z_domain(domain)%sfilter_mat(dim,dim))
    cable%Y_domain(domain)%dim=dim
    ALLOCATE(cable%Y_domain(domain)%sfilter_mat(dim,dim))
 
    cable%L_domain(domain)%mat(1,1)=LD 
    cable%Z_domain(domain)%sfilter_mat(1,1)=cable%L_domain(domain)%mat(1,1)*jw
  
    cable%C_domain(domain)%mat(1,1)=CD
    cable%Y_domain(domain)%sfilter_mat(1,1)=YD
           
! DOMAIN 2 of this sub-cable: Set the parameters for the internal common mode domain  
    domain=domain+1
    
    cable%n_internal_conductors_in_domain(domain)=2
    
! The number of modes in the internal common mode domain is 2

    dim=1
    cable%L_domain(domain)%dim=dim
    ALLOCATE(cable%L_domain(domain)%mat(dim,dim))
    cable%C_domain(domain)%dim=dim
    ALLOCATE(cable%C_domain(domain)%mat(dim,dim))
  
    cable%Z_domain(domain)%dim=dim
    ALLOCATE(cable%Z_domain(domain)%sfilter_mat(dim,dim))
    cable%Y_domain(domain)%dim=dim
    ALLOCATE(cable%Y_domain(domain)%sfilter_mat(dim,dim))

    cable%L_domain(domain)%mat(1,1)=LC 
    cable%Z_domain(domain)%sfilter_mat(1,1)=cable%L_domain(domain)%mat(1,1)*jw
  
    cable%C_domain(domain)%mat(1,1)=CC
    cable%Y_domain(domain)%sfilter_mat(1,1)=YC
  
  end do ! next inner cable 
   
! now set the parameters for the domain consisting of the outer shield and shields of the inner cables
  domain=domain+1
      
  cable%n_internal_conductors_in_domain(domain)=5
! The number of modes in this domain is 4
  dim=4
  cable%L_domain(domain)%dim=dim
  ALLOCATE(cable%L_domain(domain)%mat(dim,dim))
  cable%C_domain(domain)%dim=dim
  ALLOCATE(cable%C_domain(domain)%mat(dim,dim))
  
  cable%Z_domain(domain)%dim=dim
  ALLOCATE(cable%Z_domain(domain)%sfilter_mat(dim,dim))
  cable%Y_domain(domain)%dim=dim
  ALLOCATE(cable%Y_domain(domain)%sfilter_mat(dim,dim))

! get the parameters for this domain

  d=stpr       ! radius on which the centres of the 4 shielded twisted pairs sit
  epsr=evaluate_Sfilter_high_frequency_limit(epsr2)    ! 

  if (use_laplace) then
  
! allocate memory for the PUL parameter solver interface

    if(verbose) write(*,*)'Domain:',domain
    if(verbose) write(*,*)'Allocating PUL data structure for shielded twisted pairs'
    n_conductors=5
    
    CALL allocate_and_reset_PUL_data(PUL,n_conductors)

    PUL%shape(1:n_conductors)=circle
    
    PUL%x(1)=stpr
    PUL%y(1)=0.0
    PUL%r(1)=rs
    PUL%rd(1)=rd2
    PUL%epsr(1)=epsr2 ! permittivity of frequency dependent dielectric medium surrounding conductors
    
    PUL%x(2)=0.0
    PUL%y(2)=stpr
    PUL%r(2)=rs
    PUL%rd(2)=rd2
    PUL%epsr(2)=epsr2 ! permittivity of frequency dependent dielectric medium surrounding conductors
    
    PUL%x(3)=-stpr
    PUL%y(3)=0.0
    PUL%r(3)=rs
    PUL%rd(3)=rd2
    PUL%epsr(3)=epsr2 ! permittivity of frequency dependent dielectric medium surrounding conductors
    
    PUL%x(4)=0.0
    PUL%y(4)=-stpr
    PUL%r(4)=rs
    PUL%rd(4)=rd2
    PUL%epsr(4)=epsr2 ! permittivity of frequency dependent dielectric medium surrounding conductors

    PUL%epsr_background =1d0 ! permittivity of homogeneous medium surrounding the insulated shields (air)

! no ground plane
    PUL%ground_plane_present=.FALSE.
      
! add overshield i.e. the twinax shield
    PUL%overshield_present=.TRUE.
    PUL%overshield_x = 0d0          ! shield is centred at the origin in this calculation
    PUL%overshield_y = 0d0
    PUL%overshield_r = rs2           ! twisted pair shield radius
    
    CALL PUL_LC_Laplace(PUL,cable%cable_name,cable%Y_fit_model_order,cable%Y_fit_freq_spec,domain) 
       
    cable%L_domain(domain)%mat(:,:)=PUL%L%mat(:,:)
    cable%C_domain(domain)%mat(:,:)=PUL%C%mat(:,:)
    
    cable%Z_domain(domain)%sfilter_mat(:,:)=PUL%Zfilter%sfilter_mat(:,:)
    cable%Y_domain(domain)%sfilter_mat(:,:)=PUL%Yfilter%sfilter_mat(:,:)
    
  else

! See C.R. Paul, 1st edition, equation 3.67a,b ! Theory_Manual_Eqn 2.27, 2.28
! self inductance
    L11=(mu0/(2d0*pi))*log( (rs2**2-d**2)/(rs2*rs) )
  
! adjacent conductors, theta=90degrees, cos(theta)=0
    ctheta=0
    L12=(mu0/(2d0*pi))*log( (d/(rs2)) * &
         sqrt(  ((d*d)**2+rs2**4-2d0*d*d*rs2*rs2*ctheta)/((d*d)**2+d**4-2d0*(d**4)*ctheta)  ) )
  
! opposite conductors, theta=180 degrees, cos(theta)=-1
    ctheta=-1d0
    L13=(mu0/(2d0*pi))*log( (d/(rs2)) *  &
        sqrt(  ((d*d)**2+rs2**4-2d0*d*d*rs2*rs2*ctheta)/((d*d)**2+d**4-2d0*(d**4)*ctheta)  ) )

    cable%L_domain(domain)%mat(1,1)=L11
    cable%L_domain(domain)%mat(1,2)=L12
    cable%L_domain(domain)%mat(1,3)=L13
    cable%L_domain(domain)%mat(1,4)=L12
  
    cable%L_domain(domain)%mat(2,1)=L12
    cable%L_domain(domain)%mat(2,2)=L11
    cable%L_domain(domain)%mat(2,3)=L12
    cable%L_domain(domain)%mat(2,4)=L13
  
    cable%L_domain(domain)%mat(3,1)=L13
    cable%L_domain(domain)%mat(3,2)=L12
    cable%L_domain(domain)%mat(3,3)=L11
    cable%L_domain(domain)%mat(3,4)=L12
  
    cable%L_domain(domain)%mat(4,1)=L12
    cable%L_domain(domain)%mat(4,2)=L13
    cable%L_domain(domain)%mat(4,3)=L12
    cable%L_domain(domain)%mat(4,4)=L11
    
! calculate the capacitance matrix from the inverse of the inductance matrix *eps0*epsr*mu0

    ierr=0   ! set ierr=0 on input to matrix inverse to cause the program to stop if we have a singular matrix
    CALL  dinvert_Gauss_Jordan(cable%L_domain(domain)%mat,4,cable%C_domain(domain)%mat,4,ierr)
    cable%C_domain(domain)%mat(:,:)=eps0*epsr*mu0*cable%C_domain(domain)%mat(:,:)  
   
    CALL Z_Y_from_L_C(cable%L_domain(domain),cable%C_domain(domain),cable%Z_domain(domain),cable%Y_domain(domain))
      
  end if    

  if (use_laplace) CALL deallocate_PUL_data(PUL)  ! deallocate the PUL data structure
  
! Set the domain decomposition matrices ! Theory_Manual_Eqn 6.13, 6.14

! The dimension of the domain transformation matrices is 14
  dim=14
  cable%MI%dim=dim
  ALLOCATE(cable%MI%mat(dim,dim))
  cable%MV%dim=dim
  ALLOCATE(cable%MV%mat(dim,dim))
  
  cable%MI%mat(1:dim,1:dim)=0d0
  cable%MV%mat(1:dim,1:dim)=0d0

! domain decomposition for the 4 sheilded twisted pair cables  
  do inner_cable=1,4
  
    conductor_1=2*(inner_cable-1)+1
    conductor_2=conductor_1+1
    
    diff=2*(inner_cable-1)+1
    com=diff+1
    inner_shield=8+inner_cable
    
    cable%MI%mat(diff,conductor_1)=0.5d0
    cable%MI%mat(diff,conductor_2)=-0.5d0
    
    cable%MI%mat(com,conductor_1)=1d0
    cable%MI%mat(com,conductor_2)=1d0
    
    cable%MV%mat(diff,conductor_1)=1d0
    cable%MV%mat(diff,conductor_2)=-1d0
    
    cable%MV%mat(com,conductor_1)=0.5d0
    cable%MV%mat(com,conductor_2)=0.5d0
    cable%MV%mat(com,inner_shield)=-1d0
    
  end do
  
! domain decomposition for the domain between outer shield and twisted pair shields

  outer_shield=13

  do inner_cable=1,4
  
    conductor_1=2*(inner_cable-1)+1
    conductor_2=conductor_1+1
    inner_shield=8+inner_cable
    
    cable%MI%mat(inner_shield,conductor_1)=1d0
    cable%MI%mat(inner_shield,conductor_2)=1d0
    cable%MI%mat(inner_shield,inner_shield)=1d0
    
    cable%MV%mat(inner_shield,inner_shield)=1d0
    cable%MV%mat(inner_shield,outer_shield)=-1d0
    
  end do
  
! domain decomposition for the external domain conductor.

  cable%MI%mat(13,1:13)=1d0
  
  cable%MV%mat(13,13)=1d0
  cable%MV%mat(13,14)=-1d0

! domain decomposition for the reference domain conductor. Outer shield current is equal to the sum of all other currents
  cable%MI%mat(14,1:14)=1d0
  
  cable%MV%mat(14,14)=1d0

! Set the local reference conductor numbering  
  ALLOCATE( cable%local_reference_conductor(13) )
  cable%local_reference_conductor(1)=2      ! differential mode, reference is the second conductor
  cable%local_reference_conductor(2)=9      ! common mode, reference is the local shield conductor
  cable%local_reference_conductor(3)=3      ! differential mode, reference is the second conductor
  cable%local_reference_conductor(4)=10     ! common mode, reference is the local shield conductor
  cable%local_reference_conductor(5)=4      ! differential mode, reference is the second conductor
  cable%local_reference_conductor(6)=11     ! common mode, reference is the local shield conductor
  cable%local_reference_conductor(7)=5      ! differential mode, reference is the second conductor
  cable%local_reference_conductor(8)=12     ! common mode, reference is the local shield conductor

  cable%local_reference_conductor(9)=13     ! domain within outer shield, reference is the outer shield conductor
  cable%local_reference_conductor(10)=13    ! domain within outer shield, reference is the outer shield conductor
  cable%local_reference_conductor(11)=13    ! domain within outer shield, reference is the outer shield conductor
  cable%local_reference_conductor(12)=13    ! domain within outer shield, reference is the outer shield conductor
 
  cable%local_reference_conductor(13)=0      ! external domain conductor, reference not known

! Set the local domain information: include a reference conductor in the count
  ALLOCATE( cable%local_domain_n_conductors(1:cable%tot_n_domains) )
  cable%local_domain_n_conductors(1)=2              ! differential mode domain 
  cable%local_domain_n_conductors(2)=2              ! common mode: reference in external domain 
  cable%local_domain_n_conductors(3)=2              ! differential mode domain 
  cable%local_domain_n_conductors(4)=2              ! common mode: reference in external domain 
  cable%local_domain_n_conductors(5)=2              ! differential mode domain 
  cable%local_domain_n_conductors(6)=2              ! common mode: reference in external domain 
  cable%local_domain_n_conductors(7)=2              ! differential mode domain 
  cable%local_domain_n_conductors(8)=2              ! common mode: reference in external domain 
  cable%local_domain_n_conductors(9)=5              ! domain within outer shield
  cable%local_domain_n_conductors(10)=2             ! external domain 
  
! Set the external domain conductor and dielectric information
   
  ALLOCATE( cable%external_model(cable%n_external_conductors) )
  CALL reset_external_conductor_model(cable%external_model(1))
  cable%external_model(1)%conductor_type=circle
  cable%external_model(1)%conductor_radius=rs2
  cable%external_model(1)%dielectric_radius=rd3
  cable%external_model(1)%dielectric_epsr=epsr3 
      
! set the conductor impedance model for the four pairs of inner conductors
  do i=1,7,2
  
    cable%conductor_impedance(i)%impedance_model_type=impedance_model_type_cylindrical_with_conductivity
    cable%conductor_impedance(i)%radius=rw
    cable%conductor_impedance(i)%conductivity=sigma_w
    cable%conductor_impedance(i)%Resistance_multiplication_factor=1.5d0
    
    cable%conductor_impedance(i+1)%impedance_model_type=impedance_model_type_cylindrical_with_conductivity
    cable%conductor_impedance(i+1)%radius=rw
    cable%conductor_impedance(i+1)%conductivity=sigma_w
    cable%conductor_impedance(i+1)%Resistance_multiplication_factor=0.5d0
    
  end do

! set the conductor impedance model for the four inner shields
  
  do i=9,12
    
    cable%conductor_impedance(i)%impedance_model_type=impedance_model_type_cylindrical_shield
    cable%conductor_impedance(i)%radius=rs
    cable%conductor_impedance(i)%thickness=t1
    cable%conductor_impedance(i)%conductivity=sigma_s1
    cable%conductor_impedance(i)%ZT_filter=ZT1
    
  end do
    
! set the transfer impedance model for the outer shield conductor
  
  cable%conductor_impedance(13)%impedance_model_type=impedance_model_type_cylindrical_shield
  cable%conductor_impedance(13)%radius=rs2
  cable%conductor_impedance(13)%thickness=t2
  cable%conductor_impedance(13)%conductivity=sigma_s2
  cable%conductor_impedance(13)%ZT_filter=ZT2
  
! Deallocate all filters   
  CALL deallocate_Sfilter(epsr1)
  CALL deallocate_Sfilter(epsr2)
  CALL deallocate_Sfilter(epsr3)
  CALL deallocate_Sfilter(ZT1)
  CALL deallocate_Sfilter(ZT2)
  CALL deallocate_Sfilter(jw)
  
  ALLOCATE( cable%conductor_label(1:cable%tot_n_conductors) )
  cable%conductor_label(1)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 1 : Twisted pair 1 wire 1'
  cable%conductor_label(2)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 2 : Twisted pair 1 wire 2'
  cable%conductor_label(3)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 3 : Twisted pair 2 wire 1'
  cable%conductor_label(4)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 4 : Twisted pair 2 wire 2'
  cable%conductor_label(5)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 5 : Twisted pair 3 wire 1'
  cable%conductor_label(6)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 6 : Twisted pair 3 wire 2'
  cable%conductor_label(7)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 7 : Twisted pair 4 wire 1'
  cable%conductor_label(8)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 8 : Twisted pair 4 wire 2'
  cable%conductor_label(9)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 9 : Inner Shield 1'
  cable%conductor_label(10)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 10: Inner Shield 2'
  cable%conductor_label(11)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 11: Inner Shield 3'
  cable%conductor_label(12)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 12: Inner Shield 4'
  cable%conductor_label(13)='Cable name: '//trim(cable%cable_name)//   &
                           '. type: '//trim(cable%cable_type_string)//'. conductor 13: Outer Shield'
  
END SUBROUTINE spacewire_set_internal_domain_information
!
! NAME
!     spacewire_plot
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     plot spacewire cable 
!
! COMMENTS
!     The angle has an impact here
!     The conductor geometry must be consistent with the documentation...
!
! HISTORY
!
!     started 14/4/2016 CJS 
!
!
SUBROUTINE spacewire_plot(cable,x_offset,y_offset,theta,xmin,xmax,ymin,ymax)

USE type_specifications
USE general_module
USE constants

IMPLICIT NONE

! variables passed to subroutine

  type(cable_specification_type),intent(IN)    :: cable
  
  real(dp),intent(IN) :: x_offset,y_offset,theta
  real(dp),intent(INOUT) ::  xmin,xmax,ymin,ymax

! local variables

  real(dp) :: x,y
  real(dp) :: rw,rd
  real(dp) :: s,rstp
  real(dp) :: xstp,ystp
  real(dp) :: isr         ! inner shield radius
  real(dp) :: idr         ! inner dielectric radius
  real(dp) :: osr         ! outer shield radius
  real(dp) :: odr         ! outer dielectric radius
  
  integer :: inner_cable

! START

  rw=cable%parameters(1)    ! inner conductor radius
  rd=cable%parameters(2)    ! inner conductor dielectric radius
  s=cable%parameters(3)     ! inner conductor separation
  isr=cable%parameters(4)   ! inner shield radius
  idr=cable%parameters(6)   ! inner dielectric radius
  rstp=cable%parameters(7)  ! inner conductor separation
  osr=cable%parameters(8)   ! outer shield radius
  odr=cable%parameters(10)   ! outer dielectric radius

  do inner_cable=1,4

! calculate the centre for this inner cable

    xstp=x_offset+rstp*sin((inner_cable-1)*pi/2d0-theta)
    ystp=y_offset+rstp*cos((inner_cable-1)*pi/2d0-theta)
  
! plot inner conductor, 1
    x=xstp+(s/2d0)*sin(-theta)
    y=ystp+(s/2d0)*cos(-theta)
    CALL write_circle(x,y,rw,conductor_geometry_file_unit,xmin,xmax,ymin,ymax)

! plot inner conductor, 2
    x=xstp-(s/2d0)*sin(-theta)
    y=ystp-(s/2d0)*cos(-theta)
    CALL write_circle(x,y,rw,conductor_geometry_file_unit,xmin,xmax,ymin,ymax)
  
! plot inner conductor dielectric, 1
    x=xstp+(s/2d0)*sin(-theta)
    y=ystp+(s/2d0)*cos(-theta)
    CALL write_circle(x,y,rd,dielectric_geometry_file_unit,xmin,xmax,ymin,ymax)

! plot inner conductor dielectric, 2
    x=xstp-(s/2d0)*sin(-theta)
    y=ystp-(s/2d0)*cos(-theta)
    CALL write_circle(x,y,rd,dielectric_geometry_file_unit,xmin,xmax,ymin,ymax)

! plot shield conductor
    x=xstp
    y=ystp
    CALL write_circle(x,y,isr,conductor_geometry_file_unit,xmin,xmax,ymin,ymax)

! plot circular dielectric
    x=xstp
    y=ystp

    CALL write_circle(x,y,idr,dielectric_geometry_file_unit,xmin,xmax,ymin,ymax)
  
  end do ! next inner cable
  

! plot outer shield conductor
  x=x_offset
  y=y_offset
  CALL write_circle(x,y,osr,conductor_geometry_file_unit,xmin,xmax,ymin,ymax)

! plot outer circular dielectric
  x=x_offset
  y=y_offset

  CALL write_circle(x,y,odr,dielectric_geometry_file_unit,xmin,xmax,ymin,ymax)
  
  RETURN
  
END SUBROUTINE spacewire_plot