recursive SUBROUTINE em_field_custom(ele, param, s_rel, orb, local_ref, field, calcd, err_flag, &
                                         calc_potential, use_overlap, grid_allow_s_out_of_bounds, rf_time, used_eles)
  USE bmad !, except_dummy => em_field_custom
  use parameters_bmad
!  USE bmad_interface, except_dummy => em_field_custom
  use em_field_interface
  use muon_mod

  IMPLICIT NONE

  type(ele_struct), target :: ele
  type(lat_param_struct) param
  type(coord_struct), intent(in) :: orb
  type(em_field_struct) field, field_error
  type (ele_struct), pointer :: lord
  type (ele_pointer_struct), allocatable, optional :: used_eles(:)
  real(rp), intent(in) :: s_rel
  real(rp) t_rel
  real(rp) scale
 real(rp), optional :: rf_time
  logical, optional :: calcd, err_flag, calc_potential, grid_allow_s_out_of_bounds, use_overlap
  logical  out_of_range,local_ref
  logical grid, mpoles, err
  logical :: local_ref_frame = .false.
  logical :: first = .true.
  character(32) :: r_name = 'em_field_custom'
  integer i
  integer, save :: n

t_rel = orb%t 

mpoles = .false.
grid = .false.
if(ele%key == sbend$)then
 ! first check to see if the orbit is outside the boundaries of the map
  if(abs(orb%vec(1)) < 0.069 .or. index(ele%name,'QUAD') == 0)then  !inside the map
 ! then check to see if this element has a grid map
  if ((associated(ele%grid_field)))then
     grid = .true.
    elseif(ele%field_calc == refer_to_lords$)then
         do i=1,ele%n_lord
          lord =>pointer_to_lord(ele,i)
          if(associated(lord%grid_field))then
           grid = .true.
           exit
          endif
         end do
   endif
! next check for multipoles
   if(associated(ele%b_pole_elec) .or. associated(ele%a_pole_elec))mpoles = .true.
   if(grid .or. mpoles)then
    if(grid)ele%field_calc = fieldmap$
    if(mpoles)ele%field_calc = bmad_standard$ !multipole$

!  rf quad
    if(index(ele%name,'QUAD2') /= 0 .and. t_rel < rf_quad(2)%t_max) then
      scale = rf_quad(2)%amp_h*sin(rf_quad(2)%freq_h * t_rel + rf_quad(2)%phi_h) + rf_quad(2)%amp_v*sin(rf_quad(2)%freq_v * t_rel + rf_quad(2)%phi_v)
      if(index(ele%name,'SHORT')/=0) then
                 quad_params%short_quad_field_index(2)  = quad_params_0%short_quad_field_index(2)* (1.+ scale)
                 call set_ele_real_attribute (ele, 'FIELD_INDEX',  quad_params%short_quad_field_index(2), err)
                 call set_a_quad(ele,.false.)
  
      endif
    endif
!

    call em_field_calc(ele, param, s_rel, orb, local_ref_frame, field, calcd, err_flag, &
                                    calc_potential, use_overlap, grid_allow_s_out_of_bounds, rf_time, used_eles)
    ele%field_calc = custom$
   else
    call em_field_custom_ring(ele, param, s_rel, t_rel, orb, local_ref, field, calcd, err_flag, &
                                    calc_potential, use_overlap, grid_allow_s_out_of_bounds, rf_time, used_eles)
   endif
  else ! outside map boundaries
    field%B = [0.0_rp, ele%value(B_field$), 0.0_rp]
    field%E=0
 endif
 
   field_error%B=0
   field_error%E=0
   if(param%ixx == 1)call field_errors(ele, param,s_rel, orb, field_error)
   field%B = field%B + field_error%B
   field%E = field%E + field_error%E


!  write(98,'(a,2es12.4,5es12.4,3es12.4)')ele%name, ele%s_start, s_rel,ele%s_start+s_rel,ele%value(rho$)+orb%vec(1),field%B,field_error%B

 else
  call em_field_custom_inj(ele, param, s_rel,t_rel, orb, local_ref,field,calcd, err_flag, &
                                         calc_potential, use_overlap, grid_allow_s_out_of_bounds, rf_time, used_eles) 
endif

if(index(ele%name,'QUAD2')/= 0 .and. index(ele%name,'SHORT')/=0)then
 if(first)then
!    write(197,'(9a12)')'s_rel','t_rel','Bx','By,','Bz','Ex','Ey','Ez','ele%value(field_index$)'
    n=0
    first=.false.
 endif
     n=n+1
!     if((n/100)*100 == n) write(197,'(9es12.4)')s_rel,t_rel,field%B,field%E, ele%value(field_index$)
endif

return
end

recursive SUBROUTINE em_field_custom_ring(ele, param, s_rel, t_rel, orb, local_ref, field, calcd, err_flag, &
                                         calc_potential, use_overlap, grid_allow_s_out_of_bounds, rf_time, used_eles)
  USE bmad_struct
  USE bmad_interface
  use em_field_interface, dummy=> em_field_custom_ring

  IMPLICIT NONE


  type(ele_struct) :: ele
  type(lat_param_struct) param
  type(coord_struct) orb
  type(em_field_struct) field
  type (ele_pointer_struct), allocatable, optional :: used_eles(:)
  real(rp), intent(in) :: s_rel, t_rel
  real(rp) x,y,s,t
  real(rp), optional :: rf_time
  logical  out_of_range,local_ref
  logical, optional :: calcd, err_flag, calc_potential, grid_allow_s_out_of_bounds, use_overlap

  character(32) :: r_name = 'em_field_custom'

  x = orb%vec(1)
  y = orb%vec(3)
  s = s_rel
  t = orb%t  !t_rel            
  call fields(ele,x,y,s,t,field)

  RETURN
END SUBROUTINE em_field_custom_ring