Move to stimuli.components: image_base, resolve_grating_params, draw_…

…regions

demo/illusions/grating_parameters.ipynb

@@ -31,7 +31,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "from stimuli.components.components import image_base\n",
+    "from stimuli.components import image_base\n",
     "\n",
     "visual_size = (8., 8.)\n",
     "ppd = 10\n",
@@ -56,7 +56,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "from stimuli.components.components import resolve_grating_params\n",
+    "from stimuli.components import resolve_grating_params\n",
     "help(resolve_grating_params)"
    ]
   },

stimuli/components/__init__.py

@@ -1,7 +1,71 @@
+import itertools
+import warnings
+
 import numpy as np
 
+from stimuli.utils import resolution
+
 from .components import *
-from .components import image_base
+
+
+def image_base(visual_size=None, shape=None, ppd=None, rotation=0.0):
+    """Create coordinate-arrays to serve as image base for drawing
+
+    Parameters
+    ----------
+    visual_size : Sequence[Number, Number], Number, or None (default)
+        visual size [height, width] of image, in degrees
+    ppd : Sequence[Number, Number], Number, or None (default)
+        pixels per degree [vertical, horizontal]
+    shape : Sequence[Number, Number], Number, or None (default)
+        shape [height, width] of image, in pixels
+    rotation : float, optional
+        rotation (in degrees) counterclockwise from 3 o'clock, by default 0.0
+
+    Returns
+    -------
+    dict[str, Any]
+        dict with keys:
+        "visual_size", "ppd" : resolved from input arguments,
+        "x", "y" : single axes
+        "xx", "yy": numpy.ndarray of shape,
+        "angular" : numpy.ndarray of shape, with angle (in rad) relative to center point
+                   at each pixel
+    """
+
+    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
+
+    # Image axes
+    x = np.linspace(-visual_size.width / 2, visual_size.width / 2, shape.width)
+    y = np.linspace(-visual_size.height / 2, visual_size.height / 2, shape.height)
+
+    # Linear distance image bases
+    xx, yy = np.meshgrid(x, y)
+
+    # City-block distance (frames)
+    cityblock = np.maximum(np.abs(xx), np.abs(yy))
+
+    # Radial distance
+    radial = np.sqrt(xx**2 + yy**2)
+
+    # Angular distance
+    angular = np.arctan2(xx, yy)
+    angular -= np.deg2rad(rotation + 90)
+    angular %= 2 * np.pi
+
+    return {
+        "visual_size": visual_size,
+        "ppd": ppd,
+        "shape": shape,
+        "rotation": rotation,
+        "x": x,
+        "y": y,
+        "horizontal": xx,
+        "vertical": yy,
+        "cityblock": cityblock,
+        "radial": radial,
+        "angular": angular,
+    }
 
 
 def mask_elements(
@@ -56,3 +120,188 @@ def mask_elements(
         "visual_size": base["visual_size"],
         "ppd": base["ppd"],
     }
+
+
+def resolve_grating_params(
+    length=None,
+    visual_angle=None,
+    ppd=None,
+    frequency=None,
+    n_phases=None,
+    phase_width=None,
+    period="ignore",
+):
+    """Resolve (if possible) spatial parameters for a grating
+
+    Gratings take the regular resolution parameters (shape, ppd, visual_size).
+    In addition, there has to be an additional specification
+    of the frequency, number and with of phases (e.g., bars).
+    This can be done in two ways:
+    a phase_width (in degrees) and n_phases,
+    and/or by specifying the spatial frequency of the grating (in cycles per degree)
+
+    The total shape (in pixels) and visual size (in degrees) has to match the
+    specification of the phases and their widths.
+    Thus, not all 6 parameters have to be specified, as long as the both the resolution
+    and the distribution of phases can be resolved.
+
+    Note: all phases in a grating have the same width
+
+    Parameters
+    ----------
+    length : Number, or None (default)
+        lenght of grating, in pixels
+    visual_angle :Number, or None (default)
+        visual angle of grating, in degrees
+    ppd : Number, or None (default)
+        pixels per degree along the axis of grating
+    frequency : Number, or None (default)
+        frequency of grating, in cycles per degree
+    n_phases : int, or None (default)
+        number of phases (e.g., bars), i.e., half the number of full periods
+    phase_width : Number, or None (default)
+        extend of a single phase (e.g., bar), in degrees
+    period : "full", "half", "ignore" (default)
+        whether to ensure the grating only has "full" periods,
+        half "periods", or no guarantees ("ignore")
+
+    Returns
+    -------
+    dict[str, Any]
+        dictionary with all six resolution & size parameters resolved.
+    """
+
+    if period not in ["ignore", "full", "half"]:
+        raise TypeError(f"period not understood: {period}")
+
+    # Try to resolve resolution
+    try:
+        length, visual_angle, ppd = resolution.resolve_1D(
+            length=length, visual_angle=visual_angle, ppd=ppd
+        )
+    except ValueError:
+        ppd = ppd
+        length = length
+        visual_angle = visual_angle
+
+    # Try to resolve number and width(s) of phases:
+    # Logic here is that phase_width expresses "degrees per phase",
+    # which we can convert to "phases_per_degree"
+    # phase_width = degrees_per_phase = 1 / phases_per_degree = 1 / (2*frequency)
+    if phase_width is not None:
+        phases_pd = 1 / phase_width
+        if frequency is not None and phases_pd != 2 * frequency:
+            raise ValueError(f"phase_width {phase_width} and frequency {frequency} don't match")
+    elif frequency is not None:
+        phases_pd = 2 * frequency
+    else:  # both are None:
+        phases_pd = None
+
+    # Logic here is that phase_width expresses "degrees per phase",
+    # which we can invert to phases_per_degree, analogous to ppd:
+    # n_phases = phases_per_degree * n_degrees
+    # is analogous to
+    # pix = ppd * n_degrees
+    # Thus we can resolve the number and spacing of phases also as a resolution
+    try:
+        n_phases, min_angle, phases_pd = resolution.resolve_1D(
+            length=n_phases,
+            visual_angle=visual_angle,
+            ppd=phases_pd,
+            round=False,
+        )
+    except Exception as e:
+        raise Exception("Could not resolve grating frequency, phase_width, n_phases") from e
+
+    # Ensure full/half period?
+    if period != "ignore":
+        # Round n_phases
+        if period == "full":  # n_phases has to be even
+            n_phases = np.round(n_phases / 2) * 2
+        elif period == "half":  # n_phases can be odd
+            n_phases = np.round(n_phases)
+
+        # Check if n_phases fit in length
+        if length is not None and n_phases > 0 and length % n_phases:
+            raise resolution.ResolutionError(f"Cannot fit {n_phases} phases in {length} pix")
+
+        # Recalculate phases_pd
+        n_phases, min_angle, phases_pd = resolution.resolve_1D(
+            length=n_phases,
+            visual_angle=visual_angle,
+            ppd=None,
+            round=False,
+        )
+
+    # Convert to frequency
+    old_phase_width = phase_width
+    old_frequency = frequency
+    phase_width = 1 / phases_pd
+    frequency = phases_pd / 2
+
+    if (old_phase_width is not None and phase_width != old_phase_width) or (
+        old_frequency is not None and frequency != old_frequency
+    ):
+        warnings.warn(
+            f"Adjusted frequency and phase width to ensure {period} period: {old_frequency} ->"
+            f" {frequency}, {old_phase_width} -> {phase_width}"
+        )
+
+    # Now resolve resolution
+    visual_angle = min_angle if visual_angle is None else visual_angle
+    length, visual_angle, ppd = resolution.resolve_1D(
+        length=length, visual_angle=visual_angle, ppd=ppd
+    )
+
+    # Check that frequency does not exceed Nyquist limit:
+    if frequency > (ppd / 2):
+        raise ValueError(
+            f"Grating frequency ({frequency}) should not exceed Nyquist limit {ppd/2} (ppd/2)"
+        )
+
+    # Accumulate edges of phases
+    edges = [*itertools.accumulate(itertools.repeat(phase_width, int(n_phases)))]
+    if "period" == "ignore":
+        edges += [visual_angle]
+
+    return {
+        "length": length,
+        "visual_angle": visual_angle,
+        "ppd": ppd,
+        "frequency": frequency,
+        "phase_width": phase_width,
+        "n_phases": int(n_phases),
+        "edges": edges,
+        "period": period,
+    }
+
+
+def draw_regions(mask, intensities, intensity_background=0.5):
+    """Draw image with intensities for components in mask
+
+    Parameters
+    ----------
+    mask : numpy.ndarray
+        image-array with integer-indices for each region to draw
+    intensities : Sequence[float, ...]
+        intensity value for each masked region.
+        Can specify as many intensities as number of masked regions;
+        If fewer intensities are passed than masked regions, cycles through intensities
+    intensity_background : float, optional
+        intensity value of background, by default 0.5
+
+    Returns
+    -------
+    numpy.ndarray
+        image-array, same shape as mask, with intensity assigned to each masked region
+    """
+
+    # Create background
+    img = np.ones(mask.shape) * intensity_background
+
+    # Assign intensities to masked regions
+    ints = [*itertools.islice(itertools.cycle(intensities), len(np.unique(mask)))]
+    for frame_idx, intensity in zip(np.unique(mask), ints):
+        img = np.where(mask == frame_idx, intensity, img)
+
+    return img

stimuli/components/angular.py

@@ -1,8 +1,7 @@
 import numpy as np
 
-from stimuli.components import mask_elements
+from stimuli.components import draw_regions, mask_elements, resolve_grating_params
 from stimuli.components.circular import ring
-from stimuli.components.components import draw_regions, resolve_grating_params
 from stimuli.utils import resolution
 
 __all__ = [

stimuli/components/circular.py

@@ -3,8 +3,7 @@
 
 import numpy as np
 
-from stimuli.components import mask_elements
-from stimuli.components.components import image_base, resolve_grating_params
+from stimuli.components import image_base, mask_elements, resolve_grating_params
 from stimuli.utils import resize_array, resolution
 
 __all__ = [