code: def getionimage(p, mz_value, tol=0.1, z=1, reduce_func=sum):
      
      tol = abs(tol)
      im = np.zeros((p.imzmldict["max count of pixels y"], p.imzmldict["max count of pixels x"]))
      for i, (x, y, z_) in enumerate(p.coordinates):
          if z_ == 0:
              UserWarning("z coordinate = 0 present, if you're getting blank images set getionimage(.., .., z=0)")
          if z_ == z:
              mzs, ints = map(lambda x: np.asarray(x), p.getspectrum(i))
              min_i, max_i = _bisect_spectrum(mzs, mz_value, tol)
              im[y - 1, x - 1] = reduce_func(ints[min_i:max_i+1])
      return im
docs: Get an image representation of the intensity distribution
      of the ion with specified m/z value.
  
      By default, the intensity values within the tolerance region are summed.
  
      :param p:
          the ImzMLParser (or anything else with similar attributes) for the desired dataset
      :param mz_value:
          m/z value for which the ion image shall be returned
      :param tol:
          Absolute tolerance for the m/z value, such that all ions with values
          mz_value-|tol| <= x <= mz_value+|tol| are included. Defaults to 0.1
      :param z:
          z Value if spectrogram is 3-dimensional.
      :param reduce_func:
          the bahaviour for reducing the intensities between mz_value-|tol| and mz_value+|tol| to a single value. Must
          be a function that takes a sequence as input and outputs a number. By default, the values are summed.
  
      :return:
          numpy matrix with each element representing the ion intensity in this
          pixel. Can be easily plotted with matplotlib