code: def forecast_matrix(sh_order, d_par, d_perp, bvals):
      r"""Compute the FORECAST radial matrix"""
      n_c = int((sh_order + 1) * (sh_order + 2) / 2)
      M = np.zeros((bvals.shape[0], n_c))
      counter = 0
      for l in range(0, sh_order + 1, 2):
          for m in range(-l, l + 1):
              M[:, counter] = (
                  2
                  * np.pi
                  * np.exp(-bvals * d_perp)
                  * psi_l(l, bvals * (d_par - d_perp))
              )
              counter += 1
      return M
  
imports: from warnings import warn
  import numpy
  from dipy.reconst.cache import Cache
  from dipy.reconst.multi_voxel import multi_voxel_fit
  from dipy.reconst.csdeconv import csdeconv
  from dipy.reconst.shm import real_sh_descoteaux_from_index
  from scipy.special import gamma, hyp1f1
  from dipy.core.geometry import cart2sphere
  from dipy.data import default_sphere
  from dipy.reconst.odf import OdfModel, OdfFit
  from scipy.optimize import leastsq
  from dipy.utils.optpkg import optional_package
label_desc: Use image node_id 4 for calling a global function with example usage: forecast_matrix(sh_order, d_par, d_perp, bvals) and returns: M