docs: def load_fmm_model(filename, radius=np.pi/3.0, sphere_radius=100.0):
      '''
      load_fmm_model(filename) yields the fmm model indicated by the given file name. Fmm models are
      triangle meshes that define a field value at every vertex as well as the parameters of a
      projection to the 2D cortical surface via a registartion name and projection parameters.
      The following options may be given:
        * radius (default: pi/3) specifies the radius that should be assumed by the model; fmm models
          do not specify radius values in their projection parameters, but pi/3 is usually large for
          localized models.
        * sphere_radius (default: 100) specifies the radius of the sphere that should be assumed by
          the model. Note that in Freesurfer, spheres have a radius of 100.
      '''
      if not os.path.exists(filename):
          models_path = os.path.join(library_path(), 'models')
          # we look for it in a number of ways:
          fname = next((fnm
                        for beg in ['', models_path] for end in ['', '.fmm', '.fmm.gz']
                        for fnm0 in [filename + end]
                        for fnm in [fnm0 if beg is None else os.path.join(beg, fnm0)]
                        if os.path.exists(fnm)),
                       None)
          if fname is None:
              raise ValueError('Given model/file name does not exist: %s' % filename)
          filename = fname
      if not os.path.isfile(filename):
          raise ValueError('Given filename (%s) is not a file!' % filename)
      gz = True if len(filename) > 3 and filename[-3:] == '.gz' else False
      lines = None
      with (gzip.open(filename, 'rt') if gz else open(filename, 'rt')) as f:
          lines = f.read().split('\n')
      if len(lines) < 3 or lines[0] != 'Flat Mesh Model Version: 1.0':
          raise ValueError('Given file does not contain to a valid flat mesh model!')
      n = int(lines[1].split(':')[1].strip())
      m = int(lines[2].split(':')[1].strip())
      reg = lines[3].split(':')[1].strip()
      hemi = lines[4].split(':')[1].strip().upper()
      center = list(map(float, lines[5].split(':')[1].strip().split(',')))
      onxaxis = list(map(float, lines[6].split(':')[1].strip().split(',')))
      method = lines[7].split(':')[1].strip().lower()
      tx = np.asarray(
          [list(map(float, row.split(',')))
           for row in lines[8].split(':')[1].strip(' \t[]').split(';')])
      if lines[9].startswith('AreaNames: ['):
          # we load the area names
          s = lines[9][12:-1]
          area_names = tuple(s.split(' '))
          l0 = 10        
      else:
          area_names = None
          l0 = 9
      crds = []
      for row in lines[l0:(n+l0)]:
          (left,right) = row.split(' :: ')
          crds.append(list(map(float, left.split(','))))
      crds = np.asarray([list(map(float, left.split(',')))
                         for row in lines[l0:(n+l0)]
                         for (left,right) in [row.split(' :: ')]])
      vals = np.asarray([list(map(float, right.split(',')))
                         for row in lines[l0:(n+l0)]
                         for (left,right) in [row.split(' :: ')]])
      tris = -1 + np.asarray(
          [list(map(int, row.split(',')))
           for row in lines[(n+l0):(n+m+l0)]])
      return RegisteredRetinotopyModel(
          RetinotopyMeshModel(tris, crds,
                              90-180/np.pi*vals[:,0], vals[:,1], np.asarray(vals[:,2], dtype=np.int),
                              transform=tx,
                              area_name_to_id=area_names),
          geo.MapProjection(registration=reg,
                            center=center,
                            center_right=onxaxis,
                            method=method,
                            radius=radius,
                            sphere_radius=sphere_radius,
                            chirality=hemi))