code: def _eeg_spherepot_coil(rrs, coils, sphere):
      """Calculate the EEG in the sphere model."""
      rmags, cosmags, ws, bins = _triage_coils(coils)
      n_coils = bins[-1] + 1
      del coils
  
      # Shift to the sphere model coordinates
      rrs = rrs - sphere["r0"]
  
      B = np.zeros((3 * len(rrs), n_coils))
      for ri, rr in enumerate(rrs):
          # Only process dipoles inside the innermost sphere
          if np.sqrt(np.dot(rr, rr)) >= sphere["layers"][0]["rad"]:
              continue
          # fwd_eeg_spherepot_vec
          vval_one = np.zeros((len(rmags), 3))
  
          # Make a weighted sum over the equivalence parameters
          for eq in range(sphere["nfit"]):
              # Scale the dipole position
              rd = sphere["mu"][eq] * rr
              rd2 = np.sum(rd * rd)
              rd2_inv = 1.0 / rd2
              # Go over all electrodes
              this_pos = rmags - sphere["r0"]
  
              # Scale location onto the surface of the sphere (not used)
              # if sphere['scale_pos']:
              #     pos_len = (sphere['layers'][-1]['rad'] /
              #                np.sqrt(np.sum(this_pos * this_pos, axis=1)))
              #     this_pos *= pos_len
  
              # Vector from dipole to the field point
              a_vec = this_pos - rd
  
              # Compute the dot products needed
              a = np.sqrt(np.sum(a_vec * a_vec, axis=1))
              a3 = 2.0 / (a * a * a)
              r2 = np.sum(this_pos * this_pos, axis=1)
              r = np.sqrt(r2)
              rrd = np.sum(this_pos * rd, axis=1)
              ra = r2 - rrd
              rda = rrd - rd2
  
              # The main ingredients
              F = a * (r * a + ra)
              c1 = a3 * rda + 1.0 / a - 1.0 / r
              c2 = a3 + (a + r) / (r * F)
  
              # Mix them together and scale by lambda/(rd*rd)
              m1 = c1 - c2 * rrd
              m2 = c2 * rd2
  
              vval_one += (
                  sphere["lambda"][eq]
                  * rd2_inv
                  * (
                      m1[:, np.newaxis] * rd
                      + m2[:, np.newaxis] * this_pos
                  )
              )
  
              # compute total result
              xx = vval_one * ws[:, np.newaxis]
              zz = np.array(
                  [bincount(bins, x, bins[-1] + 1) for x in xx.T]
              )
              B[3 * ri : 3 * ri + 3, :] = zz
      # finishing by scaling by 1/(4*M_PI)
      B *= 0.25 / np.pi
      return B
  
imports: from copy import deepcopy
  import numpy
  from .._fiff.constants import FIFF
  from ..bem import _import_openmeeg, _make_openmeeg_geometry
  from ..fixes import bincount, jit
  from ..parallel import parallel_func
  from ..surface import _jit_cross, _project_onto_surface
  from ..transforms import apply_trans, invert_transform
  from ..utils import _check_option, _pl, fill_doc, logger, verbose, warn
label_desc: Use image node_id 17 for calling a global function with example usage: _eeg_spherepot_coil(rrs, coils, sphere) and returns: B