shared_setup Subroutine

  subroutine shared_setup(self, tind_in, tind_out)
    use constants, only: twopi, pi
    class(peq_base_type), intent(in out) :: self
    integer, intent(in) :: tind_in, tind_out !< Theta index of inboard and outboard
    real :: R_geo, psi_N, f_n
    integer :: j
    self%psi_0 = self%eqpsi(1) ; self%psi_a = self%eqpsi(self%nr)
    self%psi_bar = (self%eqpsi - self%psi_0) / (self%psi_a - self%psi_0)
    self%aminor = 0.5 * abs(self%R_psi(self%nr, tind_out) - self%R_psi(self%nr, tind_in))
    R_geo = 0.5 * (self%R_psi(self%nr, tind_out) + self%R_psi(self%nr, tind_in))
    self%B_T = self%fp(self%nr) / R_geo ; self%B_psi = 0.0
    self%R_psi = self%R_psi / self%aminor ; self%Z_psi = self%Z_psi / self%aminor
    self%R_mag = self%R_psi(1, 1) ; self%Z_mag = self%Z_psi(1, 1)
    self%beta = 4. * twopi * 1.e-7 * self%pressure / self%B_T**2 ! not correct definition yet?
    self%beta_0 = self%beta(1)
    self%pressure = self%pressure / self%pressure(1)

    psi_N = abs(self%B_T) * self%aminor**2
    self%psi_a = self%psi_a / psi_N ; self%psi_0 = self%psi_0 / psi_N
    self%eqpsi = self%eqpsi / psi_N
    f_N = abs(self%B_T) * self%aminor
    self%fp = self%fp / f_N

    self%diam = abs(self%R_psi(:, tind_out) - self%R_psi(:, tind_in))
    self%rc = 0.5 * (self%R_psi(:, tind_out) + self%R_psi(:, tind_in))
    do j = 1, self%nt
       self%eqth(:, j) = (j - 1) * pi / float(self%nt - 1)
       self%eqpsi_2d(:, j) = self%eqpsi
    end do
    if (self%has_full_theta_range) self%eqth = 2 * self%eqth - pi
  end subroutine shared_setup