Repository: computational-psychology/stimupy (language: python)

Issue:
Enable `rotation` in `sbcs`?
At least in `generalized()`?
SBCs
Added some convenience `.sbcs.`-stimulus functions:
- `.square()`
- `.circular()`
These are based on the respective `bullseyes.`-functions, and thus ultimately o the `rings.`-functions. As such, they take `target_radius` and `surround_radius` arguments, similar to the bulleyes and rings, which allows for consistent specification of related stimuli.

Enabled `rotation` (#39) in:
- `sbcs.generalized`,
- `sbcs.square`,
- `sbcs.circular`,

Additionally, all the `sbcs.`-functions now have a corresponding `<>_two_sided()` variant.
These are now also used in several of the `papers.`-stimuli.
Unify defaults for related stimuli (e.g., `waves` and `whites`)
When stimuli are clearly related / derived, default argument values (where applicable) should be consistent across the different stimuli.

A glaring example are the `intensity_bars` of the `waves.square` (`(1,0)`) and the `whites` (`(0,1)`). These should all be aligned, across these stimuli, and all kinds of `waves` (including `gabors`).


Discussion / hints:


Should also apply to argument names, e.g. `intensities` in `waves.sine_linear` currently becomes `intensity_bars` in `gabor`.

Fix patch:
diff --git a/.github/workflows/version.yml b/.github/workflows/version.yml
index 6f0e15d6..39b2ef7c 100644
--- a/.github/workflows/version.yml
+++ b/.github/workflows/version.yml
@@ -10,6 +10,10 @@ jobs:
   bump-version:
     if: (github.event.pull_request.merged == true) || (github.event_name == 'workflow_dispatch')
     runs-on: ubuntu-latest
+    concurrency: bump-version
+    permissions:
+      id-token: write
+      contents: write
 
     steps:
       - name: Checkout repository 
@@ -18,20 +22,14 @@ jobs:
           fetch-depth: 0 # checkout full commit history
           token: ${{ secrets.GHA_Token }}
 
-      - name: Setup Python
-        uses: actions/setup-python@v4
+      - name: Python Semantic Release
+        id: versioning
+        uses: python-semantic-release/python-semantic-release@master
         with:
-          python-version: 3.11
+          github_token: $${{ secrets.GHA_TOKEN }}
 
-      - name: Configure git to be able to push to repo
-        env:
-          GH_TOKEN: ${{ secrets.GHA_Token }}
-        run: |
-          git config user.name token
-          git config user.email token@github.com
-
-      - name: Bump version, build & upload release assets using Semantic Release
-        env:
-          GH_TOKEN: ${{ secrets.GHA_Token }}
-        run: |
-          pipx run python-semantic-release publish -v DEBUG
\ No newline at end of file
+      - name: Publish package distributions to GitHub Releases
+        uses: python-semantic-release/upload-to-gh-release@main
+        if: steps.versioning.outputs.released == 'true'
+        with:
+          github_token: ${{ secrets.GHA_TOKEN }}
\ No newline at end of file
diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml
index cb14a90d..4121c6d8 100644
--- a/.pre-commit-config.yaml
+++ b/.pre-commit-config.yaml
@@ -1,19 +1,19 @@
 repos:
   - repo: https://github.com/asottile/pyupgrade
-    rev: v3.3.1
+    rev: v3.10.1
     hooks:
     -   id: pyupgrade
 
   - repo: https://github.com/psf/black
-    rev: stable
+    rev: 23.7.0
     hooks:
     - id: black
       additional_dependencies: [".[jupyter]"]
       language_version: python3
 
   - repo: https://github.com/PyCQA/flake8
-    rev: "6.0.0"
+    rev: "6.1.0"
     hooks:
       - id: flake8
         # Ignore all format-related checks as Black takes care of those.
-        args: ["--ignore", "E2,W5,E501,F403", "--select", "E,W,F", '--per-file-ignores=__init__.py:F401,F405,E402', '--exclude','tests/*']
\ No newline at end of file
+        args: ["--ignore", "E203,E501,W503,F403", "--select", "E,W,F", '--per-file-ignores=__init__.py:F401,F405,E402', '--exclude','tests/*']
\ No newline at end of file
diff --git a/docs/reference/demos/components/texts.md b/docs/reference/demos/components/texts.md
new file mode 100644
index 00000000..d05f4ad7
--- /dev/null
+++ b/docs/reference/demos/components/texts.md
@@ -0,0 +1,106 @@
+---
+jupytext:
+  formats: md:myst
+  text_representation:
+    extension: .md
+    format_name: myst
+    format_version: 0.13
+    jupytext_version: 1.14.5
+kernelspec:
+  display_name: Python 3 (ipykernel)
+  language: python
+  name: python3
+---
+
+```{tip}
+[![Binder](https://mybinder.org/badge_logo.svg)](https://mybinder.org/v2/gh/computational-psychology/stimupy/HEAD?urlpath=lab/tree/docs/reference/demos/components/texts.md)
+ to get interactivity
+```
+```{attention}
+To run locally, the code for these interactive demos requires
+a [Jupyter Notebook](https://jupyter.org/) environment,
+and the [Jupyter Widgets extension (`ipywidgets`)](https://ipywidgets.readthedocs.io/en/latest/index.html).
+```
+
+# Components - Texts
+{py:mod}`stimupy.components.texts`
+
+
+## Text
+{py:func}`stimupy.components.texts.text`
+
+```{code-cell} ipython3
+import ipywidgets as iw
+from stimupy.utils import plot_stim
+from stimupy.components.texts import text
+
+# Define widgets
+w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
+w_width = iw.IntSlider(value=10, min=1, max=20, description="width [deg]")
+w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
+
+w_int = iw.FloatSlider(value=1., min=0, max=1, description="intensity line")
+w_int_back = iw.FloatSlider(value=0., min=0, max=1, description="intensity background")
+
+w_text = iw.Textarea(value='Hello World', placeholder='Type something', description='text')
+w_fontsize = iw.IntSlider(value=30, min=1, max=60, description="fontsize")
+w_align = iw.Dropdown(value="center", options=['left', 'center', 'right'], description="align")
+w_direction = iw.Dropdown(value="ltr", options=['ltr', 'rtl'], description="direction")
+
+w_mask = iw.ToggleButton(value=False, disabled=False, description="add mask")
+
+# Layout
+b_im_size = iw.HBox([w_height, w_width, w_ppd])
+b_text = iw.HBox([w_text, w_fontsize, w_align, w_direction])
+b_intensities = iw.HBox([w_int, w_int_back])
+b_add = iw.HBox([w_mask])
+ui = iw.VBox([b_im_size, b_text, b_intensities, b_add])
+
+# Function for showing stim
+def show_text(
+    height=None,
+    width=None,
+    ppd=None,
+    str_text=None,
+    intensity_text=None,
+    intensity_background=None,
+    fontsize=None,
+    align=None,
+    direction=None,
+    add_mask=False,
+):
+    try:
+        stim = text(
+            text=str_text,
+            visual_size=(height, width),
+            ppd=ppd,
+            intensity_text=intensity_text,
+            intensity_background=intensity_background,
+            fontsize=fontsize,
+            align=align,
+            direction=direction,
+        )
+        plot_stim(stim, mask=add_mask)
+    except Exception as e:
+        raise ValueError(f"Invalid parameter combination: {e}") from None    
+
+# Set interactivity
+out = iw.interactive_output(
+    show_text,
+    {
+        "height": w_height,
+        "width": w_width,
+        "ppd": w_ppd,
+        "str_text": w_text,
+        "intensity_text": w_int,
+        "intensity_background": w_int_back,
+        "fontsize": w_fontsize,
+        "align": w_align,
+        "direction": w_direction,
+        "add_mask": w_mask,
+    },
+)
+
+# Show
+display(ui, out)
+```
\ No newline at end of file
diff --git a/docs/reference/demos/stimuli/bullseyes.md b/docs/reference/demos/stimuli/bullseyes.md
index 58dfa2a2..3d76e7cd 100644
--- a/docs/reference/demos/stimuli/bullseyes.md
+++ b/docs/reference/demos/stimuli/bullseyes.md
@@ -115,6 +115,99 @@ out = iw.interactive_output(
 display(ui, out)
 ```
 
+## Circular bullseye, generalized
+{py:func}`stimupy.stimuli.bullseyes.circular_generalized`
+
+```{code-cell} ipython3
+import ipywidgets as iw
+from stimupy.utils import plot_stim
+from stimupy.stimuli.bullseyes import circular_generalized
+
+# Define widgets
+w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
+w_width = iw.IntSlider(value=10, min=1, max=20, description="width [deg]")
+w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
+
+w_radius1 = iw.FloatSlider(value=1, min=0, max=2, description="radius1 [deg]")
+w_radius2 = iw.FloatSlider(value=2, min=1, max=3, description="radius2 [deg]")
+w_radius3 = iw.FloatSlider(value=3, min=2, max=4, description="radius2 [deg]")
+
+w_int1 = iw.FloatSlider(value=1, min=0, max=1, description="int-ring1")
+w_int2 = iw.FloatSlider(value=0.3, min=0, max=1, description="int-ring2")
+w_int3 = iw.FloatSlider(value=0.8, min=0, max=1, description="int-ring3")
+w_int_back = iw.FloatSlider(value=0., min=0, max=1, description="intensity background")
+
+w_ori = iw.Dropdown(value="center", options=['mean', 'corner', 'center'], description="origin")
+w_rot = iw.FloatSlider(value=0, min=0, max=360, description="rotation [deg]")
+w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'ring_mask'], description="add mask")
+
+w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+
+# Layout
+b_im_size = iw.HBox([w_height, w_width, w_ppd])
+b_geometry = iw.HBox([w_radius1, w_radius2, w_radius3])
+b_intensities = iw.HBox([w_int1, w_int2, w_int3, w_int_back])
+b_target = iw.HBox([w_tint])
+b_add = iw.HBox([w_ori, w_rot, w_mask])
+ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
+
+# Function for showing stim
+def show_circular_generalized(
+    height=None,
+    width=None,
+    ppd=None,
+    radius1=None,
+    radius2=None,
+    radius3=None,
+    int1=None,
+    int2=None,
+    int3=None,
+    intensity_background=None,
+    origin=None,
+    add_mask=False,
+    intensity_target=None,
+    rotation=0.0,
+):
+    try:
+        stim = circular_generalized(
+            visual_size=(height, width),
+            ppd=ppd,
+            radii=(radius1, radius2, radius3),
+            intensity_rings=(int1, int2, int3),
+            intensity_background=intensity_background,
+            origin=origin,
+            intensity_target=intensity_target,
+        )
+        plot_stim(stim, mask=add_mask)
+    except Exception as e:
+        raise ValueError(f"Invalid parameter combination: {e}") from None
+
+# Set interactivity
+out = iw.interactive_output(
+    show_circular_generalized,
+    {
+        "height": w_height,
+        "width": w_width,
+        "ppd": w_ppd,
+        "radius1": w_radius1,
+        "radius2": w_radius2,
+        "radius3": w_radius3,
+        "int1": w_int1,
+        "int2": w_int2,
+        "int3": w_int3,
+        "intensity_background": w_int_back,
+        "origin": w_ori,
+        "add_mask": w_mask,
+        "intensity_target": w_tint,
+        "rotation": w_rot,
+    },
+)
+
+# Show
+display(ui, out)
+```
+
+
 
 ## Circular, Two-sided
 {py:func}`stimupy.stimuli.bullseyes.circular_two_sided`
@@ -139,14 +232,15 @@ w_int_back = iw.FloatSlider(value=0.5, min=0, max=1, description="int background
 w_ori = iw.Dropdown(value="mean", options=['mean', 'corner', 'center'], description="origin")
 w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'ring_mask'], description="add mask")
 
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="left target int")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="right target int")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
 b_geometry = iw.HBox([w_freq, w_phase])
 b_intensities = iw.HBox([w_int1, w_int2, w_int_back])
 b_add = iw.HBox([w_ori, w_mask])
-b_target = iw.HBox([w_tint])
+b_target = iw.HBox([w_tint_l, w_tint_r])
 ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
 
 # Function for showing stim
@@ -161,7 +255,8 @@ def show_two_sided_bullseye(
     intensity_background=None,
     origin=None,
     add_mask=False,
-    intensity_target=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
 ):
     try:
         stim = circular_two_sided(
@@ -169,10 +264,10 @@ def show_two_sided_bullseye(
             ppd=ppd,
             frequency=frequency,
             phase_shift=phase_shift,
-            intensity_rings=(int1, int2),
+            intensity_rings=((int1, int2),(int2, int1)),
             intensity_background=intensity_background,
             origin=origin,
-            intensity_target=intensity_target,
+            intensity_target=(intensity_target_l, intensity_target_r),
         )
         plot_stim(stim, mask=add_mask)
     except Exception as e:
@@ -192,7 +287,8 @@ out = iw.interactive_output(
         "intensity_background": w_int_back,
         "origin": w_ori,
         "add_mask": w_mask,
-        "intensity_target": w_tint,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
     },
 )
 
@@ -409,13 +505,14 @@ w_rot = iw.FloatSlider(value=0, min=0, max=360, description="rotation [deg]")
 
 w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'frame_mask'], description="add mask")
 
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="left target int")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="right target int")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
 b_geometry = iw.HBox([w_freq, w_phase])
 b_intensities = iw.HBox([w_int1, w_int2, w_int_back])
-b_target = iw.HBox([w_tint])
+b_target = iw.HBox([w_tint_l, w_tint_r])
 b_add = iw.HBox([w_rot, w_mask])
 ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
 
@@ -430,7 +527,8 @@ def show_rectangular_two_sided(
     int2=None,
     intensity_background=None,
     add_mask=False,
-    intensity_target=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
     rotation=0.0,
 ):
     try:
@@ -439,9 +537,9 @@ def show_rectangular_two_sided(
             ppd=ppd,
             frequency=frequency,
             phase_shift=phase_shift,
-            intensity_frames=(int1, int2),
+            intensity_frames=((int1, int2),(int2, int1)),
             intensity_background=intensity_background,
-            intensity_target=intensity_target,
+            intensity_target=(intensity_target_l, intensity_target_r),
             rotation=rotation,
         )
         plot_stim(stim, mask=add_mask)
@@ -461,7 +559,8 @@ out = iw.interactive_output(
         "int2": w_int2,
         "intensity_background": w_int_back,
         "add_mask": w_mask,
-        "intensity_target": w_tint,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
         "rotation": w_rot,
     },
 )
diff --git a/docs/reference/demos/stimuli/delboeufs.md b/docs/reference/demos/stimuli/delboeufs.md
index 528e2760..aa821436 100644
--- a/docs/reference/demos/stimuli/delboeufs.md
+++ b/docs/reference/demos/stimuli/delboeufs.md
@@ -121,23 +121,24 @@ w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
 w_width = iw.IntSlider(value=20, min=1, max=40, description="width [deg]")
 w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
 
-w_outr1 = iw.FloatSlider(value=4, min=0.5, max=8, description="outer radius1 [deg]")
-w_outr2 = iw.FloatSlider(value=2.5, min=0.5, max=8, description="outer radius2 [deg]")
+w_out_rad_l = iw.FloatSlider(value=4, min=0.5, max=8, description="outer radius1 [deg]")
+w_out_rad_r = iw.FloatSlider(value=2.5, min=0.5, max=8, description="outer radius2 [deg]")
 w_outw = iw.FloatSlider(value=0, min=0, max=2, description="outer line width [deg]")
 
 w_int1 = iw.FloatSlider(value=0, min=0, max=1, description="int line")
 w_int_back = iw.FloatSlider(value=1., min=0, max=1, description="int background")
 
 w_tr = iw.FloatSlider(value=2., min=0, max=4, description="target radius [deg]")
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="int target")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="int left target")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="int right target")
 
 w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'line_mask'], description="add mask")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
-b_geometry = iw.HBox([w_outr1, w_outr2, w_outw])
+b_geometry = iw.HBox([w_out_rad_l, w_out_rad_r, w_outw])
 b_intensities = iw.HBox([w_int1, w_int_back])
-b_target = iw.HBox([w_tr, w_tint])
+b_target = iw.HBox([w_tr, w_tint_l, w_tint_r])
 b_add = iw.HBox([w_mask])
 ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
 
@@ -152,19 +153,20 @@ def show_two_sided(
     target_radius=None,
     intensity_outer_line=None,
     intensity_background=None,
-    intensity_target=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
     add_mask=False,
 ):
     try:
         stim = two_sided(
             visual_size=(height, width),
             ppd=ppd,
-            outer_radii=(outer_radius1, outer_radius2),
+            outer_radius=(outer_radius1, outer_radius2),
             outer_line_width=outer_line_width,
             target_radius=target_radius,
             intensity_outer_line=intensity_outer_line,
             intensity_background=intensity_background,
-            intensity_target=intensity_target,
+            intensity_target=(intensity_target_l, intensity_target_r),
         )
         plot_stim(stim, mask=add_mask)
     except Exception as e:
@@ -178,11 +180,12 @@ out = iw.interactive_output(
         "width": w_width,
         "ppd": w_ppd,
         "add_mask": w_mask,
-        "outer_radius1": w_outr1,
-        "outer_radius2": w_outr2,
+        "outer_radius_l": w_out_rad_l,
+        "outer_radius_r": w_out_rad_r,
         "outer_line_width": w_outw,
         "target_radius": w_tr,
-        "intensity_target": w_tint,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
         "intensity_background": w_int_back,
         "intensity_outer_line": w_int1,
     },
diff --git a/docs/reference/demos/stimuli/mueller_lyers.md b/docs/reference/demos/stimuli/mueller_lyers.md
index 260c2cb3..4f3ad9e6 100644
--- a/docs/reference/demos/stimuli/mueller_lyers.md
+++ b/docs/reference/demos/stimuli/mueller_lyers.md
@@ -126,22 +126,27 @@ w_width = iw.IntSlider(value=10, min=1, max=20, description="width [deg]")
 w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
 
 w_outl = iw.FloatSlider(value=1, min=0.1, max=2, description="o-line length [deg]")
-w_outa = iw.IntSlider(value=45, min=-180, max=180, description="o-line angle [deg]")
+w_outa_l = iw.IntSlider(value=45, min=-180, max=180, description="o-line angle left [deg]")
+w_outa_r = iw.IntSlider(value=45, min=-180, max=180, description="o-line angle right [deg]")
+
 w_w = iw.FloatSlider(value=0, min=0, max=0.5, description="o-line width [deg]")
 
 w_int1 = iw.FloatSlider(value=1, min=0, max=1, description="int line")
 w_int_back = iw.FloatSlider(value=0, min=0, max=1, description="int background")
 
-w_tl = iw.FloatSlider(value=2.5, min=0, max=5, description="target length [deg]")
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="int target")
+w_tl_l = iw.FloatSlider(value=2.5, min=0, max=5, description="target length left [deg]")
+w_tl_r = iw.FloatSlider(value=2.5, min=0, max=5, description="target length right [deg]")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="int target left")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="int target right")
+
 
 w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'line_mask'], description="add mask")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
-b_geometry = iw.HBox([w_outl, w_outa, w_w])
+b_geometry = iw.HBox([w_outl, w_outa_l, w_outa_r, w_w])
 b_intensities = iw.HBox([w_int1, w_int_back])
-b_target = iw.HBox([w_tl, w_tint])
+b_target = iw.HBox([w_tl_l, w_tl_r, w_tint_l, w_tint_r])
 b_add = iw.HBox([w_mask])
 ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
 
@@ -152,11 +157,14 @@ def show_two_sided(
     ppd=None,
     outer_lines_length=None,
     line_width=None,
-    outer_lines_angle=None,
-    target_length=None,
+    outer_lines_angle_l=None,
+    outer_lines_angle_r=None,
+    target_length_l=None,
+    target_length_r=None,
     intensity_outer_line=None,
     intensity_background=None,
-    intensity_target=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
     add_mask=False,
 ):
     try:
@@ -164,12 +172,12 @@ def show_two_sided(
             visual_size=(height, width),
             ppd=ppd,
             outer_lines_length=outer_lines_length,
-            outer_lines_angle=outer_lines_angle,
+            outer_lines_angle=(outer_lines_angle_l, outer_lines_angle_r),
             line_width=line_width,
-            target_length=target_length,
+            target_length=(target_length_l, target_length_r),
             intensity_outer_lines=intensity_outer_line,
             intensity_background=intensity_background,
-            intensity_target=intensity_target,
+            intensity_target=(intensity_target_l, intensity_target_r),
         )
         plot_stim(stim, mask=add_mask)
     except Exception as e:
@@ -185,9 +193,12 @@ out = iw.interactive_output(
         "add_mask": w_mask,
         "outer_lines_length": w_outl,
         "line_width": w_w,
-        "outer_lines_angle": w_outa,
-        "target_length": w_tl,
-        "intensity_target": w_tint,
+        "outer_lines_angle_l": w_outa_l,
+        "outer_lines_angle_r": w_outa_r,
+        "target_length_l": w_tl_l,
+        "target_length_r": w_tl_r,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
         "intensity_background": w_int_back,
         "intensity_outer_line": w_int1,
     },
diff --git a/docs/reference/demos/stimuli/pinwheels.md b/docs/reference/demos/stimuli/pinwheels.md
index b1fc4a75..cd296edd 100644
--- a/docs/reference/demos/stimuli/pinwheels.md
+++ b/docs/reference/demos/stimuli/pinwheels.md
@@ -53,7 +53,7 @@ w_tcent = iw.FloatSlider(value=2.5, min=0, max=5, description="target center")
 w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
 
 w_ori = iw.Dropdown(value="mean", options=['mean', 'corner', 'center'], description="origin")
-w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'wedge_mask'], description="add mask")
+w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'segment_mask', 'circle_mask'], description="add mask")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
diff --git a/docs/reference/demos/stimuli/rings.md b/docs/reference/demos/stimuli/rings.md
index 7defd9d6..4a8fa5a6 100644
--- a/docs/reference/demos/stimuli/rings.md
+++ b/docs/reference/demos/stimuli/rings.md
@@ -119,6 +119,99 @@ out = iw.interactive_output(
 display(ui, out)
 ```
 
+## Circular, generalized
+{py:func}`stimupy.stimuli.rings.circular_generalized`
+
+```{code-cell} ipython3
+import ipywidgets as iw
+from stimupy.utils import plot_stim
+from stimupy.stimuli.rings import circular_generalized
+
+# Define widgets
+w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
+w_width = iw.IntSlider(value=10, min=1, max=20, description="width [deg]")
+w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
+
+w_radius1 = iw.FloatSlider(value=1, min=0, max=2, description="radius1 [deg]")
+w_radius2 = iw.FloatSlider(value=2, min=1, max=3, description="radius2 [deg]")
+w_radius3 = iw.FloatSlider(value=3, min=2, max=4, description="radius2 [deg]")
+
+w_int1 = iw.FloatSlider(value=1, min=0, max=1, description="int-ring1")
+w_int2 = iw.FloatSlider(value=0.3, min=0, max=1, description="int-ring2")
+w_int3 = iw.FloatSlider(value=0.8, min=0, max=1, description="int-ring3")
+w_int_back = iw.FloatSlider(value=0., min=0, max=1, description="intensity background")
+
+w_ori = iw.Dropdown(value="center", options=['mean', 'corner', 'center'], description="origin")
+w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'ring_mask'], description="add mask")
+
+w_tidx = iw.IntSlider(value=1, min=1, max=4, description="target idx")
+w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+
+# Layout
+b_im_size = iw.HBox([w_height, w_width, w_ppd])
+b_geometry = iw.HBox([w_radius1, w_radius2, w_radius3])
+b_intensities = iw.HBox([w_int1, w_int2, w_int3, w_int_back])
+b_target = iw.HBox([w_tidx, w_tint])
+b_add = iw.HBox([w_ori, w_mask])
+ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
+
+# Function for showing stim
+def show_circular_generalized(
+    height=None,
+    width=None,
+    ppd=None,
+    radius1=None,
+    radius2=None,
+    radius3=None,
+    int1=None,
+    int2=None,
+    int3=None,
+    intensity_background=None,
+    origin=None,
+    add_mask=False,
+    intensity_target=None,
+    target_indices=None,
+):
+    try:
+        stim = circular_generalized(
+            visual_size=(height, width),
+            ppd=ppd,
+            radii=(radius1, radius2, radius3),
+            intensity_rings=(int1, int2, int3),
+            intensity_background=intensity_background,
+            origin=origin,
+            intensity_target=intensity_target,
+            target_indices=target_indices,
+        )
+        plot_stim(stim, mask=add_mask)
+    except Exception as e:
+        raise ValueError(f"Invalid parameter combination: {e}") from None
+
+# Set interactivity
+out = iw.interactive_output(
+    show_circular_generalized,
+    {
+        "height": w_height,
+        "width": w_width,
+        "ppd": w_ppd,
+        "radius1": w_radius1,
+        "radius2": w_radius2,
+        "radius3": w_radius3,
+        "int1": w_int1,
+        "int2": w_int2,
+        "int3": w_int3,
+        "intensity_background": w_int_back,
+        "origin": w_ori,
+        "add_mask": w_mask,
+        "target_indices": w_tidx,
+        "intensity_target": w_tint,
+    },
+)
+
+# Show
+display(ui, out)
+```
+
 ## Two-sided rings
 {py:func}`stimupy.stimuli.rings.circular_two_sided`
 
@@ -140,17 +233,20 @@ w_int2 = iw.FloatSlider(value=0, min=0, max=1, description="int-ring2")
 w_int_back = iw.FloatSlider(value=0.5, min=0, max=1, description="int background")
 
 w_ori = iw.Dropdown(value="mean", options=['mean', 'corner', 'center'], description="origin")
+w_clip = iw.ToggleButton(value=False, disabled=False, description="clip")
 w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'ring_mask'], description="add mask")
 
 w_tidx = iw.IntSlider(value=2, min=0, max=10, description="target idx")
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="left target int")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="right target int")
+
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
 b_geometry = iw.HBox([w_freq, w_phase])
 b_intensities = iw.HBox([w_int1, w_int2, w_int_back])
-b_add = iw.HBox([w_ori, w_mask])
-b_target = iw.HBox([w_tidx, w_tint])
+b_add = iw.HBox([w_ori, w_mask, w_clip])
+b_target = iw.HBox([w_tidx, w_tint_l, w_tint_r])
 ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
 
 # Function for showing stim
@@ -166,7 +262,9 @@ def show_circular_two_sided(
     origin=None,
     add_mask=False,
     target_indices=None,
-    intensity_target=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
+    clip=True,
 ):
     try:
         stim = circular_two_sided(
@@ -174,11 +272,12 @@ def show_circular_two_sided(
             ppd=ppd,
             frequency=frequency,
             phase_shift=phase_shift,
-            intensity_rings=(int1, int2),
+            intensity_rings=((int1, int2),(int2, int1)),
             intensity_background=intensity_background,
             origin=origin,
             target_indices=target_indices,
-            intensity_target=intensity_target,
+            intensity_target=(intensity_target_l, intensity_target_r),
+            clip=clip,
         )
         plot_stim(stim, mask=add_mask)
     except Exception as e:
@@ -199,7 +298,9 @@ out = iw.interactive_output(
         "origin": w_ori,
         "add_mask": w_mask,
         "target_indices": w_tidx,
-        "intensity_target": w_tint,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
+        "clip": w_clip,
     },
 )
 
@@ -412,17 +513,19 @@ w_int1 = iw.FloatSlider(value=1, min=0, max=1, description="int-ring1")
 w_int2 = iw.FloatSlider(value=0, min=0, max=1, description="int-ring2")
 w_int_back = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity background")
 
+w_clip = iw.ToggleButton(value=False, disabled=False, description="clip")
 w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'frame_mask'], description="add mask")
 
 w_tidx = iw.IntSlider(value=1, min=0, max=10, description="target idx")
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="left target int")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="right target int")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
 b_geometry = iw.HBox([w_freq, w_phase])
 b_intensities = iw.HBox([w_int1, w_int2, w_int_back])
-b_target = iw.HBox([w_tidx, w_tint])
-b_add = iw.HBox([w_mask])
+b_target = iw.HBox([w_tidx, w_tint_l, w_tint_r])
+b_add = iw.HBox([w_clip, w_mask])
 ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
 
 # Function for showing stim
@@ -437,7 +540,9 @@ def show_rectangular_two_sided(
     intensity_background=None,
     add_mask=False,
     target_indices=None,
-    intensity_target=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
+    clip=None,
 ):
     try:
         stim = rectangular_two_sided(
@@ -445,10 +550,11 @@ def show_rectangular_two_sided(
             ppd=ppd,
             frequency=frequency,
             phase_shift=phase_shift,
-            intensity_frames=(int1, int2),
+            intensity_frames=((int1, int2),(int2, int1)),
             intensity_background=intensity_background,
             target_indices=target_indices,
-            intensity_target=intensity_target,
+            intensity_target=(intensity_target_l, intensity_target_r),
+            clip=clip,
         )
         plot_stim(stim, mask=add_mask)
     except Exception as e:
@@ -467,8 +573,10 @@ out = iw.interactive_output(
         "int2": w_int2,
         "intensity_background": w_int_back,
         "add_mask": w_mask,
-        "intensity_target": w_tint,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
         "target_indices": w_tidx,
+        "clip": w_clip,
     },
 )
 
diff --git a/docs/reference/demos/stimuli/sbcs.md b/docs/reference/demos/stimuli/sbcs.md
index 434380b0..318b39a1 100644
--- a/docs/reference/demos/stimuli/sbcs.md
+++ b/docs/reference/demos/stimuli/sbcs.md
@@ -73,6 +73,7 @@ def show_generalized(
     intensity_background=None,
     intensity_target=None,
     add_mask=False,
+    rotation=0.0,
 ):
     try:
         stim = generalized(
@@ -82,6 +83,7 @@ def show_generalized(
             target_position=(target_y,target_x),
             intensity_background=intensity_background,
             intensity_target=intensity_target,
+            rotation=rotation,
         )
         plot_stim(stim, mask=add_mask)
     except Exception as e:
@@ -101,6 +103,7 @@ out = iw.interactive_output(
         "intensity_background": w_int_back,
         "intensity_target": w_tint,
         "add_mask": w_mask,
+        "rotation": w_rot,
     },
 )
 
@@ -179,13 +182,180 @@ out = iw.interactive_output(
 display(ui, out)
 ```
 
-## Two sided
-{py:func}`stimupy.stimuli.sbcs.two_sided`
+## Square
+{py:func}`stimupy.stimuli.sbcs.square`
 
 ```{code-cell} ipython3
 import ipywidgets as iw
 from stimupy.utils import plot_stim
-from stimupy.stimuli.sbcs import two_sided
+from stimupy.stimuli.sbcs import square
+
+# Define widgets
+w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
+w_width = iw.IntSlider(value=10, min=1, max=20, description="width [deg]")
+w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
+
+w_radius_target = iw.FloatSlider(value=1, min=0, max=2, description="target radius [deg]")
+w_radius_surround = iw.FloatSlider(value=2, min=1, max=3, description="surround radius [deg]")
+
+w_int_surround = iw.FloatSlider(value=1, min=0, max=1, description="int surround")
+w_int_back = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity background")
+
+w_ori = iw.Dropdown(value="center", options=['mean', 'corner', 'center'], description="origin")
+w_rot = iw.FloatSlider(value=0, min=0, max=360, description="rotation [deg]")
+w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'frame_mask'], description="add mask")
+
+w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+
+# Layout
+b_im_size = iw.HBox([w_height, w_width, w_ppd])
+b_geometry = iw.HBox([w_radius_target, w_radius_surround])
+b_intensities = iw.HBox([w_int_surround, w_int_back])
+b_target = iw.HBox([w_tint])
+b_add = iw.HBox([w_ori, w_rot, w_mask])
+ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
+
+# Function for showing stim
+def show_square(
+    height=None,
+    width=None,
+    ppd=None,
+    target_radius=None,
+    surround_radius=None,
+    intensity_surround=None,
+    intensity_background=None,
+    origin=None,
+    add_mask=False,
+    intensity_target=None,
+    rotation=0.0,
+):
+    try:
+        stim = square(
+            visual_size=(height, width),
+            ppd=ppd,
+            target_radius=target_radius,
+            surround_radius=surround_radius,
+            intensity_surround=intensity_surround,
+            intensity_background=intensity_background,
+            origin=origin,
+            intensity_target=intensity_target,
+            rotation=rotation,
+        )
+        plot_stim(stim, mask=add_mask)
+    except Exception as e:
+        raise ValueError(f"Invalid parameter combination: {e}") from None
+
+# Set interactivity
+out = iw.interactive_output(
+    show_square,
+    {
+        "height": w_height,
+        "width": w_width,
+        "ppd": w_ppd,
+        "target_radius": w_radius_target,
+        "surround_radius": w_radius_surround,
+        "intensity_surround": w_int_surround,
+        "intensity_background": w_int_back,
+        "origin": w_ori,
+        "add_mask": w_mask,
+        "intensity_target": w_tint,
+        "rotation": w_rot,
+    },
+)
+
+# Show
+display(ui, out)
+```
+
+
+## Circular
+{py:func}`stimupy.stimuli.sbcs.circular`
+
+```{code-cell} ipython3
+import ipywidgets as iw
+from stimupy.utils import plot_stim
+from stimupy.stimuli.sbcs import circular
+
+# Define widgets
+w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
+w_width = iw.IntSlider(value=10, min=1, max=20, description="width [deg]")
+w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
+
+w_radius_target = iw.FloatSlider(value=1, min=0, max=2, description="target radius [deg]")
+w_radius_surround = iw.FloatSlider(value=2, min=1, max=3, description="surround radius [deg]")
+
+w_int_surround = iw.FloatSlider(value=1, min=0, max=1, description="int surround")
+w_int_back = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity background")
+
+w_ori = iw.Dropdown(value="center", options=['mean', 'corner', 'center'], description="origin")
+w_mask = iw.Dropdown(value=None, options=[None, 'target_mask', 'ring_mask'], description="add mask")
+
+w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+
+# Layout
+b_im_size = iw.HBox([w_height, w_width, w_ppd])
+b_geometry = iw.HBox([w_radius_target, w_radius_surround])
+b_intensities = iw.HBox([w_int_surround, w_int_back])
+b_target = iw.HBox([w_tint])
+b_add = iw.HBox([w_ori, w_mask])
+ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_target, b_add])
+
+# Function for showing stim
+def show_circular(
+    height=None,
+    width=None,
+    ppd=None,
+    target_radius=None,
+    surround_radius=None,
+    intensity_surround=None,
+    intensity_background=None,
+    origin=None,
+    add_mask=False,
+    intensity_target=None,
+):
+    try:
+        stim = circular(
+            visual_size=(height, width),
+            ppd=ppd,
+            target_radius=target_radius,
+            surround_radius=surround_radius,
+            intensity_surround=intensity_surround,
+            intensity_background=intensity_background,
+            origin=origin,
+            intensity_target=intensity_target,
+        )
+        plot_stim(stim, mask=add_mask)
+    except Exception as e:
+        raise ValueError(f"Invalid parameter combination: {e}") from None
+
+# Set interactivity
+out = iw.interactive_output(
+    show_circular,
+    {
+        "height": w_height,
+        "width": w_width,
+        "ppd": w_ppd,
+        "target_radius": w_radius_target,
+        "surround_radius": w_radius_surround,
+        "intensity_surround": w_int_surround,
+        "intensity_background": w_int_back,
+        "origin": w_ori,
+        "add_mask": w_mask,
+        "intensity_target": w_tint,
+    },
+)
+
+# Show
+display(ui, out)
+```
+
+## Basic, Two sided
+{py:func}`stimupy.stimuli.sbcs.basic_two_sided`
+
+```{code-cell} ipython3
+import ipywidgets as iw
+from stimupy.utils import plot_stim
+from stimupy.stimuli.sbcs import basic_two_sided
 
 # Define widgets
 w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
@@ -195,7 +365,8 @@ w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
 w_t_height = iw.IntSlider(value=3, min=1, max=6, description="target height [deg]")
 w_t_width = iw.IntSlider(value=3, min=1, max=6, description="target width [deg]")
 
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target left")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target right")
 w_int_back_l = iw.FloatSlider(value=0., min=0, max=1, description="intensity left background")
 w_int_back_r = iw.FloatSlider(value=1., min=0, max=1, description="intensity right background")
 
@@ -205,8 +376,8 @@ w_mask = iw.ToggleButton(value=False, disabled=False, description="add mask")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
-b_t_size = iw.HBox([w_t_height, w_t_width, w_rot])
-b_intensities = iw.HBox([w_tint, w_int_back_l, w_int_back_r])
+b_t_size = iw.HBox([w_t_height, w_t_width])
+b_intensities = iw.HBox([w_tint_l, w_tint_r, w_int_back_l, w_int_back_r])
 b_add = iw.HBox([w_mask])
 ui = iw.VBox([b_im_size, b_t_size, b_intensities, b_add])
 
@@ -219,16 +390,17 @@ def show_two_sided(
     target_width=None,
     intensity_background_l=None,
     intensity_background_r=None,
-    intensity_target=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
     add_mask=False,
 ):
     try:
-        stim = two_sided(
+        stim = basic_two_sided(
             visual_size=(height, width),
             ppd=ppd,
-            target_size=(target_height,target_width),
-            intensity_backgrounds=(intensity_background_l, intensity_background_r),
-            intensity_target=intensity_target,
+            target_size=((target_height,target_width),(target_height,target_width)),
+            intensity_background=(intensity_background_l, intensity_background_r),
+            intensity_target=(intensity_target_l, intensity_target_r),
         )
         plot_stim(stim, mask=add_mask)
     except Exception as e:
@@ -245,7 +417,8 @@ out = iw.interactive_output(
         "target_width": w_t_width,
         "intensity_background_l": w_int_back_l,
         "intensity_background_r": w_int_back_r,
-        "intensity_target": w_tint,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
         "add_mask": w_mask,
     },
 )
@@ -254,6 +427,299 @@ out = iw.interactive_output(
 display(ui, out)
 ```
 
+## Generalized, Two sided
+{py:func}`stimupy.stimuli.sbcs.generalized_two_sided`
+
+```{code-cell} ipython3
+import ipywidgets as iw
+from stimupy.utils import plot_stim
+from stimupy.stimuli.sbcs import generalized_two_sided
+
+# Define widgets
+w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
+w_width = iw.IntSlider(value=20, min=1, max=20, description="width [deg]")
+w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
+
+w_t_height_l = iw.IntSlider(value=3, min=1, max=6, description="target height left [deg]")
+w_t_width_l = iw.IntSlider(value=3, min=1, max=6, description="target width right [deg]")
+w_t_height_r = iw.IntSlider(value=3, min=1, max=6, description="target height left [deg]")
+w_t_width_r = iw.IntSlider(value=3, min=1, max=6, description="target width right [deg]")
+
+w_t_y_l = iw.FloatSlider(value=5.0, min=0.0, max=20, description="y left target [deg]")
+w_t_x_l = iw.FloatSlider(value=5.0, min=0.0, max=20, description="x left target [deg]") 
+w_t_y_r = iw.FloatSlider(value=3.0, min=0.0, max=20, description="y right target [deg]")
+w_t_x_r = iw.FloatSlider(value=3.0, min=0.0, max=20, description="x right target [deg]")
+
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target left")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target right")
+w_int_back_l = iw.FloatSlider(value=0., min=0, max=1, description="intensity left background")
+w_int_back_r = iw.FloatSlider(value=1., min=0, max=1, description="intensity right background")
+
+
+w_mask = iw.ToggleButton(value=False, disabled=False, description="add mask")
+
+
+# Layout
+b_im_size = iw.HBox([w_height, w_width, w_ppd])
+b_t_size_l = iw.HBox([w_t_height_l, w_t_width_l])
+b_t_size_r = iw.HBox([w_t_height_r, w_t_width_r])
+b_t_pos_l = iw.HBox([w_t_y_l, w_t_x_l])
+b_t_pos_r = iw.HBox([w_t_y_r, w_t_x_r])
+b_intensities_l = iw.HBox([w_tint_l, w_int_back_l])
+b_intensities_r = iw.HBox([w_tint_r, w_int_back_r])
+
+b_left = iw.VBox([b_t_size_l, b_t_pos_l, b_intensities_l])
+b_right = iw.VBox([b_t_size_r, b_t_pos_r, b_intensities_r])
+
+b_add = iw.HBox([w_mask])
+ui = iw.VBox([b_im_size, b_left, b_right, b_add])
+
+# Function for showing stim
+def show_gen_two_sided(
+    height=None,
+    width=None,
+    ppd=None,
+    target_height_l=None,
+    target_width_l=None,
+    target_height_r=None,
+    target_width_r=None,
+    target_pos_x_l=None,
+    target_pos_y_l=None,
+    target_pos_x_r=None,
+    target_pos_y_r=None,
+    intensity_background_l=None,
+    intensity_background_r=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
+    add_mask=False,
+):
+    try:
+        stim = generalized_two_sided(
+            visual_size=(height, width),
+            ppd=ppd,
+            target_size=((target_height_l,target_width_l),(target_height_r,target_width_r)),
+            target_position=((target_pos_y_l, target_pos_x_l), (target_pos_y_r, target_pos_x_r)),
+            intensity_background=(intensity_background_l, intensity_background_r),
+            intensity_target=(intensity_target_l, intensity_target_r),
+        )
+        plot_stim(stim, mask=add_mask)
+    except Exception as e:
+        raise ValueError(f"Invalid parameter combination: {e}") from None
+
+# Set interactivity
+out = iw.interactive_output(
+    show_gen_two_sided,
+    {
+        "height": w_height,
+        "width": w_width,
+        "ppd": w_ppd,
+        "target_height_l": w_t_height_l,
+        "target_width_l": w_t_width_l,
+        "target_height_r": w_t_height_r,
+        "target_width_r": w_t_width_r,
+        "target_pos_x_l": w_t_x_l,
+        "target_pos_y_l": w_t_y_l,
+        "target_pos_x_r": w_t_x_r,
+        "target_pos_y_r": w_t_y_r,
+        "intensity_background_l": w_int_back_l,
+        "intensity_background_r": w_int_back_r,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
+        "add_mask": w_mask,
+    },
+)
+
+# Show
+display(ui, out)
+```
+
+## Square, Two sided
+{py:func}`stimupy.stimuli.sbcs.square_two_sided`
+
+```{code-cell} ipython3
+import ipywidgets as iw
+from stimupy.utils import plot_stim
+from stimupy.stimuli.sbcs import square_two_sided
+
+# Define widgets
+w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
+w_width = iw.IntSlider(value=10, min=1, max=20, description="width [deg]")
+w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
+
+w_t_radius_l = iw.IntSlider(value=1, min=1, max=6, description="target radius left [deg]")
+w_t_radius_r = iw.IntSlider(value=1, min=1, max=6, description="target radius right [deg]")
+w_radius_surround_l = iw.IntSlider(value=2, min=1, max=6, description="surround radius left [deg]")
+w_radius_surround_r = iw.IntSlider(value=2, min=1, max=6, description="surround radius right [deg]")
+
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target left")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target right")
+w_int_surround_l = iw.FloatSlider(value=0., min=0, max=1, description="intensity left surround")
+w_int_surround_r = iw.FloatSlider(value=1., min=0, max=1, description="intensity right surround")
+w_int_background = iw.FloatSlider(value=.5, min=0, max=1, description="intenisty background")
+
+w_rot_l = iw.FloatSlider(value=0, min=0, max=360, description="rotation left [deg]")
+w_rot_r = iw.FloatSlider(value=0, min=0, max=360, description="rotation right [deg]")
+w_mask = iw.ToggleButton(value=False, disabled=False, description="add mask")
+
+
+# Layout
+b_im_size = iw.HBox([w_height, w_width, w_ppd])
+b_t_size = iw.HBox([w_t_radius_l, w_t_radius_l])
+b_s_size = iw.HBox([w_radius_surround_l, w_radius_surround_r])
+b_intensities = iw.HBox([w_tint_l, w_tint_r, w_int_back_l, w_int_back_r])
+b_add = iw.HBox([w_rot_l, w_rot_r, w_mask])
+ui = iw.VBox([b_im_size, b_t_size, b_s_size, b_intensities, b_add])
+
+# Function for showing stim
+def show_square_two_sided(
+    height=None,
+    width=None,
+    ppd=None,
+    target_radius_l=None,
+    target_radius_r=None,
+    surround_radius_l=None,
+    surround_radius_r=None,
+    intensity_surround_l=None,
+    intensity_surround_r=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
+    intensity_background=None,
+    add_mask=False,
+    rotation_l=None,
+    rotation_r=None,
+):
+    try:
+        stim = square_two_sided(
+            visual_size=(height, width),
+            ppd=ppd,
+            target_radius=(target_radius_l,target_radius_r),
+            surround_radius=(surround_radius_l, surround_radius_r),
+            intensity_surround=(intensity_surround_l, intensity_surround_r),
+            intensity_target=(intensity_target_l, intensity_target_r),
+            intensity_background=intensity_background,
+            rotation=(rotation_l, rotation_r),
+        )
+        plot_stim(stim, mask=add_mask)
+    except Exception as e:
+        raise ValueError(f"Invalid parameter combination: {e}") from None
+
+# Set interactivity
+out = iw.interactive_output(
+    show_square_two_sided,
+    {
+        "height": w_height,
+        "width": w_width,
+        "ppd": w_ppd,
+        "target_radius_l": w_t_radius_l,
+        "target_radius_r": w_t_radius_r,
+        "surround_radius_l": w_radius_surround_l,
+        "surround_radius_r": w_radius_surround_r,
+        "intensity_surround_l": w_int_surround_l,
+        "intensity_surround_r": w_int_surround_r,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
+        "intensity_background": w_int_background,
+        "rotation_l": w_rot_l,
+        "rotation_r": w_rot_r,
+        "add_mask": w_mask,
+    },
+)
+
+# Show
+display(ui, out)
+```
+
+## Circular, Two sided
+{py:func}`stimupy.stimuli.sbcs.circular_two_sided`
+
+```{code-cell} ipython3
+import ipywidgets as iw
+from stimupy.utils import plot_stim
+from stimupy.stimuli.sbcs import circular_two_sided
+
+# Define widgets
+w_height = iw.IntSlider(value=10, min=1, max=20, description="height [deg]")
+w_width = iw.IntSlider(value=10, min=1, max=20, description="width [deg]")
+w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
+
+w_t_radius_l = iw.IntSlider(value=1, min=1, max=6, description="target radius left [deg]")
+w_t_radius_r = iw.IntSlider(value=1, min=1, max=6, description="target radius right [deg]")
+w_radius_surround_l = iw.IntSlider(value=2, min=1, max=6, description="surround radius left [deg]")
+w_radius_surround_r = iw.IntSlider(value=2, min=1, max=6, description="surround radius right [deg]")
+
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target left")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="intensity target right")
+w_int_surround_l = iw.FloatSlider(value=0., min=0, max=1, description="intensity left surround")
+w_int_surround_r = iw.FloatSlider(value=1., min=0, max=1, description="intensity right surround")
+w_int_background = iw.FloatSlider(value=.5, min=0, max=1, description="intenisty background")
+
+
+w_mask = iw.ToggleButton(value=False, disabled=False, description="add mask")
+
+
+# Layout
+b_im_size = iw.HBox([w_height, w_width, w_ppd])
+b_t_size = iw.HBox([w_t_radius_l, w_t_radius_l])
+b_s_size = iw.HBox([w_radius_surround_l, w_radius_surround_r])
+b_intensities = iw.HBox([w_tint_l, w_tint_r, w_int_back_l, w_int_back_r])
+b_add = iw.HBox([w_mask])
+ui = iw.VBox([b_im_size, b_t_size, b_s_size, b_intensities, b_add])
+
+# Function for showing stim
+def show_circular_two_sided(
+    height=None,
+    width=None,
+    ppd=None,
+    target_radius_l=None,
+    target_radius_r=None,
+    surround_radius_l=None,
+    surround_radius_r=None,
+    intensity_surround_l=None,
+    intensity_surround_r=None,
+    intensity_target_l=None,
+    intensity_target_r=None,
+    intensity_background=None,
+    add_mask=False,
+):
+    try:
+        stim = circular_two_sided(
+            visual_size=(height, width),
+            ppd=ppd,
+            target_radius=(target_radius_l,target_radius_r),
+            surround_radius=(surround_radius_l, surround_radius_r),
+            intensity_surround=(intensity_surround_l, intensity_surround_r),
+            intensity_target=(intensity_target_l, intensity_target_r),
+            intensity_background=intensity_background,
+        )
+        plot_stim(stim, mask=add_mask)
+    except Exception as e:
+        raise ValueError(f"Invalid parameter combination: {e}") from None
+
+# Set interactivity
+out = iw.interactive_output(
+    show_circular_two_sided,
+    {
+        "height": w_height,
+        "width": w_width,
+        "ppd": w_ppd,
+        "target_radius_l": w_t_radius_l,
+        "target_radius_r": w_t_radius_r,
+        "surround_radius_l": w_radius_surround_l,
+        "surround_radius_r": w_radius_surround_r,
+        "intensity_surround_l": w_int_surround_l,
+        "intensity_surround_r": w_int_surround_r,
+        "intensity_target_l": w_tint_l,
+        "intensity_target_r": w_tint_r,
+        "intensity_background": w_int_background,
+        "add_mask": w_mask,
+    },
+)
+
+# Show
+display(ui, out)
+```
+
+
 ## With dots
 {py:func}`stimupy.stimuli.sbcs.with_dots`
 
@@ -410,7 +876,7 @@ def show_two_sided_with_dots(
             n_dots=(ndotsy, ndotsx),
             dot_radius=dot_radius,
             distance=dot_distance,
-            intensity_backgrounds=(intensity_background_l, intensity_background_r),
+            intensity_background=(intensity_background_l, intensity_background_r),
             intensity_target=intensity_target,
             intensity_dots=(intensity_dots_l, intensity_dots_r),
         )
@@ -600,7 +1066,7 @@ def show_two_sided_dotted(
             n_dots=(ndotsy, ndotsx),
             dot_radius=dot_radius,
             distance=dot_distance,
-            intensity_backgrounds=(intensity_background_l, intensity_background_r),
+            intensity_background=(intensity_background_l, intensity_background_r),
             intensity_target=intensity_target,
             intensity_dots=(intensity_dots_l, intensity_dots_r),
         )
diff --git a/docs/reference/demos/stimuli/todorovics.md b/docs/reference/demos/stimuli/todorovics.md
index def7415f..f6737d5a 100644
--- a/docs/reference/demos/stimuli/todorovics.md
+++ b/docs/reference/demos/stimuli/todorovics.md
@@ -251,7 +251,8 @@ w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
 
 w_theight = iw.FloatSlider(value=2, min=0, max=4, description="target height [deg]")
 w_twidth = iw.FloatSlider(value=2, min=0, max=4, description="target width [deg]")
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="target int")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="left target int")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="right target int")
 
 w_cheight = iw.FloatSlider(value=2, min=0, max=4, description="cover height [deg]")
 w_cwidth = iw.FloatSlider(value=2, min=0, max=4, description="cover width [deg]")
@@ -267,7 +268,7 @@ w_mask = iw.Dropdown(value=None, options=[None, 'target_mask'], description="add
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
-b_tsize = iw.HBox([w_theight, w_twidth, w_tint])
+b_tsize = iw.HBox([w_theight, w_twidth, w_tint_l, w_tint_r])
 b_csize = iw.HBox([w_cheight, w_cwidth, w_cint1, w_cint2])
 b_c1pos = iw.HBox([w_coffx, w_coffy])
 b_intensities = iw.HBox([w_int_back])
@@ -281,7 +282,8 @@ def show_two_sided_rectangle(
     ppd=None,
     theight=None,
     twidth=None,
-    tint=None,
+    tint_l=None,
+    tint_r=None,
     cheight=None,
     cwidth=None,
     cint1=None,
@@ -296,11 +298,11 @@ def show_two_sided_rectangle(
         stim = rectangle_two_sided(
             visual_size=(height, width),
             ppd=ppd,
-            target_size=(theight, twidth),
-            covers_size=(cheight, cwidth),
-            covers_offset=(coffy, coffx),
-            intensity_backgrounds=(intback1, intback2),
-            intensity_target=tint,
+            target_size=((theight, twidth),(theight, twidth)),
+            covers_size=((cheight, cwidth),(cheight, cwidth)),
+            covers_offset=((coffy, coffx),(coffy, coffx)),
+            intensity_background=(intback1, intback2),
+            intensity_target=(tint_l, tint_r),
             intensity_covers=(cint1, cint2),
         )
         plot_stim(stim, mask=add_mask)
@@ -316,7 +318,8 @@ out = iw.interactive_output(
         "ppd": w_ppd,
         "theight": w_theight,
         "twidth": w_twidth,
-        "tint": w_tint,
+        "tint_l": w_tint_l,
+        "tint_r": w_tint_r,
         "cheight": w_cheight,
         "cwidth": w_cwidth,
         "cint1": w_cint1,
@@ -539,20 +542,22 @@ w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
 w_theight = iw.FloatSlider(value=4, min=0, max=8, description="cross height [deg]")
 w_twidth = iw.FloatSlider(value=4, min=0, max=8, description="cross width [deg]")
 w_tthick = iw.FloatSlider(value=2, min=0, max=4, description="cross thickness [deg]")
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="cross int")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="left cross int")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="right cross int")
+
 
 w_cheight = iw.FloatSlider(value=2, min=0, max=4, description="cover height [deg]")
 w_cwidth = iw.FloatSlider(value=2, min=0, max=4, description="cover width [deg]")
 w_cint1 = iw.FloatSlider(value=1, min=0, max=1, description="cover int1")
 w_cint2 = iw.FloatSlider(value=0, min=0, max=1, description="cover int2")
 
-w_int_back1 = iw.FloatSlider(value=0., min=0, max=1, description="int back1")
-w_int_back2 = iw.FloatSlider(value=1., min=0, max=1, description="int back12")
+w_int_back_l = iw.FloatSlider(value=0., min=0, max=1, description="int back1")
+w_int_back_r = iw.FloatSlider(value=1., min=0, max=1, description="int back12")
 w_mask = iw.Dropdown(value=None, options=[None, 'target_mask'], description="add mask")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
-b_tsize = iw.HBox([w_theight, w_twidth, w_tthick, w_tint])
+b_tsize = iw.HBox([w_theight, w_twidth, w_tthick, w_tint_l, w_tint_r])
 b_csize = iw.HBox([w_cheight, w_cwidth, w_cint1, w_cint2])
 b_intensities = iw.HBox([w_int_back1, w_int_back2])
 b_add = iw.HBox([w_mask])
@@ -565,25 +570,26 @@ def show_two_sided_cross(
     ppd=None,
     theight=None,
     twidth=None,
-    tint=None,
+    tint_l=None,
+    tint_r=None,
     tthick=None,
     cheight=None,
     cwidth=None,
     cint1=None,
     cint2=None,
-    intback1=None,
-    intback2=None,
+    int_back_l=None,
+    int_back_r=None,
     add_mask=False,
 ):
     try:
         stim = cross_two_sided(
             visual_size=(height, width),
             ppd=ppd,
-            cross_size=(theight, twidth),
+            cross_size=((theight, twidth),(theight, twidth)),
             cross_thickness=tthick,
-            covers_size=(cheight, cwidth),
-            intensity_backgrounds=(intback1, intback2),
-            intensity_target=tint,
+            covers_size=((cheight, cwidth),(cheight, cwidth)),
+            intensity_background=(int_back_l, int_back_r),
+            intensity_target=(tint_l, tint_r),
             intensity_covers=(cint1, cint2),
         )
         plot_stim(stim, mask=add_mask)
@@ -599,14 +605,15 @@ out = iw.interactive_output(
         "ppd": w_ppd,
         "theight": w_theight,
         "twidth": w_twidth,
-        "tint": w_tint,
+        "tint_l": w_tint_l,
+        "tint_r": w_tint_r,
         "tthick": w_tthick,
         "cheight": w_cheight,
         "cwidth": w_cwidth,
         "cint1": w_cint1,
         "cint2": w_cint2,
-        "intback1": w_int_back1,
-        "intback2": w_int_back2,
+        "int_back_l": w_int_back_l,
+        "int_back_r": w_int_back_r,
         "add_mask": w_mask,
     },
 )
@@ -710,18 +717,19 @@ w_ppd = iw.IntSlider(value=20, min=1, max=40, description="ppd")
 w_theight = iw.FloatSlider(value=4, min=0, max=8, description="cross height [deg]")
 w_twidth = iw.FloatSlider(value=4, min=0, max=8, description="cross width [deg]")
 w_tthick = iw.FloatSlider(value=2, min=0, max=4, description="cross thickness [deg]")
-w_tint = iw.FloatSlider(value=0.5, min=0, max=1, description="cross int")
+w_tint_l = iw.FloatSlider(value=0.5, min=0, max=1, description="left cross int")
+w_tint_r = iw.FloatSlider(value=0.5, min=0, max=1, description="right cross int")
 
-w_cint1 = iw.FloatSlider(value=1, min=0, max=1, description="cover1 int")
-w_cint2 = iw.FloatSlider(value=0, min=0, max=1, description="cover2 int")
+w_cint_l = iw.FloatSlider(value=1, min=0, max=1, description="cover1 int")
+w_cint_r = iw.FloatSlider(value=0, min=0, max=1, description="cover2 int")
 
-w_int_back1 = iw.FloatSlider(value=0., min=0, max=1, description="int back1")
-w_int_back2 = iw.FloatSlider(value=1., min=0, max=1, description="int back2")
+w_int_back_l = iw.FloatSlider(value=0., min=0, max=1, description="int back1")
+w_int_back_r = iw.FloatSlider(value=1., min=0, max=1, description="int back2")
 w_mask = iw.Dropdown(value=None, options=[None, 'target_mask'], description="add mask")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
-b_tsize = iw.HBox([w_theight, w_twidth, w_tthick, w_tint])
+b_tsize = iw.HBox([w_theight, w_twidth, w_tthick, w_tint_l, w_tint_r])
 b_csize = iw.HBox([w_cint1, w_cint2])
 b_intensities = iw.HBox([w_int_back1, w_int_back2])
 b_add = iw.HBox([w_mask])
@@ -734,23 +742,24 @@ def show_two_sided_equal(
     ppd=None,
     theight=None,
     twidth=None,
-    tint=None,
+    tint_l=None,
+    tint_r=None,
     tthick=None,
-    cint1=None,
-    cint2=None,
-    intback1=None,
-    intback2=None,
+    cint_l=None,
+    cint_r=None,
+    int_back_l=None,
+    int_back_r=None,
     add_mask=False,
 ):
     try:
         stim = equal_two_sided(
             visual_size=(height, width),
             ppd=ppd,
-            cross_size=(theight, twidth),
+            cross_size=((theight, twidth),(theight, twidth)),
             cross_thickness=tthick,
-            intensity_backgrounds=(intback1, intback2),
-            intensity_target=tint,
-            intensity_covers=(cint1, cint2),
+            intensity_background=(int_back_l, int_back_r),
+            intensity_target=(tint_l, tint_r),
+            intensity_covers=(cint_l, cint_r),
         )
         plot_stim(stim, mask=add_mask)
     except Exception as e:
@@ -765,12 +774,13 @@ out = iw.interactive_output(
         "ppd": w_ppd,
         "theight": w_theight,
         "twidth": w_twidth,
-        "tint": w_tint,
+        "tint_l": w_tint_l,
+        "tint_r": w_tint_r,
         "tthick": w_tthick,
-        "cint1": w_cint1,
-        "intback1": w_int_back1,
-        "cint2": w_cint2,
-        "intback2": w_int_back2,
+        "cint_l": w_cint1,
+        "int_back_l": w_int_back_l,
+        "cint_r": w_cint2,
+        "int_back_r": w_int_back_r,
         "add_mask": w_mask,
     },
 )
diff --git a/docs/reference/demos/stimuli/whites.md b/docs/reference/demos/stimuli/whites.md
index 5574b101..047e8710 100644
--- a/docs/reference/demos/stimuli/whites.md
+++ b/docs/reference/demos/stimuli/whites.md
@@ -508,13 +508,13 @@ w_mask = iw.ToggleButton(value=False, disabled=False, description="add mask")
 
 # Layout
 b_im_size = iw.HBox([w_height, w_width, w_ppd])
-b_geometry = iw.HBox([w_freq, w_rot])
+b_geometry = iw.HBox([w_freq])
 b_intensities = iw.HBox([w_int1, w_int2])
 b_toff = w_toff
 b_target = iw.HBox([w_tidx, w_tint, w_theights])
 b_target2 = iw.HBox([w_tidx2, w_tint2, w_theights2])
 b_stripe = iw.HBox([w_sint1, w_sint2])
-b_add = iw.HBox([w_ori, w_period, w_mask])
+b_add = iw.HBox([w_period, w_mask])
 ui = iw.VBox([b_im_size, b_geometry, b_intensities, b_toff, b_target, b_target2, b_stripe, b_add])
 
 # Function for showing stim
diff --git a/stimupy/components/__init__.py b/stimupy/components/__init__.py
index 44fa60cf..5f5ce158 100644
--- a/stimupy/components/__init__.py
+++ b/stimupy/components/__init__.py
@@ -3,6 +3,7 @@
 import numpy as np
 
 from stimupy.components import *  # angulars, edges, frames, gaussians, lines, radials, shapes, waves
+from stimupy.components import texts
 from stimupy.utils import resolution
 
 __all__ = [
@@ -20,6 +21,7 @@
     "gaussians",
     "lines",
     "shapes",
+    "texts",
     "waves",
 ]
 
diff --git a/stimupy/components/angulars.py b/stimupy/components/angulars.py
index 2684fb55..43c7cfac 100644
--- a/stimupy/components/angulars.py
+++ b/stimupy/components/angulars.py
@@ -39,7 +39,7 @@ def mask_angle(
     Returns
     -------
     dict[str, Any]
-        dict with boolean mask (key: bool_mask) for pixels falling in given angle,
+        dict with boolean mask (key: "segment_mask") for pixels falling in given angle,
         and additional params
     """
     stim = mask_regions(
@@ -51,8 +51,7 @@ def mask_angle(
         ppd=ppd,
         origin=origin,
     )
-    stim["wedge_mask"] = stim["mask"]
-    del stim["mask"]
+    stim["segment_mask"] = stim.pop("mask")
     return stim
 
 
@@ -131,16 +130,19 @@ def wedge(
         origin=origin,
     )
 
-    # Remove everything but wedge, and add additional args
-    stim["img"] = np.where(mask["wedge_mask"], stim["img"], intensity_background)
-    stim["wedge_mask"] = stim["ring_mask"] * mask["wedge_mask"]
+    # Combine masks
+    stim["wedge_mask"] = stim.pop("ring_mask") * mask["segment_mask"]
     stim["wedge_mask"] = np.where(stim["wedge_mask"] != 0, 1, 0).astype(int)
+
+    # Draw image
+    stim["img"] = np.where(stim["wedge_mask"], intensity_wedge, intensity_background)
+
+    # Repackage
     stim["angle"] = angle
     stim["radius"] = radius
     stim["inner_radius"] = inner_radius
     stim["intensity_wedge"] = intensity_wedge
     stim["intensity_background"] = intensity_background
-    del stim["ring_mask"]
     return stim
 
 
@@ -175,7 +177,7 @@ def mask_segments(
     -------
     dict[str, Any]
         dict with the stimulus (key: "img"),
-        mask with integer index for each segment (key: "wedge_mask"),
+        mask with integer index for each segment (key: "segment_mask"),
         and additional keys containing stimulus parameters
     """
     stim = mask_regions(
@@ -187,8 +189,7 @@ def mask_segments(
         ppd=ppd,
         origin=origin,
     )
-    stim["wedge_mask"] = stim["mask"]
-    del stim["mask"]
+    stim["segment_mask"] = stim.pop("mask")
     return stim
 
 
@@ -199,7 +200,7 @@ def segments(
     angles=None,
     rotation=0.0,
     intensity_background=0.5,
-    intensity_segments=(0, 1),
+    intensity_segments=(0.0, 1.0),
     origin="mean",
 ):
     """Generate mask with integer indices for sequential angular segments
@@ -230,7 +231,7 @@ def segments(
     -------
     dict[str, Any]
         dict with the stimulus (key: "img"),
-        mask with integer index for each segment (key: "wedge_mask"),
+        mask with integer index for each segment (key: "segment_mask"),
         and additional keys containing stimulus parameters
     """
     if angles is None:
@@ -251,7 +252,7 @@ def segments(
 
     # Draw image
     stim["img"] = draw_regions(
-        stim["wedge_mask"],
+        stim["segment_mask"],
         intensities=intensity_segments,
         intensity_background=intensity_background,
     )
diff --git a/stimupy/components/frames.py b/stimupy/components/frames.py
index 6bd84c5e..ce2c8a4e 100644
--- a/stimupy/components/frames.py
+++ b/stimupy/components/frames.py
@@ -49,8 +49,8 @@ def mask_frames(
         ppd=ppd,
         origin=origin,
     )
-    stim["frame_mask"] = stim["mask"]
-    del stim["mask"]
+    stim["frame_mask"] = stim.pop("mask")
+
     return stim
 
 
@@ -60,7 +60,7 @@ def frames(
     shape=None,
     radii=None,
     rotation=0.0,
-    intensity_frames=(1.0, 0.0),
+    intensity_frames=(0.0, 1.0),
     intensity_background=0.5,
     origin="mean",
 ):
@@ -78,13 +78,13 @@ def frames(
         radii of each frame, in degrees visual angle
     rotation : float, optional
         rotation (in degrees), counterclockwise, by default 0.0 (horizonal)
-    intensity_frames : Sequence[float, ...]
+    intensity_frames : Sequence[float, ...], optional
         intensity value for each frame, by default (1.0, 0.0).
         Can specify as many intensities as number of frame_widths;
         If fewer intensities are passed than frame_widhts, cycles through intensities
     intensity_background : float, optional
         intensity value of background, by default 0.5
-    origin : "corner", "mean" or "center"
+    origin : "corner", "mean" or "center", optional
         if "corner": set origin to upper left corner
         if "mean": set origin to hypothetical image center (default)
         if "center": set origin to real center (closest existing value to mean)
diff --git a/stimupy/components/gaussians.py b/stimupy/components/gaussians.py
index 95af1729..3eca3cbc 100644
--- a/stimupy/components/gaussians.py
+++ b/stimupy/components/gaussians.py
@@ -114,7 +114,8 @@ def overview(**kwargs):
 
     # fmt: off
     stimuli = {
-        "gaussians_gaussian": gaussian(**default_params, sigma=(3, 1.5), rotation=70),}
+        "gaussians_gaussian": gaussian(**default_params, sigma=(3, 1.5), rotation=70),
+    }
     # fmt: on
 
     return stimuli
diff --git a/stimupy/components/lines.py b/stimupy/components/lines.py
index 529a1cc4..f964dcd4 100644
--- a/stimupy/components/lines.py
+++ b/stimupy/components/lines.py
@@ -23,8 +23,8 @@ def line(
     line_length=None,
     line_width=0,
     rotation=0.0,
-    intensity_line=1,
-    intensity_background=0,
+    intensity_line=1.0,
+    intensity_background=0.0,
     origin="corner",
 ):
     """Draw a line
@@ -138,7 +138,7 @@ def dipole(
     line_width=0,
     line_gap=None,
     rotation=0.0,
-    intensity_lines=(0, 1),
+    intensity_lines=(0.0, 1.0),
 ):
     """Draw a two centered parallel lines
 
@@ -237,8 +237,8 @@ def ellipse(
     shape=None,
     radius=None,
     line_width=0,
-    intensity_line=1,
-    intensity_background=0,
+    intensity_line=1.0,
+    intensity_background=0.0,
 ):
     """Draw an ellipse
 
@@ -309,8 +309,8 @@ def circle(
     shape=None,
     radius=None,
     line_width=0,
-    intensity_line=1,
-    intensity_background=0,
+    intensity_line=1.0,
+    intensity_background=0.0,
 ):
     """Draw a circle given the input parameters
 
@@ -381,7 +381,8 @@ def overview(**kwargs):
         "lines_line": line(**default_params, **p, origin="center"),
         "lines_dipole": dipole(**default_params, **p, line_gap=1),
         "lines_circle": circle(**default_params, radius=3),
-        "lines_ellipse": ellipse(**default_params, radius=(3, 4)),}
+        "lines_ellipse": ellipse(**default_params, radius=(3, 4)),
+    }
     # fmt: on
 
     return stimuli
diff --git a/stimupy/components/radials.py b/stimupy/components/radials.py
index c2202454..23585e42 100644
--- a/stimupy/components/radials.py
+++ b/stimupy/components/radials.py
@@ -1,8 +1,8 @@
-import itertools
+import collections
 
 import numpy as np
 
-from stimupy.components import image_base, mask_regions
+from stimupy.components import draw_regions, mask_regions
 from stimupy.utils import resolution
 
 __all__ = [
@@ -76,7 +76,7 @@ def rings(
     ppd=None,
     shape=None,
     radii=None,
-    intensity_rings=(0, 1),
+    intensity_rings=(0.0, 1.0),
     intensity_background=0.5,
     origin="mean",
 ):
@@ -92,12 +92,12 @@ def rings(
         shape [height, width] of image, in pixels
     radii : Sequence[Number]
         outer radii of rings (& disc) in degree visual angle
-    intensity_rings : Sequence[Number, ...]
-        intensity value for each ring, from inside to out.
+    intensity_rings : Sequence[Number, ...], optional
+        intensity value for each ring, from inside to out, by default (0.0, 1.0)
         If fewer intensities are passed than number of radii, cycles through intensities
-    intensity_background : float (optional)
+    intensity_background : float, optional
         value of background, by default 0.5
-    origin : "corner", "mean" or "center"
+    origin : "corner", "mean" or "center", optional
         if "corner": set origin to upper left corner
         if "mean": set origin to hypothetical image center (default)
         if "center": set origin to real center (closest existing value to mean)
@@ -109,9 +109,6 @@ def rings(
         mask with integer index for each ring (key: "ring_mask"),
         and additional keys containing stimulus parameters
     """
-    if radii is None:
-        raise ValueError("rings() missing argument 'radii' which is not 'None'")
-
     # Try to resolve resolution;
     try:
         shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
@@ -127,34 +124,32 @@ def rings(
         visual_size = [x for x in visual_size if x is not None]
         visual_size = resolution.validate_visual_size(visual_size)
 
-    if not isinstance(radii, (int, float)):
+    # Resolve radii
+    if radii is None:
+        raise ValueError("rings() missing argument 'radii' which is not 'None'")
+    if isinstance(radii, collections.abc.Sequence) or isinstance(radii, np.ndarray):
         if np.diff(radii).min() < 0:
             raise ValueError("radii need to monotonically increase")
     else:
         radii = (radii,)
-    if isinstance(intensity_rings, (int, float)):
-        ints = (intensity_rings,)
-    else:
-        ints = intensity_rings
 
     # Get masks for rings
-    params = mask_rings(radii=radii, shape=shape, visual_size=visual_size, ppd=ppd, origin=origin)
-    shape = params["shape"]
+    stim = mask_rings(radii=radii, shape=shape, visual_size=visual_size, ppd=ppd, origin=origin)
 
-    # Draw rings
-    base = image_base(shape=shape, visual_size=visual_size, origin=origin)
-    distances = base["radial"]
+    # Resolve intensities
+    if not isinstance(intensity_rings, collections.abc.Sequence):
+        intensity_rings = (intensity_rings,)
 
-    img = np.ones(shape) * intensity_background
-    ints = [*itertools.islice(itertools.cycle(ints), len(radii))]
-    for radius, intensity in zip(reversed(radii), reversed(ints)):
-        img[distances < radius] = intensity
+    # Draw rings
+    stim["img"] = draw_regions(
+        stim["ring_mask"], intensity_rings, intensity_background=intensity_background
+    )
 
     # Assemble output
-    params["intensity_rings"] = intensity_rings
-    params["radii"] = radii
-    params["intensity_background"] = intensity_background
-    return {"img": img, **params}
+    stim["intensity_rings"] = intensity_rings
+    stim["radii"] = radii
+    stim["intensity_background"] = intensity_background
+    return stim
 
 
 def disc(
@@ -180,11 +175,11 @@ def disc(
         shape [height, width] of image, in pixels
     radius : Number
         outer radius of disc in degree visual angle
-    intensity_disc : Number
+    intensity_disc : Number, optional
         intensity value of disc, by default 1.0
-    intensity_background : float (optional)
+    intensity_background : float, optional
         intensity value of background, by default 0.0
-    origin : "corner", "mean" or "center"
+    origin : "corner", "mean" or "center", optional
         if "corner": set origin to upper left corner
         if "mean": set origin to hypothetical image center (default)
         if "center": set origin to real center (closest existing value to mean)
@@ -199,26 +194,26 @@ def disc(
     if radius is None:
         raise ValueError("disc() missing argument 'radius' which is not 'None'")
 
-    radius = np.array(radius)
-    intensity = np.array(intensity_disc)
+    if isinstance(radius, np.ndarray):
+        try:
+            radius = float(radius)
+        except TypeError as e:
+            raise ValueError("Can only pass 1 radius") from e
+    elif isinstance(radius, collections.abc.Sequence):
+        if len(radius) != 1:
+            raise ValueError("Can only pass 1 radius")
 
-    if radius.size != 1:
-        raise ValueError("Can only pass 1 radius")
-    if intensity.size != 1:
-        raise ValueError("Can only pass 1 intensity")
-
-    stim = ring(
-        radii=np.insert(radius, 0, 0),
-        intensity_ring=intensity,
+    stim = rings(
+        radii=radius,
+        intensity_rings=intensity_disc,
         visual_size=visual_size,
         ppd=ppd,
         intensity_background=intensity_background,
         shape=shape,
         origin=origin,
     )
-    stim["radius"] = radius
-    stim["intensity_disc"] = intensity_disc
-    del stim["radii"], stim["intensity_rings"]
+    stim["radius"] = stim.pop("radii")
+    stim["intensity_disc"] = stim.pop("intensity_rings")
     return stim
 
 
@@ -243,11 +238,11 @@ def ring(
         shape [height, width] of image, in pixels
     radii : Sequence[Number, Number]
         inner and outer radius of ring in degree visual angle
-    intensity_ring : Number
+    intensity_ring : Number, optional
         intensity value of ring, by default 1.0
-    intensity_background : float (optional)
+    intensity_background : float, optional
         intensity value of background, by default 0.0
-    origin : "corner", "mean" or "center"
+    origin : "corner", "mean" or "center", optional
         if "corner": set origin to upper left corner
         if "mean": set origin to hypothetical image center (default)
         if "center": set origin to real center (closest existing value to mean)
diff --git a/stimupy/components/shapes.py b/stimupy/components/shapes.py
index 4d96503f..2d9dfbf2 100644
--- a/stimupy/components/shapes.py
+++ b/stimupy/components/shapes.py
@@ -624,7 +624,8 @@ def overview(**kwargs):
         "shapes_disc": disc(**default_params, radius=3),
         "shapes_ring": ring(**default_params, radii=(1, 2)),
         "shapes_annulus": annulus(**default_params, radii=(1, 2)),
-        "shapes_wedge": wedge(**default_params, angle=30, radius=4),}
+        "shapes_wedge": wedge(**default_params, angle=30, radius=4),
+    }
     # fmt: on
 
     return stimuli
diff --git a/stimupy/components/texts.py b/stimupy/components/texts.py
new file mode 100644
index 00000000..5d8dc3d3
--- /dev/null
+++ b/stimupy/components/texts.py
@@ -0,0 +1,145 @@
+import numpy as np
+from PIL import Image, ImageDraw, ImageFont
+
+from stimupy.utils import pad_dict_to_visual_size, resolution
+
+__all__ = [
+    "text",
+]
+
+
+def text(
+    text,
+    visual_size=None,
+    ppd=None,
+    shape=None,
+    intensity_text=0.0,
+    intensity_background=0.5,
+    fontsize=36,
+    align="center",
+    direction="ltr",
+):
+    """Draw given text into a (numpy) image-array
+
+    If no shape is provided / can be resolved,
+    tightly fits the bounding box of the drawn text.
+
+    Parameters
+    ----------
+    text : str
+        Text to draw
+    visual_size : Sequence[Number, Number], Number, or None (default)
+        visual size [height, width] of image, in degrees visual angle
+    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
+    intensity_text : float, optional
+        intensity of text in range (0.0; 1.0), by default 0.0
+    intensity_background : float, optional
+        intensity value of background in range (0.0; 1.0), by default 0.5
+    fontsize : int, optional
+        font size, by default 36
+    align : "left", "center" (default), "right"
+        alignment of text, by default "center"
+    direction : "ltr" (default), "rtl"
+        text drawing direction, left-to-right, right-to-left
+
+    Returns
+    -------
+    dict[str, Any]
+        dict with the stimulus (key: "img"),
+        mask with integer index for the text (key: "text_mask"),
+        and additional keys containing stimulus parameters
+    """
+
+    # Try to resolve resolution
+    try:
+        shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
+    except resolution.TooManyUnknownsError:
+        shape = resolution.validate_shape(shape)
+        visual_size = resolution.validate_visual_size(visual_size)
+        ppd = resolution.validate_ppd(ppd)
+
+    # Get font
+    font = ImageDraw.ImageDraw.font
+    if not font:
+        try:
+            # Not all machines will have Arial installed...
+            font = ImageFont.truetype(
+                "arial.ttf",
+                fontsize,
+                encoding="unic",
+            )
+        except OSError:
+            font = ImageFont.load_default()
+
+    # Determine dimensions of total text
+    n_lines = len(text.split("\n"))
+    max_length = 0
+    for line in text.split("\n"):
+        max_length = max(int(font.getlength(line)), max_length)
+    _, top, _, bottom = font.getbbox(text)
+    text_width = max_length
+    text_height = int(top + bottom) * n_lines
+    text_shape = (text_height, text_width)
+
+    # Instantiate grayscale image of correct shape (in pixels)
+    img = Image.new("L", (text_width, text_height), 0)
+    draw = ImageDraw.Draw(img)
+
+    # Draw text into this image
+    draw.text((0, 0), text, fill=1, font=font, align=align, direction=direction)
+
+    # Turn into mask-array
+    mask = np.array(img)
+    img = np.where(mask, intensity_text, intensity_background)
+
+    # Package as dict
+    stim = {
+        "img": img,
+        "text_mask": mask,
+        "text_shape": text_shape,
+    }
+
+    # Resolve resolution
+    if not shape or shape == (None, None):
+        shape = text_shape
+    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
+
+    # Pad
+    stim = pad_dict_to_visual_size(stim, visual_size=visual_size, ppd=ppd, pad_value=0.5)
+
+    # Output
+    return stim
+
+
+def overview(**kwargs):
+    """Generate example stimuli from this module
+
+    Returns
+    -------
+    stims : dict
+        dict with all stimuli containing individual stimulus dicts.
+    """
+    default_params = {
+        "visual_size": (10, 10),
+        "ppd": 32,
+    }
+    default_params.update(kwargs)
+
+    # fmt: off
+    stimuli = {
+        "text(), single line": text(text="hello world", **default_params),
+        "text(), multiline": text(text="hello\nworld", **default_params)
+    }
+    # fmt: on
+
+    return stimuli
+
+
+if __name__ == "__main__":
+    from stimupy.utils import plot_stimuli
+
+    stims = overview()
+    plot_stimuli(stims, mask=False, save=None)
diff --git a/stimupy/logos.py b/stimupy/logos.py
index 7d39de86..9457a1b2 100644
--- a/stimupy/logos.py
+++ b/stimupy/logos.py
@@ -112,7 +112,7 @@ def qrcode(ppd=128):
         embeded_image=im_logo,
     )
 
-    ## Make "eyes" circular
+    # Make "eyes" circular
     eye_size = 70  # default
     quiet_zone = 40  # default
 
diff --git a/stimupy/noises/binaries.py b/stimupy/noises/binaries.py
index 9c881eaa..0ad10ccf 100644
--- a/stimupy/noises/binaries.py
+++ b/stimupy/noises/binaries.py
@@ -71,7 +71,8 @@ def overview(**kwargs):
 
     # fmt: off
     stimuli = {
-        "binaries_binary": binary(**default_params),}
+        "binaries_binary": binary(**default_params),
+    }
     # fmt: on
 
     return stimuli
diff --git a/stimupy/noises/narrowbands.py b/stimupy/noises/narrowbands.py
index d8e9d733..9921f084 100644
--- a/stimupy/noises/narrowbands.py
+++ b/stimupy/noises/narrowbands.py
@@ -1,7 +1,7 @@
 import numpy as np
 
 from stimupy.noises import pseudo_white_spectrum
-from stimupy.utils import resolution, bandpass
+from stimupy.utils import bandpass, resolution
 from stimupy.utils.contrast_conversions import adapt_intensity_range
 
 __all__ = [
@@ -107,7 +107,8 @@ def overview(**kwargs):
     # fmt: off
     stimuli = {
         "narrowbands_narrowband3": narrowband(**default_params, center_frequency=3, bandwidth=1),
-        "narrowbands_narrowband9": narrowband(**default_params, center_frequency=9, bandwidth=1),}
+        "narrowbands_narrowband9": narrowband(**default_params, center_frequency=9, bandwidth=1),
+    }
     # fmt: on
 
     return stimuli
diff --git a/stimupy/noises/naturals.py b/stimupy/noises/naturals.py
index 89e53b60..bff07c18 100644
--- a/stimupy/noises/naturals.py
+++ b/stimupy/noises/naturals.py
@@ -183,7 +183,8 @@ def overview(**kwargs):
     stimuli = {
         "naturals_one-over-f": one_over_f(**default_params, exponent=0.5),
         "naturals_pink": pink(**default_params),
-        "naturals_brown": brown(**default_params),}
+        "naturals_brown": brown(**default_params),
+    }
     # fmt: on
 
     return stimuli
diff --git a/stimupy/noises/whites.py b/stimupy/noises/whites.py
index f58b6453..84120862 100644
--- a/stimupy/noises/whites.py
+++ b/stimupy/noises/whites.py
@@ -83,7 +83,8 @@ def overview(**kwargs):
 
     # fmt: off
     stimuli = {
-        "whites_white": white(**default_params),}
+        "whites_white": white(**default_params),
+    }
     # fmt: on
 
     return stimuli
diff --git a/stimupy/papers/RHS2007.py b/stimupy/papers/RHS2007.py
index 9fd8363f..2ed486c5 100644
--- a/stimupy/papers/RHS2007.py
+++ b/stimupy/papers/RHS2007.py
@@ -137,9 +137,9 @@ def WE_thick(ppd=PPD, pad=True):
     params = {
         "ppd": ppd,
         "frequency": 4.0 / width,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
         "intensity_target": v2,
-        "target_indices": (2, 5),
+        "target_indices": (3, 6),
         "target_heights": 4.0,
         "period": "even",
     }
@@ -198,9 +198,9 @@ def WE_thin_wide(ppd=PPD, pad=True):
     params = {
         "ppd": ppd,
         "frequency": 8.0 / width,
-        "intensity_bars": (v3, v1),
+        "intensity_bars": (v1, v3),
         "intensity_target": v2,
-        "target_indices": (3, 12),
+        "target_indices": (4, 13),
         "target_heights": 2.0,
         "period": "even",
     }
@@ -250,9 +250,9 @@ def WE_dual(ppd=PPD, pad=True):
     params = {
         "ppd": ppd,
         "frequency": 4.0 / width,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
         "intensity_target": v2,
-        "target_indices": (2, 5),
+        "target_indices": (3, 6),
         "target_heights": 2.0,
         "period": "even",
     }
@@ -320,10 +320,10 @@ def WE_anderson(ppd=PPD, pad=True):
     params = {
         "ppd": ppd,
         "frequency": 8.0 / width,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
         "intensity_target": v2,
-        "target_indices_top": (5,),
-        "target_indices_bottom": (10,),
+        "target_indices_top": (6,),
+        "target_indices_bottom": (11,),
         "target_center_offset": height / 10.0,
         "target_height": height / 5.0,
         "intensity_stripes": (v1, v3),
@@ -389,10 +389,10 @@ def WE_howe(ppd=PPD, pad=True):
     params = {
         "ppd": ppd,
         "frequency": 8.0 / width,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
         "intensity_target": v2,
-        "target_indices_top": (5,),
-        "target_indices_bottom": (10,),
+        "target_indices_top": (6,),
+        "target_indices_bottom": (11,),
         "target_center_offset": height / 5.0,
         "target_height": height / 5.0,
         "intensity_stripes": (v1, v3),
@@ -495,7 +495,7 @@ def WE_radial_thick_small(ppd=PPD, pad=True):
         "target_indices": (1, 8),
         "target_width": 4.0,
         "target_center": 4.5,
-        "intensity_segments": (v3, v1),
+        "intensity_segments": (v1, v3),
         "intensity_background": v2,
         "intensity_target": v2,
     }
@@ -550,7 +550,7 @@ def WE_radial_thick(ppd=PPD, pad=True):
         "target_indices": (1, 10),
         "target_width": 4.0,
         "target_center": 6.0,
-        "intensity_segments": (v3, v1),
+        "intensity_segments": (v1, v3),
         "intensity_background": v2,
         "intensity_target": v2,
     }
@@ -605,7 +605,7 @@ def WE_radial_thin_small(ppd=PPD, pad=True):
         "target_indices": (1, 14),
         "target_width": 2.0,
         "target_center": 4.0,
-        "intensity_segments": (v3, v1),
+        "intensity_segments": (v1, v3),
         "intensity_background": v2,
         "intensity_target": v2,
     }
@@ -660,7 +660,7 @@ def WE_radial_thin(ppd=PPD, pad=True):
         "target_indices": (1, 22),
         "target_width": 2.0,
         "target_center": 6.0,
-        "intensity_segments": (v3, v1),
+        "intensity_segments": (v1, v3),
         "intensity_background": v2,
         "intensity_target": v2,
     }
@@ -712,7 +712,7 @@ def WE_circular1(ppd=PPD, pad=True):
     params = {
         "ppd": ppd,
         "frequency": 8.0 / height,
-        "target_indices": 4,
+        "target_indices": 5,
         "intensity_background": v2,
         "intensity_target": v2,
         "clip": True,
@@ -720,13 +720,13 @@ def WE_circular1(ppd=PPD, pad=True):
 
     stim1 = stimupy.whites.radial(
         visual_size=(height, width),
-        intensity_rings=(v3, v1),
+        intensity_rings=(v1, v3),
         **params,
     )
 
     stim2 = stimupy.whites.radial(
         visual_size=(height, width),
-        intensity_rings=(v1, v3),
+        intensity_rings=(v3, v1),
         **params,
     )
     stim = stack_dicts(stim1, stim2)
@@ -773,7 +773,7 @@ def WE_circular05(ppd=PPD, pad=True):
     params = {
         "ppd": ppd,
         "frequency": 16.0 / height,
-        "target_indices": 10,
+        "target_indices": 11,
         "intensity_background": v2,
         "intensity_target": v2,
         "clip": True,
@@ -781,13 +781,13 @@ def WE_circular05(ppd=PPD, pad=True):
 
     stim1 = stimupy.whites.radial(
         visual_size=(height, width),
-        intensity_rings=(v3, v1),
+        intensity_rings=(v1, v3),
         **params,
     )
 
     stim2 = stimupy.whites.radial(
         visual_size=(height, width),
-        intensity_rings=(v1, v3),
+        intensity_rings=(v3, v1),
         **params,
     )
     stim = stack_dicts(stim1, stim2)
@@ -834,7 +834,7 @@ def WE_circular025(ppd=PPD, pad=True):
     params = {
         "ppd": ppd,
         "frequency": 32.0 / height,
-        "target_indices": 22,
+        "target_indices": 23,
         "intensity_background": v2,
         "intensity_target": v2,
         "clip": True,
@@ -842,13 +842,13 @@ def WE_circular025(ppd=PPD, pad=True):
 
     stim1 = stimupy.whites.radial(
         visual_size=(height, width),
-        intensity_rings=(v3, v1),
+        intensity_rings=(v1, v3),
         **params,
     )
 
     stim2 = stimupy.whites.radial(
         visual_size=(height, width),
-        intensity_rings=(v1, v3),
+        intensity_rings=(v3, v1),
         **params,
     )
     stim = stack_dicts(stim1, stim2)
@@ -955,23 +955,13 @@ def sbc_large(ppd=PPD, pad=True):
         Vision Research, 39, 4361-4377.
     """
 
-    params = {
-        "ppd": ppd,
-        "target_size": 3.0,
-        "intensity_target": v2,
-    }
-
-    stim1 = stimupy.sbcs.basic(
-        visual_size=(13.0, 15.5),
-        intensity_background=0.0,
-        **params,
-    )
-    stim2 = stimupy.sbcs.basic(
-        visual_size=(13.0, 15.5),
-        intensity_background=1.0,
-        **params,
+    stim = stimupy.sbcs.basic_two_sided(
+        visual_size=(13.0, 31.0),
+        ppd=ppd,
+        target_size=3.0,
+        intensity_target=v2,
+        intensity_background=(0.0, 1.0),
     )
-    stim = stack_dicts(stim1, stim2)
 
     if pad:
         stim = pad_dict_to_visual_size(stim, VISEXTENT, ppd, pad_value=v2)
@@ -980,13 +970,13 @@ def sbc_large(ppd=PPD, pad=True):
         "effect_strength": 11.35,
     }
 
-    params.update(
+    stim.update(
         visual_size=np.array(stim["img"].shape) / ppd,
         shape=stim["img"].shape,
         intensity_range=(v1, v3),
         experimental_data=experimental_data,
     )
-    return {**stim, **params}
+    return stim
 
 
 def sbc_small(ppd=PPD, pad=True):
@@ -1015,23 +1005,13 @@ def sbc_small(ppd=PPD, pad=True):
         Vision Research, 39, 4361-4377.
     """
 
-    params = {
-        "ppd": ppd,
-        "target_size": 1.0,
-        "intensity_target": v2,
-    }
-
-    stim1 = stimupy.sbcs.basic(
-        visual_size=(13.0, 15.5),
-        intensity_background=0.0,
-        **params,
-    )
-    stim2 = stimupy.sbcs.basic(
-        visual_size=(13.0, 15.5),
-        intensity_background=1.0,
-        **params,
+    stim = stimupy.sbcs.basic_two_sided(
+        visual_size=(13.0, 31.0),
+        ppd=ppd,
+        target_size=1.0,
+        intensity_background=(0.0, 1.0),
+        intensity_target=v2,
     )
-    stim = stack_dicts(stim1, stim2)
 
     if pad:
         stim = pad_dict_to_visual_size(stim, VISEXTENT, ppd, pad_value=v2)
@@ -1040,13 +1020,13 @@ def sbc_small(ppd=PPD, pad=True):
         "effect_strength": 19.78,
     }
 
-    params.update(
+    stim.update(
         visual_size=np.array(stim["img"].shape) / ppd,
         shape=stim["img"].shape,
         intensity_range=(v1, v3),
         experimental_data=experimental_data,
     )
-    return {**stim, **params}
+    return stim
 
 
 def todorovic_equal(ppd=PPD, pad=True):
@@ -1645,38 +1625,22 @@ def todorovic_benary1_2(ppd=PPD, pad=True):
         Perception, 26, 379-395.
     """
 
-    params = {
-        "ppd": ppd,
-        "L_width": 2.5,
-        "target_size": 2.5,
-        "target_type": ("t", "t"),
-        "target_rotation": (0.0, 180.0),
-        "target_x": (2.5, 26.0),
-        "target_y": (4.0, 6.5),
-        "intensity_background": v3,
-        "intensity_cross": v1,
-        "intensity_target": v2,
-    }
+    # Generate full stim
+    stim = todorovic_benary1_2_3_4(ppd=ppd, pad=pad)
 
-    stim = stimupy.benarys.todorovic_generalized(
-        visual_size=(13.0, 31.0),
-        **params,
-    )
+    # Remove unused targets
+    stim["img"] = np.where(stim["target_mask"] == 3, v3, stim["img"])
+    stim["img"] = np.where(stim["target_mask"] == 4, v1, stim["img"])
+    stim["target_mask"] = np.where(stim["target_mask"] > 2, 0, stim["target_mask"])
+    stim["target_type"] = stim["target_type"][:2]
+    stim["target_rotation"] = stim["target_rotation"][:2]
+    stim["target_x"] = stim["target_x"][:2]
+    stim["target_y"] = stim["target_y"][:2]
 
-    if pad:
-        stim = pad_dict_to_visual_size(stim, VISEXTENT, ppd, pad_value=v2)
-
-    experimental_data = {
-        "effect_strength": 11.96,
-    }
+    # Add experimental data
+    stim["effect_strength"] = 11.96
 
-    params.update(
-        visual_size=np.array(stim["img"].shape) / ppd,
-        shape=stim["img"].shape,
-        intensity_range=(v1, v3),
-        experimental_data=experimental_data,
-    )
-    return {**stim, **params}
+    return stim
 
 
 def todorovic_benary3_4(ppd=PPD, pad=True):
@@ -1707,38 +1671,25 @@ def todorovic_benary3_4(ppd=PPD, pad=True):
         Perception, 26, 379-395.
     """
 
-    params = {
-        "ppd": ppd,
-        "L_width": 2.5,
-        "target_size": 2.5,
-        "target_type": ("t", "t"),
-        "target_rotation": (45.0, 225.0),
-        "target_x": (9.5, 18.0),
-        "target_y": (6.5, 6.5 - np.sqrt(12.5) / 2.0 + 1 / ppd),
-        "intensity_background": v3,
-        "intensity_cross": v1,
-        "intensity_target": v2,
-    }
+    # Generate full stim
+    stim = todorovic_benary1_2_3_4(ppd=ppd, pad=pad)
 
-    stim = stimupy.benarys.todorovic_generalized(
-        visual_size=(13.0, 31.0),
-        **params,
+    # Remove unused targets
+    stim["img"] = np.where(stim["target_mask"] == 1, v3, stim["img"])
+    stim["img"] = np.where(stim["target_mask"] == 2, v1, stim["img"])
+    stim["target_mask"] = np.where(stim["target_mask"] < 3, 0, stim["target_mask"])
+    stim["target_mask"] = np.where(
+        stim["target_mask"], stim["target_mask"] - 2, stim["target_mask"]
     )
+    stim["target_type"] = stim["target_type"][2:]
+    stim["target_rotation"] = stim["target_rotation"][2:]
+    stim["target_x"] = stim["target_x"][2:]
+    stim["target_y"] = stim["target_y"][2:]
 
-    if pad:
-        stim = pad_dict_to_visual_size(stim, VISEXTENT, ppd, pad_value=v2)
-
-    experimental_data = {
-        "effect_strength": 9.55,
-    }
+    # Add experimental data
+    stim["effect_strength"] = 9.55
 
-    params.update(
-        visual_size=np.array(stim["img"].shape) / ppd,
-        shape=stim["img"].shape,
-        intensity_range=(v1, v3),
-        experimental_data=experimental_data,
-    )
-    return {**stim, **params}
+    return stim
 
 
 def todorovic_benary1_2_3_4(ppd=PPD, pad=True):
diff --git a/stimupy/papers/domijan2015.py b/stimupy/papers/domijan2015.py
index b4e8a8ef..7afb4784 100644
--- a/stimupy/papers/domijan2015.py
+++ b/stimupy/papers/domijan2015.py
@@ -362,18 +362,18 @@ def grating(visual_size=VSIZES["grating"], ppd=PPD, shape=SHAPES["grating"]):
         "visual_size": single_vissize,
         "ppd": ppd,
         "n_bars": 9,
-        "target_indices": (4,),
+        "target_indices": (5,),
         "bar_width": 1.0 * visual_resize,
     }
 
     stim1 = stimupy.waves.square_linear(
         **params,
-        intensity_bars=(v3, v1),
+        intensity_bars=(v1, v3),
         intensity_target=v2,
     )
     stim2 = stimupy.waves.square_linear(
         **params,
-        intensity_bars=(v1, v3),
+        intensity_bars=(v3, v1),
         intensity_target=v2,
     )
 
@@ -564,30 +564,17 @@ def simultaneous_brightness_contrast(
     shape, visual_size, ppd, visual_resize = resolve(
         shape, visual_size, ppd, VSIZES["simultaneous_brightness_contrast"]
     )
-    ppd = ppd[0]
 
-    params = {
-        "visual_size": visual_size[0],
-        "ppd": ppd,
-        "target_size": (2.1 * visual_resize, 2.1 * visual_resize),
-        "target_position": (3.8 * visual_resize, 3.8 * visual_resize),
-    }
-
-    stim1 = stimupy.sbcs.generalized(
-        **params,
-        intensity_background=v3,
+    stim = stimupy.sbcs.generalized_two_sided(
+        visual_size=visual_size,
+        ppd=ppd,
+        target_size=(2.1 * visual_resize, 2.1 * visual_resize),
+        target_position=(3.8 * visual_resize, 3.8 * visual_resize),
+        intensity_background=(v3, v1),
         intensity_target=v2,
     )
-    stim2 = stimupy.sbcs.generalized(
-        **params,
-        intensity_background=v1,
-        intensity_target=v2,
-    )
-
-    # Stacking
-    stim = stack_dicts(stim1, stim2)
 
-    params.update(
+    stim.update(
         original_shape=SHAPES["simultaneous_brightness_contrast"],
         original_ppd=PPD,
         original_visual_size=VSIZES["simultaneous_brightness_contrast"],
@@ -596,7 +583,7 @@ def simultaneous_brightness_contrast(
         visual_size=resolution.visual_size_from_shape_ppd(stim["img"].shape, ppd),
         shape=stim["img"].shape,
     )
-    return {**stim, **params}
+    return stim
 
 
 def white(visual_size=VSIZES["white"], ppd=PPD, pad=PAD, shape=SHAPES["white"]):
@@ -636,14 +623,14 @@ def white(visual_size=VSIZES["white"], ppd=PPD, pad=PAD, shape=SHAPES["white"]):
         "visual_size": visual_size,
         "ppd": ppd,
         "frequency": 4.0 / visual_size[1],
-        "target_indices": (2, 5),
+        "target_indices": (3, 6),
         "target_heights": 2.1 * visual_resize,
         "period": "even",
     }
 
     stim = stimupy.whites.white(
         **params,
-        intensity_bars=(v1, v3),
+        intensity_bars=(v3, v1),
         intensity_target=v2,
     )
 
@@ -1080,8 +1067,8 @@ def white_yazdanbakhsh(
         "visual_size": visual_size,
         "ppd": ppd,
         "frequency": 4.0 / visual_size[1],
-        "target_indices_top": (2,),
-        "target_indices_bottom": (5,),
+        "target_indices_top": (3,),
+        "target_indices_bottom": (6,),
         "target_center_offset": 0.0,
         "target_heights": visual_size[0] / 4.0,
         "gap_size": visual_size[0] / 10.0,
@@ -1090,7 +1077,7 @@ def white_yazdanbakhsh(
 
     stim = stimupy.whites.yazdanbakhsh(
         **params,
-        intensity_bars=(v1, v3),
+        intensity_bars=(v3, v1),
         intensity_target=v2,
     )
 
@@ -1165,8 +1152,8 @@ def white_anderson(
         "visual_size": visual_size,
         "ppd": ppd,
         "frequency": 5.0 / visual_size[1],
-        "target_indices_top": (2,),
-        "target_indices_bottom": (7,),
+        "target_indices_top": (3,),
+        "target_indices_bottom": (8,),
         "target_center_offset": visual_size[0] / 10.0,
         "target_height": visual_size[0] / 5.0,
         "stripe_center_offset": visual_size[0] / 5.0,
@@ -1176,7 +1163,7 @@ def white_anderson(
 
     stim = stimupy.whites.anderson(
         **params,
-        intensity_bars=(v3, v1),
+        intensity_bars=(v1, v3),
         intensity_target=v2,
         intensity_stripes=(v1, v3),
     )
@@ -1249,8 +1236,8 @@ def white_howe(visual_size=VSIZES["white_howe"], ppd=PPD, shape=SHAPES["white_ho
         "visual_size": visual_size,
         "ppd": ppd,
         "frequency": 5.0 / visual_size[1],
-        "target_indices_top": (2,),
-        "target_indices_bottom": (7,),
+        "target_indices_top": (3,),
+        "target_indices_bottom": (8,),
         "target_center_offset": visual_size[0] / 5.0,
         "target_height": visual_size[0] / 5.0,
         "period": "even",
@@ -1258,7 +1245,7 @@ def white_howe(visual_size=VSIZES["white_howe"], ppd=PPD, shape=SHAPES["white_ho
 
     stim = stimupy.whites.howe(
         **params,
-        intensity_bars=(v3, v1),
+        intensity_bars=(v1, v3),
         intensity_target=v2,
         intensity_stripes=(v1, v3),
     )
diff --git a/stimupy/papers/modelfest.py b/stimupy/papers/modelfest.py
index d39b95ad..ac5db9a8 100644
--- a/stimupy/papers/modelfest.py
+++ b/stimupy/papers/modelfest.py
@@ -44,7 +44,7 @@
 from stimupy.noises.binaries import binary as binary_noise
 from stimupy.stimuli.gabors import gabor
 from stimupy.stimuli.plaids import gabors as plaid
-from stimupy.utils import pad_dict_to_shape, resize_dict, roll_dict, stack_dicts
+from stimupy.utils import pad_dict_to_shape, resize_dict, resolution, roll_dict, stack_dicts
 
 __all__ = [
     "GaborPatch1",
@@ -95,6 +95,7 @@
 # Default values:
 PPD = 120  # pixel size of 0.5 arcmin
 PPD2 = 60  # pixel size of 1 arcmin (pixel replication)
+VISUAL_SIZE = resolution.visual_size_from_shape_ppd(shape=256, ppd=PPD)
 
 mean_lum = 0.5
 
@@ -151,7 +152,7 @@ def GaborPatch1(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 1.12,
         "sigma": 0.5,
@@ -200,7 +201,7 @@ def GaborPatch2(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 2,
         "sigma": 0.5,
@@ -249,7 +250,7 @@ def GaborPatch3(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 2.83,
         "sigma": 0.5,
@@ -298,7 +299,7 @@ def GaborPatch4(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": 0.5,
@@ -347,7 +348,7 @@ def GaborPatch5(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 5.66,
         "sigma": 0.5,
@@ -396,7 +397,7 @@ def GaborPatch6(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 8,
         "sigma": 0.5,
@@ -445,7 +446,7 @@ def GaborPatch7(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 11.3,
         "sigma": 0.5,
@@ -494,7 +495,7 @@ def GaborPatch8(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 16,
         "sigma": 0.5,
@@ -543,7 +544,7 @@ def GaborPatch9(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 22.6,
         "sigma": 0.5,
@@ -592,7 +593,7 @@ def GaborPatch10(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 30,
         "sigma": 0.5,
@@ -641,7 +642,7 @@ def GaborPatch11(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 2,
         "sigma": 0.28,
@@ -690,7 +691,7 @@ def GaborPatch12(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": 0.14,
@@ -739,7 +740,7 @@ def GaborPatch13(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 8,
         "sigma": 0.07,
@@ -788,7 +789,7 @@ def GaborPatch14(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 16,
         "sigma": 0.035,
@@ -837,7 +838,7 @@ def ElongatedGabor15(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": (0.28, 0.5),
@@ -886,7 +887,7 @@ def ElongatedGabor16(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 8,
         "sigma": (0.14, 0.5),
@@ -935,7 +936,7 @@ def ElongatedGabor17(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 16,
         "sigma": (0.07, 0.5),
@@ -984,7 +985,7 @@ def Baguette18(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": (0.14, 0.28),
@@ -1033,7 +1034,7 @@ def Baguette19(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": (0.14, 0.5),
@@ -1082,7 +1083,7 @@ def Baguette20(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": (0.28, 0.14),
@@ -1131,7 +1132,7 @@ def Baguette21(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": (0.5, 0.14),
@@ -1180,7 +1181,7 @@ def Subthreshold22(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "sigma": 0.5,
         "rotation": -90,
@@ -1235,7 +1236,7 @@ def Subthreshold23(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "sigma": 0.5,
         "rotation": -90,
@@ -1290,7 +1291,7 @@ def Subthreshold24(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "sigma": 0.5,
         "rotation": -90,
@@ -1345,7 +1346,7 @@ def Subthreshold25(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "sigma": 0.5,
         "rotation": -90,
@@ -1399,7 +1400,7 @@ def Gaussians26(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "sigma": 0.5,
         "intensity_max": 1,
@@ -1446,7 +1447,7 @@ def Gaussians27(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "sigma": 8.43 / 60,
         "intensity_max": 1,
@@ -1493,7 +1494,7 @@ def Gaussians28(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "sigma": 2.106 / 60,
         "intensity_max": 1,
@@ -1540,7 +1541,7 @@ def Gaussians29(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "sigma": 1.05 / 60,
         "intensity_max": 1,
@@ -1586,7 +1587,7 @@ def Edge30(ppd=PPD):
         https://doi.org/10.1117/12.348473
     """
 
-    stim = gaussian_edge(visual_size=256 / PPD, ppd=ppd, rotation=-90, sigma=0.5)
+    stim = gaussian_edge(visual_size=VISUAL_SIZE, ppd=ppd, rotation=-90, sigma=0.5)
 
     v = 117
     experimental_data = {
@@ -1626,14 +1627,14 @@ def Line31(ppd=PPD):
     """
 
     stim = lines.line(
-        visual_size=256 / PPD,
+        visual_size=VISUAL_SIZE,
         ppd=ppd,
-        line_length=256 / PPD,
-        line_width=0.5 / 60,
+        line_length=VISUAL_SIZE.width,
+        line_width=resolution.visual_angle_from_length_ppd(length=0.5, ppd=PPD2),
         rotation=90,
         intensity_background=0.5,
     )
-    window = gaussians.gaussian(visual_size=256 / PPD, ppd=ppd, sigma=0.5)
+    window = gaussians.gaussian(visual_size=VISUAL_SIZE, ppd=ppd, sigma=0.5)
 
     img = (stim["img"] - 0.5) * window["img"] + 0.5
     stim["img"] = img
@@ -1678,15 +1679,15 @@ def Dipole32(ppd=PPD):
     """
 
     stim = lines.dipole(
-        visual_size=256 / PPD,
+        visual_size=VISUAL_SIZE,
         ppd=ppd,
-        line_length=256 / PPD,
-        line_width=1 / PPD,
-        line_gap=2 / PPD,
+        line_length=VISUAL_SIZE.width,
+        line_width=resolution.visual_angle_from_length_ppd(length=1, ppd=PPD),
+        line_gap=resolution.visual_angle_from_length_ppd(length=2, ppd=PPD),
         rotation=270,
     )
     stim = roll_dict(stim, -1, axes=0)
-    window = gaussians.gaussian(visual_size=256 / PPD, ppd=ppd, sigma=0.5)
+    window = gaussians.gaussian(visual_size=VISUAL_SIZE, ppd=ppd, sigma=0.5)
 
     img = (stim["img"] - 0.5) * window["img"] + 0.5
     stim["img"] = img
@@ -1731,7 +1732,7 @@ def GaborString33(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD / 6,
+        "visual_size": VISUAL_SIZE.width / 6,
         "ppd": ppd,
         "frequency": 8,
         "sigma": 0.07,
@@ -1793,7 +1794,7 @@ def GaborString34(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD / 6,
+        "visual_size": VISUAL_SIZE.width / 6,
         "ppd": ppd,
         "frequency": 8,
         "sigma": 0.07,
@@ -1854,9 +1855,9 @@ def Noise35_random(ppd=PPD):
         https://doi.org/10.1117/12.348473
     """
 
-    stim = binary_noise(visual_size=256 / PPD, ppd=ppd / 2, rms_contrast=1)
+    stim = binary_noise(visual_size=VISUAL_SIZE, ppd=ppd / 2, rms_contrast=1)
     stim = resize_dict(stim, (2, 2))
-    window = gaussians.gaussian(visual_size=256 / PPD, ppd=ppd, sigma=0.5)
+    window = gaussians.gaussian(visual_size=VISUAL_SIZE, ppd=ppd, sigma=0.5)
 
     img = stim["img"] * window["img"]
     stim["img"] = img / 2 + 0.5
@@ -1894,7 +1895,7 @@ def Noise35(ppd=PPD):
 
     stim = {
         "img": img,
-        "visual_size": (256 / PPD, 256 / PPD),
+        "visual_size": VISUAL_SIZE,
         "shape": img.shape,
         "ppd": PPD,
         "intensity_range": (img.min(), img.max()),
@@ -1940,7 +1941,7 @@ def Orientation36(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": 0.14,
@@ -1990,7 +1991,7 @@ def Orientation37(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": 0.14,
@@ -2040,7 +2041,7 @@ def Plaids38(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": 0.14,
@@ -2093,7 +2094,7 @@ def Plaids39(ppd=PPD):
     """
 
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 4,
         "sigma": 0.14,
@@ -2143,7 +2144,7 @@ def Disk40(ppd=PPD):
     """
 
     stim = shapes.disc(
-        visual_size=256 / PPD, ppd=ppd, radius=0.125, origin="center", intensity_background=0.5
+        visual_size=VISUAL_SIZE, ppd=ppd, radius=0.1249, origin="center", intensity_background=0.5
     )
     stim = roll_dict(stim, (-2, -2), axes=(0, 1))
 
@@ -2184,8 +2185,8 @@ def Bessel41(ppd=PPD):
         https://doi.org/10.1117/12.348473
     """
 
-    stim = bessel(visual_size=256 / PPD, ppd=ppd, frequency=4, origin="center")
-    window = gaussians.gaussian(visual_size=256 / PPD, ppd=ppd, sigma=0.5)
+    stim = bessel(visual_size=VISUAL_SIZE, ppd=ppd, frequency=4, origin="center")
+    window = gaussians.gaussian(visual_size=VISUAL_SIZE, ppd=ppd, sigma=0.5)
     stim["sigma"] = window["sigma"]
 
     # Apply Gaussian windows to Bessel
@@ -2231,7 +2232,7 @@ def Checkerboard42(ppd=PPD):
         https://doi.org/10.1117/12.348473
     """
     params = {
-        "visual_size": 256 / PPD,
+        "visual_size": VISUAL_SIZE,
         "ppd": ppd,
         "frequency": 2.81,
         "rotation": -45,
@@ -2240,7 +2241,7 @@ def Checkerboard42(ppd=PPD):
     }
 
     stim = checkerboards.checkerboard(**params)
-    window = gaussians.gaussian(visual_size=256 / PPD, ppd=ppd, sigma=0.5)
+    window = gaussians.gaussian(visual_size=VISUAL_SIZE, ppd=ppd, sigma=0.5)
 
     img = (stim["img"] - 0.5) * window["img"] + 0.5
     stim["img"] = img
@@ -2284,11 +2285,11 @@ def NaturalScene43(ppd=PPD):
 
     # Read natural image from Modelfest
     img = read_tif(Path(__file__).parents[0] / "modelfest_natural_scene.tif")
-    img = img / 255
+    img /= 255
 
     stim = {
         "img": img,
-        "visual_size": (256 / PPD, 256 / PPD),
+        "visual_size": VISUAL_SIZE,
         "shape": img.shape,
         "ppd": PPD,
         "intensity_range": (img.min(), img.max()),
@@ -2313,41 +2314,8 @@ def read_tif(filename):
     return img
 
 
-def compare(o1, s1, filename):
-    import matplotlib.pyplot as plt
-
-    o1 = o1 / 255
-    s1 = s1 / 1
-    vmin, vmax = 0, 1
-
-    plt.figure(figsize=(20, 6))
-    plt.subplot(141), plt.imshow(o1, vmin=vmin, vmax=vmax), plt.title("Original")
-    plt.subplot(142), plt.imshow(s1, vmin=vmin, vmax=vmax), plt.title("Our implementation")
-    plt.subplot(143), plt.imshow(o1 - s1, vmin=-vmax, vmax=vmax), plt.title(
-        "Difference"
-    ), plt.colorbar()
-    plt.subplot(144), plt.plot(o1[:, 128], label="Original")
-    plt.plot(s1[:, 128], label="Our implementation"), plt.legend()
-    plt.savefig("./comparisons/" + filename)
-    plt.close()
-
-
-def compare_all():
-    import os
-
-    if not os.path.exists("./comparisons/"):
-        os.makedirs("./comparisons/")
-
-    for stim_name in __all__:
-        func = globals()[stim_name]
-        o1 = read_tif(str(Path(__file__).parents[0]) + "/modelfest/" + stim_name + ".tif")
-        s1 = func()["img"]
-        compare(o1, s1, stim_name + ".png")
-
-
 if __name__ == "__main__":
     from stimupy.utils import plot_stimuli
 
     stims = gen_all(skip=True)
     plot_stimuli(stims, mask=False, units="visual_size")
-    # compare_all()
diff --git a/stimupy/papers/murray2020.py b/stimupy/papers/murray2020.py
index 230d2dbb..0de8beb0 100644
--- a/stimupy/papers/murray2020.py
+++ b/stimupy/papers/murray2020.py
@@ -731,11 +731,11 @@ def white(ppd=PPD):
         "ppd": ppd,
         "visual_size": (8.0, 8.0),
         "frequency": 4.0 / 8.0,
-        "target_indices_top": (1, 3, 5),
-        "target_indices_bottom": (2, 4, 6),
+        "target_indices_top": (2, 4, 6),
+        "target_indices_bottom": (3, 5, 7),
         "target_center_offset": 2,
         "target_heights": 2,
-        "intensity_bars": (70, 17.5),
+        "intensity_bars": (17.5, 70),
         "intensity_target": 35.0,
         "period": "even",
     }
diff --git a/stimupy/papers/white1981.py b/stimupy/papers/white1981.py
index 7bc1d435..33bb4048 100644
--- a/stimupy/papers/white1981.py
+++ b/stimupy/papers/white1981.py
@@ -202,9 +202,9 @@ def grating_white_white(ppd=PPD):
         "bar_width": bar_width,
         "rotation": ROTATION,
         "intensity_background": v3,
-        "intensity_bars": (v3, v2),
+        "intensity_bars": (v2, v3),
         "intensity_target": v2,
-        "target_indices": (1, 3, 5, 7),
+        "target_indices": (2, 4, 6, 8),
     }
     stim = stimupy.gratings.on_uniform(**params)
     stim["target_mask"] = np.where(stim["target_mask"] != 0, 1, 0)
@@ -245,9 +245,9 @@ def grating_white_black(ppd=PPD):
         "bar_width": bar_width,
         "rotation": ROTATION,
         "intensity_background": v1,
-        "intensity_bars": (v3, v2),
+        "intensity_bars": (v2, v3),
         "intensity_target": v2,
-        "target_indices": (1, 3, 5, 7),
+        "target_indices": (2, 4, 6, 8),
     }
     stim = stimupy.gratings.on_uniform(**params)
     stim["target_mask"] = np.where(stim["target_mask"] != 0, 1, 0)
@@ -288,9 +288,9 @@ def grating_black_white(ppd=PPD):
         "bar_width": bar_width,
         "rotation": ROTATION,
         "intensity_background": v3,
-        "intensity_bars": (v1, v2),
+        "intensity_bars": (v2, v1),
         "intensity_target": v2,
-        "target_indices": (1, 3, 5, 7),
+        "target_indices": (2, 4, 6, 8),
     }
     stim = stimupy.gratings.on_uniform(**params)
     stim["target_mask"] = np.where(stim["target_mask"] != 0, 1, 0)
@@ -331,9 +331,9 @@ def grating_black_black(ppd=PPD):
         "bar_width": bar_width,
         "rotation": ROTATION,
         "intensity_background": v1,
-        "intensity_bars": (v1, v2),
+        "intensity_bars": (v2, v1),
         "intensity_target": v2,
-        "target_indices": (1, 3, 5, 7),
+        "target_indices": (2, 4, 6, 8),
     }
     stim = stimupy.gratings.on_uniform(**params)
     stim["target_mask"] = np.where(stim["target_mask"] != 0, 1, 0)
@@ -372,16 +372,16 @@ def grating_white_in(ppd=PPD):
         "n_bars": 9,
         "ppd": ppd,
         "rotation": ROTATION,
-        "intensity_bars": (v3, v1),
+        "intensity_bars": (v1, v3),
         "intensity_target": v2,
-        "target_indices": (1, 3, 5, 7),
+        "target_indices": (2, 4, 6, 8),
     }
     large_params = {
         "bar_width": bar_width,
         "visual_size": 51 * bar_width,
         "ppd": ppd,
         "rotation": ROTATION,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
     }
     stim = stimupy.gratings.on_grating(
         small_grating_params=small_params,
@@ -424,16 +424,16 @@ def grating_black_in(ppd=PPD):
         "n_bars": 9,
         "ppd": ppd,
         "rotation": ROTATION,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
         "intensity_target": v2,
-        "target_indices": (1, 3, 5, 7),
+        "target_indices": (2, 4, 6, 8),
     }
     large_params = {
         "bar_width": bar_width,
         "visual_size": 51 * bar_width,
         "ppd": ppd,
         "rotation": ROTATION,
-        "intensity_bars": (v3, v1),
+        "intensity_bars": (v1, v3),
     }
     stim = stimupy.gratings.on_grating(
         small_grating_params=small_params,
@@ -476,16 +476,16 @@ def grating_white_out(ppd=PPD):
         "visual_size": (50 * bar_width, 9 * bar_width),
         "ppd": ppd,
         "rotation": -ROTATION,
-        "intensity_bars": (v3, v1),
+        "intensity_bars": (v1, v3),
         "intensity_target": v2,
-        "target_indices": (23, 25, 27, 29),
+        "target_indices": (24, 26, 28, 30),
     }
     large_params = {
         "bar_width": bar_width,
         "visual_size": 50 * bar_width,
         "ppd": ppd,
         "rotation": -ROTATION,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
     }
     stim = stimupy.gratings.on_grating(
         small_grating_params=small_params,
@@ -528,16 +528,16 @@ def grating_black_out(ppd=PPD):
         "visual_size": (50 * bar_width, 9 * bar_width),
         "ppd": ppd,
         "rotation": -ROTATION,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
         "intensity_target": v2,
-        "target_indices": (23, 25, 27, 29),
+        "target_indices": (24, 26, 28, 30),
     }
     large_params = {
         "bar_width": bar_width,
         "visual_size": 50 * bar_width,
         "ppd": ppd,
         "rotation": -ROTATION,
-        "intensity_bars": (v3, v1),
+        "intensity_bars": (v1, v3),
     }
     stim = stimupy.gratings.on_grating(
         small_grating_params=small_params,
@@ -580,16 +580,16 @@ def grating_white_orthogonal(ppd=PPD):
         "bar_width": bar_width,
         "ppd": ppd,
         "rotation": ROTATION + 90,
-        "intensity_bars": (v3, v2),
+        "intensity_bars": (v2, v3),
         "intensity_target": v2,
-        "target_indices": (1, 3, 5, 7),
+        "target_indices": (2, 4, 6, 8),
     }
     large_params = {
         "bar_width": bar_width,
         "visual_size": 51 * bar_width,
         "ppd": ppd,
         "rotation": ROTATION,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
     }
     stim = stimupy.gratings.on_grating_masked(
         small_grating_params=small_params,
@@ -634,16 +634,16 @@ def grating_black_orthogonal(ppd=PPD):
         "bar_width": bar_width,
         "ppd": ppd,
         "rotation": ROTATION + 90,
-        "intensity_bars": (v1, v3),
+        "intensity_bars": (v3, v1),
         "intensity_target": v2,
-        "target_indices": (1, 3, 5, 7),
+        "target_indices": (2, 4, 6, 8),
     }
     large_params = {
         "bar_width": bar_width,
         "visual_size": 51 * bar_width,
         "ppd": ppd,
         "rotation": ROTATION,
-        "intensity_bars": (v3, v1),
+        "intensity_bars": (v1, v3),
     }
     stim = stimupy.gratings.on_grating_masked(
         small_grating_params=small_params,
diff --git a/stimupy/papers/white1985.py b/stimupy/papers/white1985.py
index 51cbae7b..eb600ea0 100644
--- a/stimupy/papers/white1985.py
+++ b/stimupy/papers/white1985.py
@@ -57,7 +57,7 @@
 START_PHASE2 = 0
 
 COMMON_PARAMS = {
-    "intensity_bars": (v1, v3),
+    "intensity_bars": (v3, v1),
     "intensity_target": v2,
     "distance_metric": "vertical",
 }
diff --git a/stimupy/stimuli/__init__.py b/stimupy/stimuli/__init__.py
index 22787463..e3925e7c 100644
--- a/stimupy/stimuli/__init__.py
+++ b/stimupy/stimuli/__init__.py
@@ -1,6 +1,13 @@
+import itertools
+
+import numpy as np
+
+from stimupy.components import draw_regions
 from stimupy.stimuli import *
 
 __all__ = [
+    "mask_targets",
+    "place_targets",
     "overview",
     "plot_overview",
     "benarys",
@@ -27,6 +34,103 @@
 ]
 
 
+def mask_targets(element_mask, target_indices):
+    """Indicate elements as targets
+
+    Creates a new target_mask from a mask of elements (e.g., grating bars, rings, frames, etc.),
+    by indexing these elements.
+
+    Parameters
+    ----------
+    element_mask : numpy.ndarray
+        mask with integer values for different elements / regions in a stimulus
+    target_indices : Sequence[int] or int
+        index or indices of elements to be designated as targets.
+        Index 0 should always refer to background region.
+        Indices can be negative, which results in "counting backwards"
+        from the highest index in element_mask.
+
+    Returns
+    -------
+    numpy.ndarray
+        target mask, with integer values indicating target regions,
+        in order that they appear in target_indices
+
+    Raises
+    ------
+    ValueError
+        if a target_idx is greater than any value in the element_mask
+    """
+    if target_indices is None:
+        target_indices = ()
+    if isinstance(target_indices, (int, float)):
+        target_indices = [
+            target_indices,
+        ]
+
+    target_mask = np.zeros_like(element_mask)
+    for target_idx, element_idx in enumerate(target_indices):
+        if element_idx < 0:
+            element_idx = int(element_mask.max()) + element_idx
+
+        if element_idx > element_mask.max():
+            raise ValueError("target_idx is outside stimulus")
+        target_mask = np.where(element_mask == element_idx, target_idx + 1, target_mask)
+
+    return target_mask
+
+
+def place_targets(stim, element_mask_key, target_indices, intensity_target=0.5):
+    """Place targets in stimulus
+
+    Turns image regions/elements defined by element_mask_key
+    and indicated by target_indices, into targets.
+    Targets are defined in a new target_mask, and drawn into image with intensity_target.
+
+    Parameters
+    ----------
+    stim : dict[str, Any]
+        stimulus dictionary, with at least an "img" key, and mask indicated by element_mask_key
+    element_mask_key : str
+        key of the mask in stim-dict indicating image "elements"/regions
+    target_indices : Sequence[int] or int
+        index or indices of elements to be designated as targets
+    intensity_target : float, optional
+        intensity value for target, by default 0.5
+
+    Returns
+    -------
+    dict[str, Any]
+        dict with the stimulus (key: "img") with targets placed,
+        mask with integer index for the target (key: "target_mask")
+
+    See also
+    --------
+        mask_targets, draw_regions
+    """
+    stim["target_mask"] = mask_targets(
+        element_mask=stim[element_mask_key], target_indices=target_indices
+    )
+
+    if isinstance(intensity_target, (int, float)):
+        intensity_target = [
+            intensity_target,
+        ]
+    intensity_target = itertools.cycle(intensity_target)
+
+    stim["img"] = np.where(
+        stim["target_mask"],
+        draw_regions(
+            mask=stim["target_mask"], intensities=intensity_target, intensity_background=0.0
+        ),
+        stim["img"],
+    )
+    stim["target_indices"] = target_indices
+    stim["intensity_target"] = intensity_target
+
+    return stim
+
+
 def overview(skip=False):
     """Generate example stimuli from this module
 
@@ -37,7 +141,7 @@ def overview(skip=False):
     """
     stimuli = {}
     for stimmodule_name in __all__:
-        if stimmodule_name in ["overview", "plot_overview"]:
+        if stimmodule_name in ["overview", "plot_overview", "mask_targets", "place_targets"]:
             continue
 
         print(f"Generating stimuli from {stimmodule_name}")
diff --git a/stimupy/stimuli/bullseyes.py b/stimupy/stimuli/bullseyes.py
index 1dd99481..7d2a6ee7 100644
--- a/stimupy/stimuli/bullseyes.py
+++ b/stimupy/stimuli/bullseyes.py
@@ -1,7 +1,9 @@
 from stimupy.stimuli import rings
+from stimupy.utils import make_two_sided
 
 __all__ = [
     "circular",
+    "circular_generalized",
     "circular_two_sided",
     "rectangular",
     "rectangular_generalized",
@@ -18,7 +20,7 @@ def circular(
     ring_width=None,
     phase_shift=0,
     intensity_target=0.5,
-    intensity_rings=(1.0, 0.0),
+    intensity_rings=(0.0, 1.0),
     intensity_background=0.5,
     origin="mean",
     clip=True,
@@ -105,43 +107,24 @@ def circular(
         intensity_rings=intensity_rings,
         intensity_background=intensity_background,
         intensity_target=intensity_target,
-        target_indices=0,
+        target_indices=1,
         origin=origin,
         clip=clip,
     )
     return stim
 
 
-def circular_two_sided(
+def circular_generalized(
     visual_size=None,
     ppd=None,
     shape=None,
-    frequency=None,
-    n_rings=None,
-    ring_width=None,
-    phase_shift=0,
-    intensity_target=0.5,
-    intensity_rings=(1.0, 0.0),
+    radii=None,
+    intensity_rings=(0.0, 1.0),
     intensity_background=0.5,
+    intensity_target=0.5,
     origin="mean",
 ):
-    """Two-sided circular Bullseye stimulus
-
-    Circular grating, where the target is the central disc.
-    Alias for circular_white(target_indices=0,...)
-
-    Specification of the number of rings, and their width can be done in two ways:
-    a ring_width (in degrees) and n_rings, and/or by specifying the spatial frequency
-    of a circular grating (in cycles per degree)
-
-    The total shape (in pixels) and visual size (in degrees) has to match the
-    specification of the rings and their widths.
-    Thus, not all 6 parameters (visual_size, ppd, shape, frequency, ring_width, n_rings)
-    have to be specified, as long as both the resolution, and the distribution of rings,
-    can be resolved.
-
-    Note: all rings in a grating have the same width -- if more control is required
-    see disc_and_rings
+    """Draw sequential set of circular rings with specified radii, with central target
 
     Parameters
     ----------
@@ -151,21 +134,18 @@ def circular_two_sided(
         pixels per degree [vertical, horizontal]
     shape : Sequence[Number, Number], Number, or None (default)
         shape [height, width] of image, in pixels
-    frequency : Number, or None (default)
-        spatial frequency of circular grating, in cycles per degree
-    n_rings : int, or None (default)
-        number of rings
-    ring_width : Number, or None (default)
-        width of a single ring, in degrees
-    phase_shift : float
-        phase shift of grating in degrees
-    intensity_target : float (optional)
-        intensity value of target ring(s), by default 0.5
-    intensity_rings : Sequence[Number, ...]
-        intensity value for each ring, from inside to out, by default [1,0]
-        If fewer intensities are passed than number of radii, cycles through intensities
+    radii : Sequence[Number] or None (default)
+        radii of each ring, in degrees visual angle
+    rotation : float, optional
+        rotation (in degrees), counterclockwise, by default 0.0 (horizonal)
+    intensity_rings : Sequence[float, float]
+        intensities of rings, by default (1.0, 0.0)
     intensity_background : float (optional)
         intensity value of background, by default 0.5
+    intensity_target : float, or Sequence[float, ...], optional
+        intensity value for each target, by default 0.5.
+        Can specify as many intensities as number of target_indices;
+        If fewer intensities are passed than target_indices, cycles through intensities
     origin : "corner", "mean" or "center"
         if "corner": set origin to upper left corner
         if "mean": set origin to hypothetical image center (default)
@@ -175,43 +155,29 @@ def circular_two_sided(
     -------
     dict[str, Any]
         dict with the stimulus (key: "img"),
-        mask with integer index for each target (key: "target_mask"),
+        mask with integer index for each frame (key: "target_mask"),
         and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Bindman, D., & Chubb, C. (2004).
-        Brightness assimilation in bullseye displays.
-        Vision Research, 44, 309-319.
-        https://doi.org/10.1016/S0042-6989(03)00430-9
-    Hong, S. W., and Shevell, S. K. (2004).
-        Brightness contrast and assimilation from patterned inducing backgrounds.
-        Vision Research, 44, 35-43.
-        https://doi.org/10.1016/j.visres.2003.07.010
-    Howe, P. D. L. (2005).
-        White's effect:
-        removing the junctions but preserving the strength of the illusion.
-        Perception, 34, 557-564.
-        https://doi.org/10.1068/p5414
     """
-
-    stim = rings.circular_two_sided(
+    stim = rings.circular_generalized(
+        radii=radii,
         visual_size=visual_size,
         ppd=ppd,
         shape=shape,
-        frequency=frequency,
-        n_rings=n_rings,
-        ring_width=ring_width,
-        phase_shift=phase_shift,
         intensity_rings=intensity_rings,
         intensity_background=intensity_background,
+        target_indices=1,
         intensity_target=intensity_target,
-        target_indices=0,
         origin=origin,
     )
     return stim
 
 
+circular_two_sided = make_two_sided(
+    circular,
+    two_sided_params=("intensity_rings", "intensity_target"),
+)
+
+
 def rectangular(
     visual_size=None,
     ppd=None,
@@ -221,7 +187,7 @@ def rectangular(
     frame_width=None,
     rotation=0.0,
     phase_shift=0,
-    intensity_frames=(1.0, 0.0),
+    intensity_frames=(0.0, 1.0),
     intensity_background=0.5,
     intensity_target=0.5,
     origin="mean",
@@ -289,7 +255,7 @@ def rectangular(
         rotation=rotation,
         phase_shift=phase_shift,
         intensity_frames=intensity_frames,
-        target_indices=0,
+        target_indices=1,
         intensity_target=intensity_target,
         origin=origin,
         clip=clip,
@@ -304,7 +270,7 @@ def rectangular_generalized(
     shape=None,
     radii=None,
     rotation=0.0,
-    intensity_frames=(1.0, 0.0),
+    intensity_frames=(0.0, 1.0),
     intensity_background=0.5,
     intensity_target=0.5,
     origin="mean",
@@ -360,90 +326,10 @@ def rectangular_generalized(
     return stim
 
 
-def rectangular_two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    frequency=None,
-    n_frames=None,
-    frame_width=None,
-    rotation=0.0,
-    phase_shift=0,
-    intensity_target=0.5,
-    intensity_frames=(1.0, 0.0),
-    intensity_background=0.5,
-    origin="mean",
-):
-    """Two-sided square "bullseye", i.e., set of rings with target in center
-
-    Essentially frames(target_indices=1)
-
-    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
-    frequency : Number, or None (default)
-        spatial frequency of grating, in cycles per degree visual angle
-    n_frames : int, or None (default)
-        number of frames in the grating
-    frame_width : Number, or None (default)
-        width of a single frame, in degrees visual angle
-    rotation : float, optional
-        rotation (in degrees), counterclockwise, by default 0.0 (horizonal)
-    phase_shift : float
-        phase shift of grating in degrees
-    period : "full", "half", "ignore" (default)
-        whether to ensure the grating only has "full" periods,
-        half "periods", or no guarantees ("ignore")
-    intensity_target : float, or Sequence[float, ...], optional
-        intensity value for each target, by default 0.5.
-        Can specify as many intensities as number of target_indices;
-        If fewer intensities are passed than target_indices, cycles through intensities
-    intensity_frames : Sequence[float, float]
-        min and max intensity of square-wave, by default (0.0, 1.0)
-    intensity_background : float (optional)
-        intensity value of background, by default 0.5
-    origin : "corner", "mean" or "center"
-        if "corner": set origin to upper left corner
-        if "mean": set origin to hypothetical image center (default)
-        if "center": set origin to real center (closest existing value to mean)
-
-
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for each target (key: "target_mask"),
-        and additional keys containing stimulus parameters
-
-    References
-    ----------
-        Bindman, D., & Chubb, C. (2004).
-        Brightness assimilation in bullseye displays.
-        Vision Research, 44, 309-319.
-        https://doi.org/10.1016/S0042-6989(03)00430-9
-    """
-
-    stim = rings.rectangular_two_sided(
-        visual_size=visual_size,
-        ppd=ppd,
-        shape=shape,
-        frequency=frequency,
-        n_frames=n_frames,
-        frame_width=frame_width,
-        rotation=rotation,
-        phase_shift=phase_shift,
-        intensity_frames=intensity_frames,
-        intensity_background=intensity_background,
-        intensity_target=intensity_target,
-        target_indices=0,
-        origin=origin,
-    )
-    return stim
+rectangular_two_sided = make_two_sided(
+    rectangular,
+    two_sided_params=("intensity_frames", "intensity_target"),
+)
 
 
 def overview(**kwargs):
@@ -463,9 +349,9 @@ def overview(**kwargs):
     # fmt: off
     stimuli = {
         "bullseyes_circular": circular(**default_params, frequency=1.0, clip=True),
-        "bullseyes_circular_2sided": circular_two_sided(**default_params, frequency=1.0),
+        "bullseyes_circular_2sided": circular_two_sided(**default_params, frequency=1.0, intensity_rings=((0.0, 1.0), (1.0, 0.0))),
         "bullseyes_rectangular": rectangular(**default_params, frequency=1.0, clip=True),
-        "bullseyes_rectangular_2sided": rectangular_two_sided(**default_params, frequency=1.0),
+        "bullseyes_rectangular_2sided": rectangular_two_sided(**default_params, frequency=1.0, intensity_frames=((0.0, 1.0), (1.0, 0.0))),
     }
     # fmt: on
 
diff --git a/stimupy/stimuli/checkerboards.py b/stimupy/stimuli/checkerboards.py
index 9a966ba4..bce13ec4 100644
--- a/stimupy/stimuli/checkerboards.py
+++ b/stimupy/stimuli/checkerboards.py
@@ -72,7 +72,7 @@ def extend_target_idx(target_index, offsets=[(-1, 0), (0, -1), (0, 0), (0, 1), (
     return extended_idc
 
 
-def add_targets(checkerboard_stim, target_indices, extend_targets=False, intensity_target=0.5):
+def add_targets(checkerboard_stim, target_indices=(), extend_targets=False, intensity_target=0.5):
     """Add targets to a checkerboard stimulus
 
     Parameters
@@ -94,7 +94,10 @@ def add_targets(checkerboard_stim, target_indices, extend_targets=False, intensi
         mask with integer index for each target (key: "target_mask"),
         and additional keys containing stimulus parameters
     """
+    if target_indices is None:
+        target_indices = ()
     mask = np.zeros(checkerboard_stim["shape"])
+
     for i, target in enumerate(target_indices):
         if extend_targets:
             target_idc = extend_target_idx(target)
@@ -116,7 +119,7 @@ def checkerboard(
     frequency=None,
     board_shape=None,
     check_visual_size=None,
-    target_indices=None,
+    target_indices=(),
     extend_targets=False,
     period="ignore",
     rotation=0.0,
@@ -298,15 +301,13 @@ def checkerboard(
     }
 
     # Add targets
-    if target_indices is not None:
-        stim = add_targets(
-            stim,
-            target_indices=target_indices,
-            extend_targets=extend_targets,
-            intensity_target=intensity_target,
-        )
-    else:
-        stim["target_mask"] = np.zeros(stim["shape"])
+    stim = add_targets(
+        stim,
+        target_indices=target_indices,
+        extend_targets=extend_targets,
+        intensity_target=intensity_target,
+    )
+
     return stim
 
 
diff --git a/stimupy/stimuli/cubes.py b/stimupy/stimuli/cubes.py
index eeca6e26..ce9a5699 100644
--- a/stimupy/stimuli/cubes.py
+++ b/stimupy/stimuli/cubes.py
@@ -13,7 +13,7 @@ def varying_cells(
     cell_lengths=None,
     cell_thickness=None,
     cell_spacing=None,
-    target_indices=None,
+    target_indices=(),
     intensity_background=0.0,
     intensity_cells=1.0,
     intensity_target=0.5,
@@ -154,7 +154,7 @@ def cube(
     ppd=None,
     shape=None,
     n_cells=None,
-    target_indices=None,
+    target_indices=(),
     cell_thickness=None,
     cell_spacing=None,
     intensity_background=0.0,
@@ -318,8 +318,8 @@ def overview(**kwargs):
 
     # fmt: off
     stimuli = {
-        "cubes_regular": cube(**default_params, visual_size=10, n_cells=5, target_indices=(1,2)),
-        "cubes_variable": varying_cells(**default_params, cell_lengths=(2,4,2), target_indices=1),
+        "cubes_regular": cube(**default_params, visual_size=10, n_cells=5, target_indices=(1, 2)),
+        "cubes_variable": varying_cells(**default_params, cell_lengths=(2, 4, 2), target_indices=1),
     }
     # fmt: on
 
diff --git a/stimupy/stimuli/delboeufs.py b/stimupy/stimuli/delboeufs.py
index 079f4b44..a9b9638b 100644
--- a/stimupy/stimuli/delboeufs.py
+++ b/stimupy/stimuli/delboeufs.py
@@ -1,6 +1,6 @@
 from stimupy.components import lines
 from stimupy.components.shapes import disc
-from stimupy.utils import resolution, stack_dicts
+from stimupy.utils import make_two_sided
 
 __all__ = ["delboeuf", "two_sided"]
 
@@ -84,101 +84,23 @@ def delboeuf(
     inner["line_mask"] = outer["line_mask"]
     inner["outer_radius"] = outer_radius
     inner["target_radius"] = target_radius
-    inner["target_mask"] = inner["ring_mask"]
+    inner["target_mask"] = inner.pop("ring_mask")
     inner["outer_line_width"] = outer["line_width"]
     inner["intensity_outer_line"] = intensity_outer_line
     inner["intensity_target"] = intensity_target
-    del inner["ring_mask"]
     return inner
 
 
-def two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    outer_radii=None,
-    outer_line_width=0,
-    target_radius=None,
-    intensity_outer_line=0.0,
-    intensity_target=0.0,
-    intensity_background=1.0,
-):
-    """Two-sided Delboeuf's (1865) stimulus
-
-    Parameters
-    ----------
-    visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, 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 grating, in pixels
-    outer_radii : Sequence[Number, Number], or None (default)
-        radius of outer circle
-    outer_line_width : Number
-        line width of outer circle in degrees visual angle
-        if 0 (default), set line width to 1 px
-    target_radius : Number or None (default)
-        radius of target circle
-    intensity_outer_line : Number
-        intensity value of outer circle line (default: 0)
-    intensity_target : Number
-        intensity value of target (default: 0)
-    intensity_background : Number
-        intensity value of background (default: 1)
-
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for the target (key: "target_mask"),
-        and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Delboeuf, F. J. (1865).
-        Note sur certaines illusions d'optique:
-        Essai d'une théorie psychophysique de la maniere
-        dont l'oeil apprécie les distances et les angles.
-        Bulletins de l'Académie Royale des Sciences, Lettres et
-        Beaux-arts de Belgique, 19, 195-216.
-    """
-    if outer_radii is None:
-        raise ValueError("two_sided() missing argument 'outer_radii' which is not 'None'")
-    if target_radius is None:
-        raise ValueError("delboeuf() missing argument 'target_radius' which is not 'None'")
-
-    # Resolve resolution
-    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-
-    stim1 = delboeuf(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        outer_radius=outer_radii[0],
-        outer_line_width=outer_line_width,
-        target_radius=target_radius,
-        intensity_outer_line=intensity_outer_line,
-        intensity_target=intensity_target,
-        intensity_background=intensity_background,
-    )
-
-    stim2 = delboeuf(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        outer_radius=outer_radii[1],
-        outer_line_width=outer_line_width,
-        target_radius=target_radius,
-        intensity_outer_line=intensity_outer_line,
-        intensity_target=intensity_target,
-        intensity_background=intensity_background,
-    )
-
-    stim = stack_dicts(stim1, stim2)
-    stim["shape"] = shape
-    stim["visual_size"] = visual_size
-    del stim["outer_radius"]
-    stim["outer_radii"] = outer_radii
-    return stim
+two_sided = make_two_sided(
+    delboeuf,
+    two_sided_params=(
+        "outer_radius",
+        "outer_line_width",
+        "target_radius",
+        "intensity_target",
+        "intensity_outer_line",
+    ),
+)
 
 
 def overview(**kwargs):
@@ -198,7 +120,7 @@ def overview(**kwargs):
     # fmt: off
     stimuli = {
         "delboeuf": delboeuf(**default_params, target_radius=1, outer_radius=4),
-        "delboeuf_2sided": two_sided(**default_params, target_radius=1, outer_radii=(2, 1.1)),
+        "delboeuf_2sided": two_sided(**default_params, target_radius=1, outer_radius=(2, 1.1)),
     }
     # fmt: on
 
diff --git a/stimupy/stimuli/gabors.py b/stimupy/stimuli/gabors.py
index f39af859..bf1bb81b 100644
--- a/stimupy/stimuli/gabors.py
+++ b/stimupy/stimuli/gabors.py
@@ -14,7 +14,7 @@ def gabor(
     period="ignore",
     rotation=0.0,
     phase_shift=0,
-    intensity_bars=(0.0, 1.0),
+    intensities=(0.0, 1.0),
     origin="center",
     round_phase_width=False,
     sigma=None,
@@ -45,7 +45,7 @@ def gabor(
         rotation (in degrees), counterclockwise, by default 0.0 (horizonal)
     phase_shift : float
         phase shift of grating in degrees
-    intensity_bars : Sequence[float, ...]
+    intensities : Sequence[float, ...]
         maximal intensity value for each bar, by default (0.0, 1.0).
     origin : "corner", "mean" or "center"
         if "corner": set origin to upper left corner
@@ -72,7 +72,7 @@ def gabor(
         period=period,
         rotation=rotation,
         phase_shift=phase_shift,
-        intensities=intensity_bars,
+        intensities=intensities,
         origin=origin,
         distance_metric="oblique",
         round_phase_width=round_phase_width,
@@ -84,7 +84,7 @@ def gabor(
         sigma=sigma,
         origin=origin,
     )
-    mean_int = (intensity_bars[0] + intensity_bars[1]) / 2
+    mean_int = (intensities[0] + intensities[1]) / 2
     stim["img"] = (stim["img"] - mean_int) * gaussian_window["img"] + mean_int
     del stim["intensities"]
 
@@ -92,7 +92,7 @@ def gabor(
         **stim,
         "sigma": sigma,
         "gaussian_mask": gaussian_window["gaussian_mask"],
-        "intensity_bars": intensity_bars,
+        "intensities": intensities,
     }
 
 
diff --git a/stimupy/stimuli/gratings.py b/stimupy/stimuli/gratings.py
index efacd7ca..ec8fe3a9 100644
--- a/stimupy/stimuli/gratings.py
+++ b/stimupy/stimuli/gratings.py
@@ -32,7 +32,7 @@ def on_uniform(
     period="ignore",
     rotation=0.0,
     phase_shift=0,
-    intensity_bars=(1.0, 0.0),
+    intensity_bars=(0.0, 1.0),
     target_indices=(),
     intensity_target=0.5,
     intensity_background=0.5,
@@ -259,7 +259,7 @@ def phase_shifted(
     period="ignore",
     distance_metric="horizontal",
     phase_shift=0,
-    intensity_bars=(1.0, 0.0),
+    intensity_bars=(0.0, 1.0),
     target_size=None,
     target_phase_shift=0,
     intensity_target=0.5,
@@ -299,7 +299,7 @@ def phase_shifted(
         rotation=rotation,
         phase_shift=0,
         period=period,
-        intensity_bars=(0, intensity_target),
+        intensity_bars=(intensity_target, 0),
         origin=origin,
         round_phase_width=True,
     )
@@ -473,7 +473,7 @@ def grating_induction_blur(
     period="ignore",
     rotation=0.0,
     phase_shift=0,
-    intensity_bars=(1.0, 0.0),
+    intensity_bars=(0.0, 1.0),
     target_width=None,
     sigma=None,
     intensity_target=0.5,
diff --git a/stimupy/stimuli/mondrians.py b/stimupy/stimuli/mondrians.py
index 7d6ac9da..a46792d7 100644
--- a/stimupy/stimuli/mondrians.py
+++ b/stimupy/stimuli/mondrians.py
@@ -150,7 +150,7 @@ def corrugated_mondrian(
     ncols=None,
     depths=0,
     intensities=(0, 1),
-    target_indices=None,
+    target_indices=(),
     intensity_background=0.5,
     intensity_target=None,
 ):
@@ -208,6 +208,9 @@ def corrugated_mondrian(
     if isinstance(depths, (float, int)):
         depths = (depths,) * nrows
 
+    if target_indices is None:
+        target_indices = ()
+
     if len(depths) != nrows:
         raise ValueError("Unclear number of rows. Check nrows, intensities and depths.")
 
@@ -259,7 +262,7 @@ def corrugated_mondrian(
             poses.append(mpos)
             intenses.append(mint)
 
-            if (target_indices is not None) and (r, c) in target_indices:
+            if (r, c) in target_indices:
                 tlist.append(target_counter)
 
             xst += mwidth_px
diff --git a/stimupy/stimuli/mueller_lyers.py b/stimupy/stimuli/mueller_lyers.py
index d5389d64..b44cb2e0 100644
--- a/stimupy/stimuli/mueller_lyers.py
+++ b/stimupy/stimuli/mueller_lyers.py
@@ -3,7 +3,7 @@
 import numpy as np
 
 from stimupy.components import lines
-from stimupy.utils import resolution, stack_dicts
+from stimupy.utils import make_two_sided, resolution
 
 __all__ = [
     "mueller_lyer",
@@ -35,18 +35,18 @@ def mueller_lyer(
         shape [height, width] of grating, in pixels
     outer_lines_length : Number
         length of outer lines in degrees visual angle
-    outer_lines_angle : Number
-        angle of outer lines in degrees. Must be between -180 and 180 degrees.
+    outer_lines_angle : Number (optional)
+        angle of outer lines in degrees, by default 45. Must be between -180 and 180 degrees.
     target_length : Number
         length of target line in degrees visual angle
-    line_width :
+    line_width : Number (optional)
         line width in degrees visual angle; if 0 (default), line width is 1 px
-    intensity_outer_lines : Number
-        intensity value of outer lines
-    intensity_target : Number
-        intensity value of target line
-    intensity_background : Number
-        intensity value of background
+    intensity_outer_lines : Number (optional)
+        intensity value of outer lines, by default 0.01
+    intensity_target : Number (optional)
+        intensity value of target line, by default 0.5
+    intensity_background : Number (optional)
+        intensity value of background, by default 0.0
 
     Returns
     -------
@@ -174,88 +174,17 @@ def mueller_lyer(
     return target_line
 
 
-def two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    outer_lines_length=None,
-    outer_lines_angle=45,
-    target_length=None,
-    line_width=0,
-    intensity_outer_lines=1.0,
-    intensity_target=0.5,
-    intensity_background=0.0,
-):
-    """Two-sided Mueller-Lyer's (1896) illusion
-
-    Parameters
-    ----------
-    visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, 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 grating, in pixels
-    outer_lines_length : Number
-        length of outer lines in degrees visual angle
-    outer_lines_angle : Number
-        angle of outer lines in degrees. Must be between -180 and 180 degrees.
-    target_length : Number
-        length of target line in degrees visual angle
-    line_width :
-        line width in degrees visual angle; if 0 (default), line width is 1 px
-    intensity_outer_lines : Number
-        intensity value of outer lines
-    intensity_target : Number
-        intensity value of target line
-    intensity_background : Number
-        intensity value of background
-
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for each target (key: "target_mask"),
-        and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Mueller-Lyer, F. (1896).
-        Zur Lehre von den optischen Taeuschungen.
-        Ueber Kontrast und Konfluxion.
-        Zeitschrift fuer Psychologie und Physiologie der Sinnesorgane, IX, 1-16.
-    """
-    # Resolve resolution
-    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-
-    stim1 = mueller_lyer(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        outer_lines_length=outer_lines_length,
-        outer_lines_angle=outer_lines_angle,
-        target_length=target_length,
-        line_width=line_width,
-        intensity_outer_lines=intensity_outer_lines,
-        intensity_target=intensity_target,
-        intensity_background=intensity_background,
-    )
-
-    stim2 = mueller_lyer(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        outer_lines_length=outer_lines_length,
-        outer_lines_angle=outer_lines_angle + 90,
-        target_length=target_length,
-        line_width=line_width,
-        intensity_outer_lines=intensity_outer_lines,
-        intensity_target=intensity_target,
-        intensity_background=intensity_background,
-    )
-
-    stim = stack_dicts(stim1, stim2)
-    stim["shape"] = shape
-    stim["visual_size"] = visual_size
-    return stim
+two_sided = make_two_sided(
+    mueller_lyer,
+    two_sided_params=(
+        "outer_lines_length",
+        "outer_lines_angle",
+        "target_length",
+        "line_width",
+        "intensity_outer_lines",
+        "intensity_target",
+    ),
+)
 
 
 def overview(**kwargs):
@@ -274,15 +203,14 @@ def overview(**kwargs):
 
     stim_params = {
         "outer_lines_length": 1,
-        "outer_lines_angle": 45,
         "target_length": 3,
         "line_width": 0.1,
     }
 
     # fmt: off
     stimuli = {
-        "mueller_lyer": mueller_lyer(**default_params, **stim_params),
-        "mueller_lyer_2sided": two_sided(**default_params, **stim_params),
+        "mueller_lyer": mueller_lyer(**default_params, **stim_params, outer_lines_angle=45),
+        "mueller_lyer_2sided": two_sided(**default_params, **stim_params, outer_lines_angle=(30, 135)),
     }
     # fmt: on
 
diff --git a/stimupy/stimuli/pinwheels.py b/stimupy/stimuli/pinwheels.py
index e4958fc4..52c5411f 100644
--- a/stimupy/stimuli/pinwheels.py
+++ b/stimupy/stimuli/pinwheels.py
@@ -2,9 +2,10 @@
 
 import numpy as np
 
+from stimupy.components import combine_masks, draw_regions
 from stimupy.components.shapes import circle
 from stimupy.components.shapes import ring as ring_shape
-from stimupy.stimuli import waves
+from stimupy.stimuli import mask_targets, waves
 
 __all__ = [
     "pinwheel",
@@ -19,10 +20,10 @@ def pinwheel(
     n_segments=None,
     segment_width=None,
     rotation=0.0,
-    target_indices=None,
+    target_indices=(),
     target_width=None,
     target_center=None,
-    intensity_segments=(1.0, 0.0),
+    intensity_segments=(0.0, 1.0),
     intensity_background=0.5,
     intensity_target=0.5,
     origin="mean",
@@ -87,6 +88,7 @@ def pinwheel(
         https://doi.org/10.1016/j.visres.2007.02.017
     """
 
+    # Draw angular grating
     stim = waves.square_angular(
         visual_size=visual_size,
         ppd=ppd,
@@ -99,76 +101,79 @@ def pinwheel(
         origin=origin,
         intensity_segments=intensity_segments,
     )
-    radius = min(stim["visual_size"]) / 2
 
-    stim["target_indices"] = target_indices
-    stim["target_center"] = target_center
-    stim["target_width"] = target_width
-    stim["intensity_target"] = intensity_target
-
-    circle_mask = circle(
+    # Mask to circular aperture
+    radius = min(stim["visual_size"]) / 2
+    stim["circle_mask"] = circle(
         visual_size=visual_size,
         ppd=ppd,
         shape=shape,
         radius=radius,
         origin=origin,
     )["circle_mask"]
+    stim["segment_mask"] = np.where(stim["circle_mask"], stim["segment_mask"], 0)
+    stim["img"] = np.where(stim["circle_mask"], stim["img"], intensity_background)
+    stim["intensity_background"] = intensity_background
 
-    stim["img"] = np.where(circle_mask, stim["img"], intensity_background)
+    # Target segment mask
+    if isinstance(target_indices, (int, float)):
+        target_indices = (target_indices,)
+    target_segment_mask = mask_targets(
+        element_mask=stim["grating_mask"], target_indices=target_indices
+    )
+    stim["target_indices"] = target_indices
 
+    # Mask ring regions
     if target_center is None:
         target_center = radius / 2
-
-    # Place target(s)
-    if isinstance(target_indices, (int)):
-        target_indices = [
-            target_indices,
-        ]
     if isinstance(target_center, (int, float)):
-        target_center = [
-            target_center,
-        ]
+        target_center = (target_center,)
+    stim["target_center"] = target_center
+    target_center = tuple(itertools.islice(itertools.cycle(target_center), len(target_indices)))
+
+    if target_width is None:
+        raise ValueError("pinwheel() missing argument 'target_width' which is not 'None'")
     if isinstance(target_width, (int, float)):
-        target_width = [
-            target_width,
-        ]
-    if isinstance(intensity_target, (int, float)):
-        intensity_target = [
-            intensity_target,
-        ]
-
-    # Initiate target mask
-    target_mask = np.zeros_like(stim["grating_mask"])
-
-    if target_indices is not None:
-        if target_width is None:
-            raise ValueError("pinwheel() missing argument 'target_width' which is not 'None'")
-
-        target_center = itertools.cycle(target_center)
-        target_width = itertools.cycle(target_width)
-        intensity_target = itertools.cycle(intensity_target)
-
-        for target_idx, (segment_idx, center, width, intensity) in enumerate(
-            zip(target_indices, target_center, target_width, intensity_target)
-        ):
-            # Draw ring
-            inner_radius = center - (width / 2)
-            outer_radius = center + (width / 2)
-            if inner_radius < 0 or outer_radius > np.min(visual_size) / 2:
-                raise ValueError("target does not fully fit into pinwheel")
-            ring_stim = ring_shape(
-                radii=[inner_radius, outer_radius],
-                intensity_ring=intensity,
-                visual_size=stim["visual_size"],
-                ppd=stim["ppd"],
-                shape=stim["shape"],
-            )
-            condition1 = stim["grating_mask"] == segment_idx
-            condition2 = ring_stim["ring_mask"] == 1
-            target_mask = np.where(condition1 & condition2, target_idx + 1, target_mask)
-            stim["img"] = np.where(target_mask == (target_idx + 1), intensity, stim["img"])
-    stim["target_mask"] = target_mask
-    stim["intensity_background"] = intensity_background
+        target_width = (target_width,)
+    stim["target_width"] = target_width
+    target_width = tuple(itertools.islice(itertools.cycle(target_width), len(target_indices)))
+
+    target_ring_masks = []
+    for target_idx, (center, width) in enumerate(zip(target_center, target_width)):
+        # Draw ring
+        inner_radius = center - (width / 2)
+        outer_radius = center + (width / 2)
+        if inner_radius < 0 or outer_radius > np.min(visual_size) / 2:
+            raise ValueError("target does not fully fit into pinwheel")
+        ring = ring_shape(
+            radii=[inner_radius, outer_radius],
+            intensity_ring=target_idx,
+            visual_size=stim["visual_size"],
+            ppd=stim["ppd"],
+            shape=stim["shape"],
+        )
+        target_ring_masks.append(ring["ring_mask"])
+
+    # Combine segment & ring masks
+    target_masks = []
+    for target_idx, ring_mask in enumerate(target_ring_masks):
+        # Find where ring intesects with target segment
+        target_mask = (target_segment_mask == target_idx + 1) & ring_mask
+        target_masks.append(target_mask)
+
+    # Combine target masks
+    if len(target_masks) > 0:
+        target_mask = combine_masks(*target_masks)
+    else:
+        target_mask = np.zeros_like(stim["img"])
+    stim["target_mask"] = target_mask.astype(int)
+
+    # Draw target(s)
+    stim["img"] = np.where(
+        target_mask, draw_regions(mask=target_mask, intensities=intensity_target), stim["img"]
+    )
+    stim["intensity_target"] = intensity_target
+
     return stim
 
 
diff --git a/stimupy/stimuli/rings.py b/stimupy/stimuli/rings.py
index bb16ffbc..d192a91e 100644
--- a/stimupy/stimuli/rings.py
+++ b/stimupy/stimuli/rings.py
@@ -1,14 +1,13 @@
-import itertools
-
-import numpy as np
-
 from stimupy.components.frames import frames
+from stimupy.components.radials import rings
+from stimupy.stimuli import place_targets
 from stimupy.stimuli.waves import square_radial as circular
 from stimupy.stimuli.waves import square_rectilinear as rectangular
-from stimupy.utils import resolution, stack_dicts
+from stimupy.utils import make_two_sided
 
 __all__ = [
     "circular",
+    "circular_generalized",
     "circular_two_sided",
     "rectangular",
     "rectangular_generalized",
@@ -16,34 +15,18 @@
 ]
 
 
-def circular_two_sided(
+def circular_generalized(
     visual_size=None,
     ppd=None,
     shape=None,
-    frequency=None,
-    n_rings=None,
-    ring_width=None,
-    phase_shift=0,
-    target_indices=None,
-    intensity_target=0.5,
-    intensity_rings=(1.0, 0.0),
+    radii=None,
+    intensity_rings=(0.0, 1.0),
     intensity_background=0.5,
+    target_indices=(),
+    intensity_target=0.5,
     origin="mean",
 ):
-    """Two-sided circular grating, with one or more target rings
-
-    Specification of the number of rings, and their width can be done in two ways:
-    a ring_width (in degrees) and n_rings, and/or by specifying the spatial frequency
-    of a circular grating (in cycles per degree)
-
-    The total shape (in pixels) and visual size (in degrees) has to match the
-    specification of the rings and their widths.
-    Thus, not all 6 parameters (visual_size, ppd, shape, frequency, ring_width, n_rings)
-    have to be specified, as long as both the resolution, and the distribution of rings,
-    can be resolved.
-
-    Note: all rings in a grating have the same width -- if more control is required
-    see disc_and_rings
+    """Sequential set of circular rings with specified radii and targets
 
     Parameters
     ----------
@@ -53,23 +36,18 @@ def circular_two_sided(
         pixels per degree [vertical, horizontal]
     shape : Sequence[Number, Number], Number, or None (default)
         shape [height, width] of image, in pixels
-    frequency : Number, or None (default)
-        spatial frequency of circular grating, in cycles per degree
-    n_rings : int, or None (default)
-        number of rings
-    ring_width : Number, or None (default)
-        width of a single ring, in degrees
-    phase_shift : float
-        phase shift of grating in degrees
-    target_indices : int or Sequence[int, ] (optional)
-        indices of target discs. If not specified, use middle ring (round down)
-    intensity_target : float (optional)
-        intensity value of target ring(s), by default 0.5
-    intensity_rings : Sequence[Number, ...]
-        intensity value for each ring, from inside to out, by default [1,0]
-        If fewer intensities are passed than number of radii, cycles through intensities
+    radii : Sequence[Number] or None (default)
+        radii of each ring, in degrees visual angle
+    intensity_rings : Sequence[float, float]
+        intensities of rings, by default (1.0, 0.0)
     intensity_background : float (optional)
         intensity value of background, by default 0.5
+    target_indices : int, or Sequence[int, ...]
+        indices rings where targets will be placed
+    intensity_target : float, or Sequence[float, ...], optional
+        intensity value for each target, by default 0.5.
+        Can specify as many intensities as number of target_indices;
+        If fewer intensities are passed than target_indices, cycles through intensities
     origin : "corner", "mean" or "center"
         if "corner": set origin to upper left corner
         if "mean": set origin to hypothetical image center (default)
@@ -79,68 +57,43 @@ def circular_two_sided(
     -------
     dict[str, Any]
         dict with the stimulus (key: "img"),
-        mask with integer index for each target (key: "target_mask"),
+        mask with integer index for each ring (key: "target_mask"),
         and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Hong, S. W., and Shevell, S. K. (2004).
-        Brightness contrast and assimilation from patterned inducing backgrounds.
-        Vision Research, 44, 35-43.
-        https://doi.org/10.1016/j.visres.2003.07.010
-    Howe, P. D. L. (2005).
-        White's effect:
-        removing the junctions but preserving the strength of the illusion.
-        Perception, 34, 557-564.
-        https://doi.org/10.1068/p5414
     """
 
-    # Resolve resolution
-    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-
-    stim1 = circular(
-        visual_size=(visual_size[0], visual_size[1] / 2),
+    # Rings component
+    stim = rings(
+        radii=radii,
+        visual_size=visual_size,
         ppd=ppd,
-        frequency=frequency,
-        n_rings=n_rings,
-        ring_width=ring_width,
-        phase_shift=phase_shift,
-        target_indices=target_indices,
-        intensity_target=intensity_target,
+        shape=shape,
         intensity_rings=intensity_rings,
         intensity_background=intensity_background,
         origin=origin,
-        clip=True,
     )
 
-    stim2 = circular(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        frequency=frequency,
-        n_rings=n_rings,
-        ring_width=ring_width,
-        phase_shift=phase_shift,
+    # Place target(s)
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="ring_mask",
         target_indices=target_indices,
         intensity_target=intensity_target,
-        intensity_rings=intensity_rings[::-1],
-        intensity_background=intensity_background,
-        origin=origin,
-        clip=True,
     )
 
-    stim = stack_dicts(stim1, stim2)
-    stim["shape"] = shape
-    stim["visual_size"] = visual_size
-    stim["intensity_background"] = intensity_background
     return stim
 
 
+circular_two_sided = make_two_sided(
+    circular, two_sided_params=("intensity_rings", "intensity_target", "target_indices")
+)
+
+
 def rectangular_generalized(
     visual_size=None,
     ppd=None,
     shape=None,
     radii=None,
-    intensity_frames=(1.0, 0.0),
+    intensity_frames=(0.0, 1.0),
     intensity_background=0.5,
     target_indices=(),
     intensity_target=0.5,
@@ -198,139 +151,22 @@ def rectangular_generalized(
     )
     stim["intensity_target"] = intensity_target
 
-    # Resolve target parameters
-    if isinstance(target_indices, (int)):
-        target_indices = [
-            target_indices,
-        ]
-    if isinstance(intensity_target, (int, float)):
-        intensity_target = [
-            intensity_target,
-        ]
-    intensity_target = itertools.cycle(intensity_target)
-
     # Place target(s)
-    targets_mask = np.zeros_like(stim["frame_mask"])
-    for target_idx, (bar_idx, intensity) in enumerate(zip(target_indices, intensity_target)):
-        targets_mask = np.where(stim["frame_mask"] == bar_idx, target_idx + 1, targets_mask)
-        stim["img"] = np.where(targets_mask == target_idx + 1, intensity, stim["img"])
-        if bar_idx > stim["frame_mask"].max():
-            raise ValueError("target_idx is outside stimulus")
-
-    # Update and return stimulus
-    stim["target_mask"] = targets_mask
-    stim["target_indices"] = target_indices
-    return stim
-
-
-def rectangular_two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    frequency=None,
-    n_frames=None,
-    frame_width=None,
-    phase_shift=0,
-    intensity_frames=(1.0, 0.0),
-    intensity_background=0.5,
-    target_indices=(),
-    intensity_target=0.5,
-    origin="mean",
-    rotation=0.0,
-):
-    """Draw set of square frames, with some frame(s) as target
-
-    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
-    frequency : Number, or None (default)
-        spatial frequency of grating, in cycles per degree visual angle
-    n_frames : int, or None (default)
-        number of frames in the grating
-    frame_width : Number, or None (default)
-        width of a single frame, in degrees visual angle
-    phase_shift : float
-        phase shift of grating in degrees
-    intensity_frames : Sequence[float, float]
-        min and max intensity of square-wave, by default (0.0, 1.0)
-    intensity_background : float (optional)
-        intensity value of background, by default 0.5
-    target_indices : int, or Sequence[int, ...]
-        indices frames where targets will be placed
-    intensity_target : float, or Sequence[float, ...], optional
-        intensity value for each target, by default 0.5.
-        Can specify as many intensities as number of target_indices;
-        If fewer intensities are passed than target_indices, cycles through intensities
-    origin : "corner", "mean" or "center"
-        if "corner": set origin to upper left corner
-        if "mean": set origin to hypothetical image center (default)
-        if "center": set origin to real center (closest existing value to mean)
-    rotation : float, optional
-        rotation (in degrees), counterclockwise, by default 0.0 (horizonal)
-
-
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for each target (key: "target_mask"),
-        and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Domijan, D. (2015).
-        A neurocomputational account
-        of the role of contour facilitation in brightness perception.
-        Frontiers in Human Neuroscience, 9, 93.
-        https://doi.org/10.3389/fnhum.2015.00093
-    """
-
-    # Resolve resolution
-    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-
-    stim1 = rectangular(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        frequency=frequency,
-        n_frames=n_frames,
-        frame_width=frame_width,
-        phase_shift=phase_shift,
-        intensity_target=intensity_target,
-        intensity_frames=intensity_frames,
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="frame_mask",
         target_indices=target_indices,
-        origin=origin,
-        clip=True,
-        intensity_background=intensity_background,
-        rotation=rotation,
-    )
-
-    stim2 = rectangular(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        frequency=frequency,
-        n_frames=n_frames,
-        frame_width=frame_width,
-        phase_shift=phase_shift,
         intensity_target=intensity_target,
-        intensity_frames=intensity_frames[::-1],
-        target_indices=target_indices,
-        origin=origin,
-        clip=True,
-        intensity_background=intensity_background,
-        rotation=rotation,
     )
 
-    stim = stack_dicts(stim1, stim2)
-    stim["shape"] = shape
-    stim["visual_size"] = visual_size
     return stim
 
 
+rectangular_two_sided = make_two_sided(
+    rectangular, two_sided_params=("intensity_frames", "intensity_target", "target_indices")
+)
+
+
 def overview(**kwargs):
     """Generate example stimuli from this module
 
@@ -350,12 +186,12 @@ def overview(**kwargs):
         "rings_circular": circular(**default_params, frequency=1.0),
         "rings_circular_with_targets": circular(**default_params, frequency=1.0, target_indices=3),
         "rings_circular_clipped": circular(**default_params, frequency=1.0, clip=True),
-        "rings_circular_2sided": circular_two_sided(**default_params, frequency=1.0),
+        "rings_circular_2sided": circular_two_sided(**default_params, frequency=1.0, intensity_rings=((0.0, 1.0), (1.0, 0.0)), clip=True, target_indices=2),
 
         "rings_rectangular": rectangular(**default_params, frequency=1.0),
         "rings_rectangular_with_targets": rectangular(**default_params, frequency=1.0, target_indices=3),
         "rings_rectangular_clipped": rectangular(**default_params, frequency=1.0, clip=True),
-        "rings_rectangular_2sided": rectangular_two_sided(**default_params, frequency=1.0),
+        "rings_rectangular_2sided": rectangular_two_sided(**default_params, frequency=1.0, intensity_frames=((0.0, 1.0), (1.0, 0.0)), clip=True, target_indices=2),
     }
     # fmt: on
 
diff --git a/stimupy/stimuli/sbcs.py b/stimupy/stimuli/sbcs.py
index 958bab27..c1110e19 100644
--- a/stimupy/stimuli/sbcs.py
+++ b/stimupy/stimuli/sbcs.py
@@ -1,15 +1,20 @@
 import numpy as np
 
 from stimupy.components.shapes import disc, rectangle
-from stimupy.utils import pad_by_visual_size, pad_to_shape, resolution, stack_dicts
+from stimupy.stimuli import bullseyes
+from stimupy.utils import make_two_sided, pad_by_visual_size, pad_to_shape, resolution
 
 __all__ = [
     "generalized",
     "basic",
-    "two_sided",
+    "square",
+    "circular",
     "with_dots",
-    "with_dots_two_sided",
     "dotted",
+    "basic_two_sided",
+    "square_two_sided",
+    "circular_two_sided",
+    "with_dots_two_sided",
     "dotted_two_sided",
 ]
 
@@ -22,26 +27,30 @@ def generalized(
     target_position=None,
     intensity_background=0.0,
     intensity_target=0.5,
+    rotation=0.0,
 ):
-    """Simultaneous contrast stimulus with free target placement.
+    """Simultaneous contrast stimulus with free target placement
 
     Parameters
     ----------
     visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, in degrees
+        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 grating, in pixels
+        shape [height, width] of image, in pixels
     target_size : float or (float, float)
-        size of the target in degree visual angle (height, width)
+        size [height, width] of the target, in degrees visual angle
     target_position : float or (float, float)
         position of the target in degree visual angle (height, width);
         if None, place target centrally
-    intensity_background : float
-        intensity value for background
-    intensity_target : float
-        intensity value for target
+    intensity_background : float, optional
+        intensity value for background, by default 0.0
+    intensity_target : float, optional
+        intensity value for target, by default 0.5
+    rotation : float, optional
+        rotation (in degrees), counterclockwise, by default 0.0 (horizonal)
+
 
     Returns
     -------
@@ -66,18 +75,14 @@ def generalized(
         rectangle_position=target_position,
         intensity_background=intensity_background,
         intensity_rectangle=intensity_target,
+        rotation=rotation,
     )
 
-    stim["target_mask"] = stim["rectangle_mask"]
-    stim["target_size"] = stim["rectangle_size"]
-    stim["target_position"] = stim["rectangle_position"]
-    stim["intensity_target"] = stim["intensity_rectangle"]
-    del (
-        stim["rectangle_mask"],
-        stim["rectangle_size"],
-        stim["rectangle_position"],
-        stim["intensity_rectangle"],
-    )
+    stim["target_mask"] = stim.pop("rectangle_mask")
+    stim["target_size"] = stim.pop("rectangle_size")
+    stim["target_position"] = stim.pop("rectangle_position")
+    stim["intensity_target"] = stim.pop("intensity_rectangle")
+
     return stim
 
 
@@ -89,22 +94,22 @@ def basic(
     intensity_background=0.0,
     intensity_target=0.5,
 ):
-    """Simultaneous contrast stimulus with central target.
+    """Simultaneous contrast stimulus with central target
 
     Parameters
     ----------
     visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, in degrees
+        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 grating, in pixels
+        shape [height, width] of image, in pixels
     target_size : float or (float, float)
-        size of the target in degree visual angle (height, width)
-    intensity_background : float
-        intensity value for background
-    intensity_target : float
-        intensity value for target
+        size [height, width] of the target, in degrees visual angle
+    intensity_background : float, optional
+        intensity value for background, by default 0.0
+    intensity_target : float, optional
+        intensity value for target, by default 0.5
 
     Returns
     -------
@@ -124,6 +129,7 @@ def basic(
     stim = generalized(
         visual_size=visual_size,
         ppd=ppd,
+        shape=shape,
         target_size=target_size,
         target_position=None,
         intensity_background=intensity_background,
@@ -132,72 +138,123 @@ def basic(
     return stim
 
 
-def two_sided(
+def square(
+    target_radius,
     visual_size=None,
     ppd=None,
     shape=None,
-    target_size=None,
-    intensity_backgrounds=(0.0, 1.0),
+    surround_radius=None,
+    rotation=0.0,
+    intensity_surround=0.0,
+    intensity_background=0.5,
     intensity_target=0.5,
+    origin="mean",
 ):
-    """Two-sided simultaneous contrast display with central targets.
+    """Simultaneous contrast stimulus with square target and surround field
 
     Parameters
     ----------
     visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, in degrees
+        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 grating, in pixels
-    target_size : float or (float, float)
-        size of the target in degree visual angle (height, width)
-    intensity_background : Sequence[Number, Number]
-        intensity values for backgrounds
-    intensity_target : float
-        intensity value for target
+        shape [height, width] of image, in pixels
+    target_radius : float
+        radius of target, in degrees visual angle
+    surround_radius : float
+        radius of surround context field, in degrees visual angle
+    rotation : float, optional
+        rotation (in degrees), counterclockwise, by default 0.0 (horizonal)
+    intensity_surrond : float, optional
+        intensity of surround context field, by default 0.0
+    intensity_background : float, optional
+        intensity value of background, by default 0.5
+    intensity_target : float, or Sequence[float, ...], optional
+        intensity value for each target, by default 0.5.
+    origin : "corner", "mean" or "center"
+        if "corner": set origin to upper left corner
+        if "mean": set origin to hypothetical image center (default)
+        if "center": set origin to real center (closest existing value to mean)
 
     Returns
     -------
     dict[str, Any]
         dict with the stimulus (key: "img"),
-        mask with integer index for the target (key: "target_mask"),
+        mask with integer index for each frame (key: "target_mask"),
         and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Chevreul, M. (1855).
-        The principle of harmony and contrast of colors.
     """
-    if target_size is None:
-        raise ValueError("two_sided() missing argument 'target_size' which is not 'None'")
-
-    # Resolve resolution
-    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-
-    stim1 = basic(
-        visual_size=(visual_size[0], visual_size[1] / 2),
+    stim = bullseyes.rectangular_generalized(
+        radii=(target_radius, surround_radius),
+        visual_size=visual_size,
         ppd=ppd,
-        target_size=target_size,
-        intensity_background=intensity_backgrounds[0],
+        shape=shape,
+        rotation=rotation,
+        intensity_frames=(0.0, intensity_surround),
+        intensity_background=intensity_background,
         intensity_target=intensity_target,
+        origin=origin,
     )
+    return stim
 
-    stim2 = basic(
-        visual_size=(visual_size[0], visual_size[1] / 2),
+
+def circular(
+    target_radius,
+    visual_size=None,
+    ppd=None,
+    shape=None,
+    surround_radius=None,
+    intensity_surround=0.0,
+    intensity_background=0.5,
+    intensity_target=0.5,
+    origin="mean",
+):
+    """Simultaneous contrast stimulus with circular target and surround field
+
+    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
+    target_radius : float
+        radius of target, in degrees visual angle
+    surround_radius : float
+        radius of surround context field, in degrees visual angle
+    rotation : float, optional
+        rotation (in degrees), counterclockwise, by default 0.0 (horizonal)
+    intensity_surrond : float, optional
+        intensity of surround context field, by default 0.0
+    intensity_background : float (optional)
+        intensity value of background, by default 0.5
+    intensity_target : float, or Sequence[float, ...], optional
+        intensity value for each target, by default 0.5.
+        Can specify as many intensities as number of target_indices;
+        If fewer intensities are passed than target_indices, cycles through intensities
+    origin : "corner", "mean" or "center"
+        if "corner": set origin to upper left corner
+        if "mean": set origin to hypothetical image center (default)
+        if "center": set origin to real center (closest existing value to mean)
+
+    Returns
+    -------
+    dict[str, Any]
+        dict with the stimulus (key: "img"),
+        mask with integer index for each frame (key: "target_mask"),
+        and additional keys containing stimulus parameters
+    """
+    stim = bullseyes.circular_generalized(
+        radii=(target_radius, surround_radius),
+        visual_size=visual_size,
         ppd=ppd,
-        target_size=target_size,
-        intensity_background=intensity_backgrounds[1],
+        shape=shape,
+        intensity_rings=(0.0, intensity_surround),
+        intensity_background=intensity_background,
         intensity_target=intensity_target,
+        origin=origin,
     )
-
-    stim = stack_dicts(stim1, stim2)
-    del stim["intensity_background"]
-    del stim["target_position"]
-    stim["intensity_backgrounds"] = intensity_backgrounds
-    stim["target_positions"] = (stim1["target_position"], stim2["target_position"])
-    stim["shape"] = shape
-    stim["visual_size"] = visual_size
     return stim
 
 
@@ -333,96 +390,6 @@ def with_dots(
     return stim
 
 
-def with_dots_two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    n_dots=None,
-    dot_radius=None,
-    distance=None,
-    target_shape=None,
-    intensity_backgrounds=(0.0, 1.0),
-    intensity_dots=(1.0, 0.0),
-    intensity_target=0.5,
-):
-    """Two-sided simultaneous contrast stimulus with dots
-
-    Parameters
-    ----------
-    visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, 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 grating, in pixels
-    n_dots : int or (int, int)
-        stimulus size defined as the number of dots in y and x-directions
-    dot_radius : float
-        radius of dots
-    distance : float
-        distance between dots in degree visual angle
-    target_shape : int or (int, int)
-        target shape defined as the number of dots that fit into the target
-    intensity_backgrounds : Sequence[Number, Number]
-        intensity values for background
-    intensity_dots : Sequence[Number, Number]
-        intensity values for dots
-    intensity_target : float
-        intensity value for target
-
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for the target (key: "target_mask"),
-        and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Bressan, P., & Kramer, P. (2008).
-        Gating of remote effects on lightness.
-        Journal of Vision, 8(2), 16-16.
-        https://doi.org/10.1167/8.2.16
-    """
-
-    # Resolve resolution
-    try:
-        shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-        visual_size_ = (visual_size[0], visual_size[1] / 2)
-    except resolution.TooManyUnknownsError:
-        visual_size_ = None
-
-    stim1 = with_dots(
-        visual_size=visual_size_,
-        ppd=ppd,
-        n_dots=n_dots,
-        dot_radius=dot_radius,
-        distance=distance,
-        target_shape=target_shape,
-        intensity_background=intensity_backgrounds[0],
-        intensity_dots=intensity_dots[0],
-        intensity_target=intensity_target,
-    )
-
-    stim2 = with_dots(
-        visual_size=visual_size_,
-        ppd=ppd,
-        n_dots=n_dots,
-        dot_radius=dot_radius,
-        distance=distance,
-        target_shape=target_shape,
-        intensity_background=intensity_backgrounds[1],
-        intensity_dots=intensity_dots[1],
-        intensity_target=intensity_target,
-    )
-
-    stim = stack_dicts(stim1, stim2)
-    del stim["intensity_background"]
-    stim["intensity_backgrounds"] = intensity_backgrounds
-    stim["intensity_dots"] = intensity_dots
-    return stim
-
-
 def dotted(
     visual_size=None,
     ppd=None,
@@ -555,94 +522,56 @@ def dotted(
     return stim
 
 
-def dotted_two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    n_dots=None,
-    dot_radius=None,
-    distance=None,
-    target_shape=None,
-    intensity_backgrounds=(0.0, 1.0),
-    intensity_dots=(1.0, 0.0),
-    intensity_target=0.5,
-):
-    """Two-sided dotted simultaneous contrast stimulus
+generalized_two_sided = make_two_sided(
+    generalized,
+    two_sided_params=(
+        "target_size",
+        "target_position",
+        "rotation",
+        "intensity_target",
+        "intensity_background",
+    ),
+)
 
-    Parameters
-    ----------
-    visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, 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 grating, in pixels
-    n_dots : int or (int, int)
-        stimulus size defined as the number of dots in y and x-directions
-    dot_radius : float
-        radius of dots
-    distance : float
-        distance between dots in degree visual angle
-    target_shape : int or (int, int)
-        target shape defined as the number of dots that fit into the target
-    intensity_background : equence[Number, Number]
-        intensity values for background
-    intensity_dots : equence[Number, Number]
-        intensity values for dots
-    intensity_target : float
-        intensity value for target
 
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for the targets (key: "target_mask"),
-        and additional keys containing stimulus parameters
+basic_two_sided = make_two_sided(
+    basic, two_sided_params=("target_size", "intensity_target", "intensity_background")
+)
 
-    References
-    ----------
-    Bressan, P., & Kramer, P. (2008).
-        Gating of remote effects on lightness.
-        Journal of Vision, 8(2), 16-16.
-        https://doi.org/10.1167/8.2.16
-    """
 
-    # Resolve resolution
-    try:
-        shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-        visual_size_ = (visual_size[0], visual_size[1] / 2)
-    except resolution.TooManyUnknownsError:
-        visual_size_ = None
+square_two_sided = make_two_sided(
+    square,
+    two_sided_params=(
+        "target_radius",
+        "surround_radius",
+        "rotation",
+        "intensity_target",
+        "intensity_surround",
+        "intensity_background",
+    ),
+)
 
-    stim1 = dotted(
-        visual_size=visual_size_,
-        ppd=ppd,
-        n_dots=n_dots,
-        dot_radius=dot_radius,
-        distance=distance,
-        target_shape=target_shape,
-        intensity_background=intensity_backgrounds[0],
-        intensity_dots=intensity_dots[0],
-        intensity_target=intensity_target,
-    )
 
-    stim2 = dotted(
-        visual_size=visual_size_,
-        ppd=ppd,
-        n_dots=n_dots,
-        dot_radius=dot_radius,
-        distance=distance,
-        target_shape=target_shape,
-        intensity_background=intensity_backgrounds[1],
-        intensity_dots=intensity_dots[1],
-        intensity_target=intensity_target,
-    )
+circular_two_sided = make_two_sided(
+    circular,
+    two_sided_params=(
+        "target_radius",
+        "surround_radius",
+        "intensity_target",
+        "intensity_surround",
+        "intensity_background",
+    ),
+)
 
-    stim = stack_dicts(stim1, stim2)
-    del stim["intensity_background"]
-    stim["intensity_backgrounds"] = intensity_backgrounds
-    stim["intensity_dots"] = intensity_dots
-    return stim
+
+with_dots_two_sided = make_two_sided(
+    with_dots, two_sided_params=("intensity_dots", "intensity_background", "intensity_target")
+)
+
+
+dotted_two_sided = make_two_sided(
+    dotted, two_sided_params=("intensity_target", "intensity_background", "intensity_dots")
+)
 
 
 def overview(**kwargs):
@@ -658,15 +587,27 @@ def overview(**kwargs):
     }
     default_params.update(kwargs)
 
+    dot_params = {
+        "n_dots": 5,
+        "dot_radius": 2,
+        "distance": 0.05,
+        "target_shape": 3,
+    }
+
     # fmt: off
     stimuli = {
-        "sbc_generalized": generalized(**default_params, visual_size=10, target_size=(3,4), target_position=(1, 2)),
+        "sbc_generalized": generalized(**default_params, visual_size=10, target_size=(3, 4), target_position=(1, 2)),
         "sbc_basic": basic(**default_params, visual_size=10, target_size=3),
-        "sbc_2sided": two_sided(**default_params, visual_size=10, target_size=2),
-        "sbc_with_dots": with_dots(**default_params, n_dots=5, dot_radius=2, distance=0.5, target_shape=3),
-        "sbc_with_dots_2sided": with_dots_two_sided(**default_params, n_dots=5, dot_radius=2, distance=0.5, target_shape=3),
-        "sbc_dotted": dotted(**default_params, n_dots=5, dot_radius=2, distance=0.5, target_shape=3),
-        "sbc_dotted_2sided": dotted_two_sided(**default_params, n_dots=5, dot_radius=2, distance=0.5, target_shape=3),
+        "sbc_square": square(**default_params, visual_size=10, target_radius=1, surround_radius=2),
+        "sbc_circular": circular(**default_params, visual_size=10, target_radius=1, surround_radius=2),
+        "sbc_with_dots": with_dots(**default_params, **dot_params),
+        "sbc_dotted": dotted(**default_params, **dot_params),
+
+        "sbc_2sided": basic_two_sided(**default_params, visual_size=10, target_size=2, intensity_background=(0.0, 1.0)),
+        "sbc_square_2sided": square_two_sided(**default_params, visual_size=10, target_radius=1, surround_radius=2, intensity_surround=(1.0, 0.0)),
+        "sbc_circular_2sided": circular_two_sided(**default_params, visual_size=10, target_radius=1, surround_radius=2, intensity_surround=(1.0, 0.0)),
+        "sbc_with_dots_2sided": with_dots_two_sided(**default_params, **dot_params, intensity_background=(0.0, 1.0), intensity_dots=(1.0, 0.0)),
+        "sbc_dotted_2sided": dotted_two_sided(**default_params, **dot_params, intensity_background=(0.0, 1.0), intensity_dots=(1.0, 0.0)),
     }
     # fmt: on
 
@@ -677,4 +618,4 @@ def overview(**kwargs):
     from stimupy.utils import plot_stimuli
 
     stims = overview()
-    plot_stimuli(stims, mask=True, save=None)
+    plot_stimuli(stims, mask=False, save=None)
diff --git a/stimupy/stimuli/todorovics.py b/stimupy/stimuli/todorovics.py
index fe1342a5..b45d6687 100644
--- a/stimupy/stimuli/todorovics.py
+++ b/stimupy/stimuli/todorovics.py
@@ -4,7 +4,7 @@
 
 from stimupy.components.shapes import cross as cross_shape
 from stimupy.components.shapes import rectangle as rectangle_shape
-from stimupy.utils import pad_dict_to_shape, resolution, stack_dicts
+from stimupy.utils import make_two_sided, pad_dict_to_shape, resolution
 
 __all__ = [
     "rectangle_generalized",
@@ -251,96 +251,17 @@ def rectangle(
     return stim
 
 
-def rectangle_two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    target_size=None,
-    covers_size=None,
-    covers_offset=None,
-    intensity_backgrounds=(0.0, 1.0),
-    intensity_target=0.5,
-    intensity_covers=(1.0, 0.0),
-):
-    """
-    Two-sided Todorovic's illusion with rectangular target in the center and four
-    rectangular covers added symmetrically around target center
-
-    Parameters
-    ----------
-    visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, 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 grating, in pixels
-    target_size : float or (float, float)
-        size of the target in degrees of visual angle (height, width)
-    covers_size : float or (float, float)
-        size of covers in degrees of visual angle (height, width)
-    covers_offset : float or (float, float)
-        distance from cover center to target center in (y, x)
-    intensity_background : Sequence[Number, Number]
-        intensity values for backgrounds
-    intensity_target : float
-        intensity value for target
-    intensity_covers : Sequence[Number, Number]
-        intensity values for covers
-
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for the target (key: "target_mask"),
-        and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Blakeslee, B., & McCourt, M. E. (1999).
-        A multiscale spatial ﬁltering account
-        of the White eﬀect, simultaneous brightness contrast and grating induction.
-        Vision Research, 39, 4361-4377.
-    Pessoa, L., Baratoff, G., Neumann, H., & Todorovic, D. (1998).
-        Lightness and junctions: variations on White's display.
-        Investigative Ophthalmology and Visual Science (Supplement), 39, S159.
-    Todorovic, D. (1997).
-        Lightness and junctions. Perception, 26, 379-395.
-    """
-
-    # Resolve resolution
-    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-
-    stim1 = rectangle(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        target_size=target_size,
-        covers_size=covers_size,
-        covers_offset=covers_offset,
-        intensity_background=intensity_backgrounds[0],
-        intensity_target=intensity_target,
-        intensity_covers=intensity_covers[0],
-    )
-
-    stim2 = rectangle(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        target_size=target_size,
-        covers_size=covers_size,
-        covers_offset=covers_offset,
-        intensity_background=intensity_backgrounds[1],
-        intensity_target=intensity_target,
-        intensity_covers=intensity_covers[1],
-    )
-
-    stim = stack_dicts(stim1, stim2)
-    del stim["intensity_background"]
-    del stim["target_position"]
-    stim["intensity_backgrounds"] = intensity_backgrounds
-    stim["intensity_covers"] = intensity_covers
-    stim["target_positions"] = (stim1["target_position"], stim2["target_position"])
-    stim["shape"] = shape
-    stim["visual_size"] = visual_size
-    return stim
+rectangle_two_sided = make_two_sided(
+    rectangle,
+    two_sided_params=(
+        "target_size",
+        "covers_size",
+        "covers_offset",
+        "intensity_background",
+        "intensity_target",
+        "intensity_covers",
+    ),
+)
 
 
 def cross_generalized(
@@ -584,92 +505,17 @@ def cross(
     return stim
 
 
-def cross_two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    cross_size=None,
-    cross_thickness=None,
-    covers_size=None,
-    intensity_backgrounds=(0.0, 1.0),
-    intensity_target=0.5,
-    intensity_covers=(1.0, 0.0),
-):
-    """Two-sided with cross target and four rectangular covers added at inner cross corners
-
-    Parameters
-    ----------
-    visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, 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 grating, in pixels
-    cross_size : float or (float, float)
-        size of target cross in visual angle
-    cross_thickness : float
-        thickness of target cross in visual angle
-    covers_size : float or (float, float)
-        size of covers in degrees of visual angle (height, width)
-    intensity_backgrounds : Sequence[Number, Number]
-        intensity values for backgrounds
-    intensity_target : float
-        intensity value for target
-    intensity_covers : Sequence[Number, Number]
-        intensity values for covers
-
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for the target (key: "target_mask"),
-        and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Blakeslee, B., & McCourt, M. E. (1999).
-        A multiscale spatial ﬁltering account
-        of the White eﬀect, simultaneous brightness contrast and grating induction.
-        Vision Research, 39, 4361-4377.
-    Pessoa, L., Baratoff, G., Neumann, H., & Todorovic, D. (1998).
-        Lightness and junctions: variations on White's display.
-        Investigative Ophthalmology and Visual Science (Supplement), 39, S159.
-    Todorovic, D. (1997).
-        Lightness and junctions. Perception, 26, 379-395.
-    """
-
-    # Resolve resolution
-    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-
-    stim1 = cross(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        cross_size=cross_size,
-        cross_thickness=cross_thickness,
-        covers_size=covers_size,
-        intensity_background=intensity_backgrounds[0],
-        intensity_target=intensity_target,
-        intensity_covers=intensity_covers[0],
-    )
-
-    stim2 = cross(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        cross_size=cross_size,
-        cross_thickness=cross_thickness,
-        covers_size=covers_size,
-        intensity_background=intensity_backgrounds[1],
-        intensity_target=intensity_target,
-        intensity_covers=intensity_covers[1],
-    )
-
-    stim = stack_dicts(stim1, stim2)
-    del stim["intensity_background"]
-    stim["intensity_backgrounds"] = intensity_backgrounds
-    stim["intensity_covers"] = intensity_covers
-    stim["shape"] = shape
-    stim["visual_size"] = visual_size
-    return stim
+cross_two_sided = make_two_sided(
+    cross,
+    two_sided_params=(
+        "cross_size",
+        "cross_thickness",
+        "covers_size",
+        "intensity_background",
+        "intensity_target",
+        "intensity_covers",
+    ),
+)
 
 
 def equal(
@@ -767,87 +613,16 @@ def equal(
     return stim
 
 
-def equal_two_sided(
-    visual_size=None,
-    ppd=None,
-    shape=None,
-    cross_size=None,
-    cross_thickness=None,
-    intensity_backgrounds=(0.0, 1.0),
-    intensity_target=0.5,
-    intensity_covers=(1.0, 0.0),
-):
-    """Two-sided with cross target and four rectangular covers added at inner cross corners
-
-    Parameters
-    ----------
-    visual_size : Sequence[Number, Number], Number, or None (default)
-        visual size [height, width] of grating, 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 grating, in pixels
-    cross_size : float or (float, float)
-        size of target cross in visual angle
-    cross_thickness : float
-        thickness of target cross in visual angle
-    intensity_background : float
-        intensity value for background
-    intensity_target : float
-        intensity value for target
-    intensity_covers : float
-        intensity value for covers
-
-    Returns
-    -------
-    dict[str, Any]
-        dict with the stimulus (key: "img"),
-        mask with integer index for the target (key: "target_mask"),
-        and additional keys containing stimulus parameters
-
-    References
-    ----------
-    Blakeslee, B., & McCourt, M. E. (1999).
-        A multiscale spatial ﬁltering account
-        of the White eﬀect, simultaneous brightness contrast and grating induction.
-        Vision Research, 39, 4361-4377.
-    Pessoa, L., Baratoff, G., Neumann, H., & Todorovic, D. (1998).
-        Lightness and junctions: variations on White's display.
-        Investigative Ophthalmology and Visual Science (Supplement), 39, S159.
-    Todorovic, D. (1997).
-        Lightness and junctions. Perception, 26, 379-395.
-    """
-
-    # Resolve resolution
-    shape, visual_size, ppd = resolution.resolve(shape=shape, visual_size=visual_size, ppd=ppd)
-
-    stim1 = equal(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        cross_size=cross_size,
-        cross_thickness=cross_thickness,
-        intensity_background=intensity_backgrounds[0],
-        intensity_target=intensity_target,
-        intensity_covers=intensity_covers[0],
-    )
-
-    stim2 = equal(
-        visual_size=(visual_size[0], visual_size[1] / 2),
-        ppd=ppd,
-        cross_size=cross_size,
-        cross_thickness=cross_thickness,
-        intensity_background=intensity_backgrounds[1],
-        intensity_target=intensity_target,
-        intensity_covers=intensity_covers[1],
-    )
-
-    stim = stack_dicts(stim1, stim2)
-    del stim["intensity_background"]
-    stim["intensity_backgrounds"] = intensity_backgrounds
-    stim["intensity_covers"] = intensity_covers
-    stim["shape"] = shape
-    stim["visual_size"] = visual_size
-    return stim
+equal_two_sided = make_two_sided(
+    equal,
+    two_sided_params=(
+        "cross_size",
+        "cross_thickness",
+        "intensity_background",
+        "intensity_target",
+        "intensity_covers",
+    ),
+)
 
 
 def overview(**kwargs):
@@ -877,12 +652,12 @@ def overview(**kwargs):
     stimuli = {
         "todorovic_rectangle": rectangle(**default_params, **rectangle_params, covers_offset=1.5),
         "todorovic_rectangle_general": rectangle_generalized(**default_params, **rectangle_params, target_position=3.5, covers_x=(2, 6), covers_y=(2, 6)),
-        "todorovic_rectangle_2sided": rectangle_two_sided(**default_params, **rectangle_params, covers_offset=1),
+        "todorovic_rectangle_2sided": rectangle_two_sided(**default_params, **rectangle_params, covers_offset=1, intensity_background=(0.0, 1.0), intensity_covers=(1.0, 0.0)),
         "todorovic_cross": cross(**default_params, **cross_params, covers_size=2),
         "todorovic_cross_general": cross_generalized(**default_params, **cross_params, covers_size=2, covers_x=(2, 6), covers_y=(2, 6)),
-        "todorovic_cross_2sided": cross_two_sided(**default_params, **cross_params, covers_size=1),
+        "todorovic_cross_2sided": cross_two_sided(**default_params, **cross_params, covers_size=1, intensity_background=(0.0, 1.0), intensity_covers=(1.0, 0.0)),
         "todorovic_equal": equal(**default_params, **cross_params,),
-        "todorovic_equal_2sided": equal_two_sided(**default_params, **cross_params,),
+        "todorovic_equal_2sided": equal_two_sided(**default_params, **cross_params, intensity_background=(0.0, 1.0), intensity_covers=(1.0, 0.0)),
     }
     # fmt: on
 
diff --git a/stimupy/stimuli/waves.py b/stimupy/stimuli/waves.py
index a1f27a88..eda677c0 100644
--- a/stimupy/stimuli/waves.py
+++ b/stimupy/stimuli/waves.py
@@ -1,9 +1,8 @@
-import itertools
-
 import numpy as np
 
-from stimupy.components import draw_regions, waves
+from stimupy.components import waves
 from stimupy.components.shapes import disc, rectangle
+from stimupy.stimuli import place_targets
 
 __all__ = [
     "sine_linear",
@@ -21,29 +20,6 @@
 ]
 
 
-def add_targets(wave_stim, target_indices, intensity_target):
-    # Create target-mask
-    if isinstance(target_indices, (int)):
-        target_indices = [target_indices]
-
-    targets_mask = np.zeros_like(wave_stim["grating_mask"])
-    for target_idx, bar_idx in enumerate(target_indices):
-        targets_mask = np.where(
-            wave_stim["grating_mask"] == (bar_idx + 1), target_idx + 1, targets_mask
-        )
-    targets_mask = targets_mask.astype(int)
-    wave_stim["target_mask"] = targets_mask
-
-    # Place target(s)
-    if isinstance(intensity_target, (int, float)):
-        intensities = [intensity_target]
-        intensities = itertools.cycle(intensities)
-    for target_idx, intensity in zip(np.unique(targets_mask[targets_mask > 0]), intensities):
-        wave_stim["img"] = np.where(targets_mask == target_idx, intensity, wave_stim["img"])
-
-    return wave_stim
-
-
 def sine_linear(
     visual_size=None,
     ppd=None,
@@ -145,12 +121,16 @@ def sine_linear(
     stim["bar_width"] = stim.pop("phase_width")
     stim.pop("distance_metric")
 
-    # Add targets(?)
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
+    # Rename mask
+    stim["bar_mask"] = stim["grating_mask"]
 
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -164,7 +144,7 @@ def square_linear(
     period="ignore",
     rotation=0.0,
     phase_shift=0,
-    intensity_bars=(1.0, 0.0),
+    intensity_bars=(0.0, 1.0),
     target_indices=(),
     intensity_target=0.5,
     origin="corner",
@@ -239,23 +219,25 @@ def square_linear(
         origin=origin,
         round_phase_width=round_phase_width,
         distance_metric="oblique",
+        intensities=intensity_bars,
     )
 
-    # Adjust intensities to passed-in values
-    stim["img"] = draw_regions(mask=stim["grating_mask"], intensities=intensity_bars)
-
     # Repackage output
     stim["n_bars"] = stim.pop("n_phases")
     stim["bar_width"] = stim.pop("phase_width")
     stim["intensity_bars"] = stim.pop("intensities")
     stim.pop("distance_metric")
 
-    # Add targets(?)
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
+    # Rename mask
+    stim["bar_mask"] = stim["grating_mask"]
 
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -358,12 +340,16 @@ def staircase_linear(
     stim["intensity_bars"] = stim.pop("intensities")
     stim.pop("distance_metric")
 
-    # Add targets(?)
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
+    # Rename mask
+    stim["bar_mask"] = stim["grating_mask"]
 
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -478,15 +464,19 @@ def sine_radial(
         )
         stim["img"] = np.where(circle["ring_mask"], stim["img"], intensity_background)
         stim["grating_mask"] = np.where(circle["ring_mask"], stim["grating_mask"], 0)
-
-    # Resolve target parameters
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
-
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
     stim["clip"] = clip
     stim["intensity_background"] = intensity_background
+
+    # Rename mask
+    stim["ring_mask"] = stim["grating_mask"]
+
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -500,7 +490,7 @@ def square_radial(
     period="ignore",
     rotation=0.0,
     phase_shift=0,
-    intensity_rings=(1.0, 0.0),
+    intensity_rings=(0.0, 1.0),
     target_indices=(),
     intensity_target=0.5,
     origin="mean",
@@ -579,11 +569,9 @@ def square_radial(
         origin=origin,
         round_phase_width=round_phase_width,
         distance_metric="radial",
+        intensities=intensity_rings,
     )
 
-    # Adjust intensities to passed-in values
-    stim["img"] = draw_regions(mask=stim["grating_mask"], intensities=intensity_rings)
-
     # Repackage output
     stim["n_rings"] = stim.pop("n_phases")
     stim["ring_width"] = stim.pop("phase_width")
@@ -601,15 +589,19 @@ def square_radial(
         )
         stim["img"] = np.where(circle["ring_mask"], stim["img"], intensity_background)
         stim["grating_mask"] = np.where(circle["ring_mask"], stim["grating_mask"], 0)
-
-    # Resolve target parameters
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
-
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
     stim["clip"] = clip
     stim["intensity_background"] = intensity_background
+
+    # Rename mask
+    stim["ring_mask"] = stim["grating_mask"]
+
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -623,7 +615,7 @@ def staircase_radial(
     period="ignore",
     rotation=0.0,
     phase_shift=0,
-    intensity_rings=(1.0, 0.0),
+    intensity_rings=(0.0, 1.0),
     target_indices=(),
     intensity_target=0.5,
     origin="mean",
@@ -725,15 +717,19 @@ def staircase_radial(
         )
         stim["img"] = np.where(circle["ring_mask"], stim["img"], intensity_background)
         stim["grating_mask"] = np.where(circle["ring_mask"], stim["grating_mask"], 0)
-
-    # Resolve target parameters
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
-
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
     stim["clip"] = clip
     stim["intensity_background"] = intensity_background
+
+    # Rename mask
+    stim["ring_mask"] = stim["grating_mask"]
+
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -856,15 +852,19 @@ def sine_rectilinear(
             )
         stim["img"] = np.where(rect["rectangle_mask"], stim["img"], intensity_background)
         stim["grating_mask"] = np.where(rect["rectangle_mask"], stim["grating_mask"], 0)
-
-    # Add targets(?)
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
-
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
     stim["clip"] = clip
     stim["intensity_background"] = intensity_background
+
+    # Rename mask
+    stim["frame_mask"] = stim["grating_mask"]
+
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -878,7 +878,7 @@ def square_rectilinear(
     period="ignore",
     rotation=0.0,
     phase_shift=0,
-    intensity_frames=(1.0, 0.0),
+    intensity_frames=(0.0, 1.0),
     target_indices=(),
     intensity_target=0.5,
     origin="mean",
@@ -958,11 +958,9 @@ def square_rectilinear(
         origin=origin,
         round_phase_width=round_phase_width,
         distance_metric="rectilinear",
+        intensities=intensity_frames,
     )
 
-    # Adjust intensities to passed-in values
-    stim["img"] = draw_regions(mask=stim["grating_mask"], intensities=intensity_frames)
-
     # Repackage output
     stim["n_frames"] = stim.pop("n_phases")
     stim["frame_width"] = stim.pop("phase_width")
@@ -988,15 +986,19 @@ def square_rectilinear(
             )
         stim["img"] = np.where(rect["rectangle_mask"], stim["img"], intensity_background)
         stim["grating_mask"] = np.where(rect["rectangle_mask"], stim["grating_mask"], 0)
-
-    # Add targets(?)
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
-
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
     stim["clip"] = clip
     stim["intensity_background"] = intensity_background
+
+    # Rename mask
+    stim["frame_mask"] = stim["grating_mask"]
+
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -1020,8 +1022,6 @@ def staircase_rectilinear(
 ):
     """Rectiinear staircase, with some frame(s) as target(s)
 
-    Parameters
-    ----------
     Parameters
     ----------
     visual_size : Sequence[Number, Number], Number, or None (default)
@@ -1122,15 +1122,19 @@ def staircase_rectilinear(
             )
         stim["img"] = np.where(rect["rectangle_mask"], stim["img"], intensity_background)
         stim["grating_mask"] = np.where(rect["rectangle_mask"], stim["grating_mask"], 0)
-
-    # Add targets(?)
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
-
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
     stim["clip"] = clip
     stim["intensity_background"] = intensity_background
+
+    # Rename mask
+    stim["frame_mask"] = stim["grating_mask"]
+
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -1227,32 +1231,17 @@ def sine_angular(
     stim["n_segments"] = stim.pop("n_phases")
     stim["segment_width"] = stim.pop("phase_width")
     stim.pop("distance_metric")
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
-
-    # Resolve target parameters
-    if target_indices is not None and target_indices != ():
-        if isinstance(target_indices, (int)):
-            target_indices = [
-                target_indices,
-            ]
-        if isinstance(intensity_target, (int, float)):
-            intensity_target = [
-                intensity_target,
-            ]
-        intensity_target = itertools.cycle(intensity_target)
-
-        # Place target(s)
-        targets_mask = np.zeros_like(stim["grating_mask"])
-        for target_idx, (segment_idx, intensity) in enumerate(
-            zip(target_indices, intensity_target)
-        ):
-            targets_mask = np.where(
-                stim["grating_mask"] == (segment_idx + 1), target_idx + 1, targets_mask
-            )
-            stim["img"] = np.where(targets_mask == (target_idx + 1), intensity, stim["img"])
-        stim["target_mask"] = targets_mask.astype(int)
 
+    # Rename mask
+    stim["segment_mask"] = stim["grating_mask"]
+
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -1266,7 +1255,7 @@ def square_angular(
     period="ignore",
     rotation=0.0,
     phase_shift=0,
-    intensity_segments=(1.0, 0.0),
+    intensity_segments=(0.0, 1.0),
     target_indices=(),
     intensity_target=0.5,
     origin="mean",
@@ -1341,42 +1330,25 @@ def square_angular(
         origin=origin,
         round_phase_width=round_phase_width,
         distance_metric="angular",
+        intensities=intensity_segments,
     )
 
-    # Adjust intensities to passed-in values
-    stim["img"] = draw_regions(mask=stim["grating_mask"], intensities=intensity_segments)
-
     # Repackage output
     stim["n_segments"] = stim.pop("n_phases")
     stim["segment_width"] = stim.pop("phase_width")
     stim["intensity_segments"] = stim.pop("intensities")
     stim.pop("distance_metric")
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
-
-    # Resolve target parameters
-    if target_indices is not None and target_indices != ():
-        if isinstance(target_indices, (int)):
-            target_indices = [
-                target_indices,
-            ]
-        if isinstance(intensity_target, (int, float)):
-            intensity_target = [
-                intensity_target,
-            ]
-        intensity_target = itertools.cycle(intensity_target)
-
-        # Place target(s)
-        targets_mask = np.zeros_like(stim["grating_mask"])
-        for target_idx, (segment_idx, intensity) in enumerate(
-            zip(target_indices, intensity_target)
-        ):
-            targets_mask = np.where(
-                stim["grating_mask"] == (segment_idx + 1), target_idx + 1, targets_mask
-            )
-            stim["img"] = np.where(targets_mask == (target_idx + 1), intensity, stim["img"])
-        stim["target_mask"] = targets_mask.astype(int)
 
+    # Rename mask
+    stim["segment_mask"] = stim["grating_mask"]
+
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
@@ -1466,6 +1438,7 @@ def staircase_angular(
         origin=origin,
         round_phase_width=round_phase_width,
         distance_metric="angular",
+        intensities=intensity_segments,
     )
 
     # Repackage output
@@ -1474,12 +1447,16 @@ def staircase_angular(
     stim["intensity_segments"] = stim.pop("intensities")
     stim.pop("distance_metric")
 
-    # Add targets(?)
-    if target_indices is not None and target_indices != ():
-        stim = add_targets(stim, target_indices=target_indices, intensity_target=intensity_target)
+    # Rename mask
+    stim["segment_mask"] = stim["grating_mask"]
 
-    stim["target_indices"] = target_indices
-    stim["intensity_target"] = intensity_target
+    # Add targets
+    stim = place_targets(
+        stim=stim,
+        element_mask_key="grating_mask",
+        target_indices=target_indices,
+        intensity_target=intensity_target,
+    )
     return stim
 
 
diff --git a/stimupy/stimuli/wedding_cakes.py b/stimupy/stimuli/wedding_cakes.py
index 813c0184..15c9a391 100644
--- a/stimupy/stimuli/wedding_cakes.py
+++ b/stimupy/stimuli/wedding_cakes.py
@@ -16,7 +16,7 @@ def wedding_cake(
     target_height=None,
     target_indices1=None,
     target_indices2=None,
-    intensity_bars=(1.0, 0.0),
+    intensity_bars=(0.0, 1.0),
     intensity_target=0.5,
 ):
     """Wedding cake stimulus
@@ -195,7 +195,7 @@ def overview(**kwargs):
     # fmt: off
     stimuli = {
         "wedding_cake": wedding_cake(**default_params,
-                                     L_size= (4, 3, 1),
+                                     L_size=(4, 3, 1),
                                      target_height=1,
                                      target_indices1=((2, 2), (2, 1)),
                                      target_indices2=((2, -1), (2, 0)),
diff --git a/stimupy/stimuli/whites.py b/stimupy/stimuli/whites.py
index daca172f..4683ce9e 100644
--- a/stimupy/stimuli/whites.py
+++ b/stimupy/stimuli/whites.py
@@ -4,6 +4,7 @@
 
 from stimupy.components import combine_masks, draw_regions
 from stimupy.components.shapes import rectangle
+from stimupy.stimuli import mask_targets
 from stimupy.stimuli.gratings import squarewave
 from stimupy.stimuli.pinwheels import pinwheel as angular
 from stimupy.stimuli.waves import square_radial as radial
@@ -113,34 +114,29 @@ def generalized(
         round_phase_width=round_phase_width,
     )
 
-    # Resolve target parameters
-    if isinstance(target_indices, (int)):
-        target_indices = (target_indices,)
-    if isinstance(intensity_target, (int, float)):
-        intensity_target = (intensity_target,)
+    # Mask target bars
+    target_bar_mask = mask_targets(
+        element_mask=stim["grating_mask"], target_indices=target_indices
+    )
+    stim["target_indices"] = target_indices
+
+    # Mask rectangular regions
     if isinstance(target_heights, (int, float)):
         target_heights = (target_heights,)
     if isinstance(target_center_offsets, (int, float)):
         target_center_offsets = (target_center_offsets,)
-
     if len(target_indices) != 0 and target_heights is None:
         raise ValueError("generalized() missing argument 'target_heights' which is not 'None'")
     if len(target_indices) == 0 and target_heights is None:
         target_heights = (0,)
-
-    intensity_target = tuple(
-        itertools.islice(itertools.cycle(intensity_target), len(target_indices))
-    )
     target_heights = tuple(itertools.islice(itertools.cycle(target_heights), len(target_indices)))
     target_center_offsets = tuple(
         itertools.islice(itertools.cycle(target_center_offsets), len(target_indices))
     )
 
-    # Place target(s)
     stim_center = stim["visual_size"].height / 2
-    target_zip = zip(target_indices, intensity_target, target_heights, target_center_offsets)
-    targets_mask = np.zeros_like(stim["grating_mask"])
-    for target_idx, (bar_idx, intensity, height, offset) in enumerate(target_zip):
+    target_rect_masks = []
+    for target_idx, (height, offset) in enumerate(zip(target_heights, target_center_offsets)):
         # Draw a stripe of target_height x stim_width, at center + offset
         rect = rectangle(
             visual_size=stim["visual_size"],
@@ -148,24 +144,34 @@ def generalized(
             shape=stim["shape"],
             rectangle_size=(height, stim["visual_size"].width),
             rectangle_position=(stim_center + offset - (height / 2), 0),
-            intensity_rectangle=intensity,
+            intensity_rectangle=target_idx,
         )
+        target_rect_masks.append(rect["rectangle_mask"])
+    stim["target_heights"] = target_heights
+    stim["target_center_offsets"] = target_center_offsets
 
+    # Combine rect & bar masks
+    target_masks = []
+    for target_idx, rect_mask in enumerate(target_rect_masks):
         # Find where strip intersects with the target bar
-        if bar_idx < 0:
-            bar_idx = int(stim["n_bars"]) + bar_idx
-        mask = (stim["grating_mask"] == bar_idx + 1) & rect["rectangle_mask"]
-        targets_mask = np.where(mask, target_idx + 1, targets_mask)
-
-        # Draw target
-        stim["img"] = np.where(targets_mask == target_idx + 1, intensity, stim["img"])
-
-    # Update and return stimulus
-    stim["target_mask"] = targets_mask.astype(int)
-    stim["target_indices"] = target_indices
+        target_mask = (target_bar_mask == target_idx + 1) & rect_mask
+        target_masks.append(target_mask)
+
+    # Combine masks
+    if len(target_masks) > 0:
+        target_mask = combine_masks(*target_masks)
+    else:
+        target_mask = np.zeros_like(stim["img"])
+    stim["target_mask"] = target_mask.astype(int)
+
+    # Draw targets
+    stim["img"] = np.where(
+        target_mask,
+        draw_regions(mask=target_mask, intensities=intensity_target, intensity_background=0.0),
+        stim["img"],
+    )
     stim["intensity_target"] = intensity_target
-    stim["target_heights"] = target_heights
-    stim["target_center_offsets"] = target_center_offsets
+
     return stim
 
 
@@ -178,7 +184,7 @@ def white(
     bar_width=None,
     period="ignore",
     rotation=0.0,
-    intensity_bars=(1.0, 0.0),
+    intensity_bars=(0.0, 1.0),
     target_indices=(),
     intensity_target=0.5,
     target_heights=None,
@@ -272,7 +278,7 @@ def white_two_rows(
     bar_width=None,
     period="ignore",
     rotation=0.0,
-    intensity_bars=(1.0, 0.0),
+    intensity_bars=(0.0, 1.0),
     intensity_target=0.5,
     target_indices_top=(),
     target_indices_bottom=(),
@@ -388,7 +394,7 @@ def anderson(
     target_indices_bottom=None,
     target_center_offset=0,
     target_height=None,
-    intensity_stripes=(1.0, 0.0),
+    intensity_stripes=(0.0, 1.0),
     stripe_center_offset=0,
     stripe_height=None,
     round_phase_width=True,
@@ -492,7 +498,7 @@ def anderson(
         if bar_idx < 0:
             bar_idx = int(stim["n_bars"]) + bar_idx
         stripe_top["rectangle_mask"] = np.where(
-            stim["grating_mask"] == bar_idx + 1, 0, stripe_top["rectangle_mask"]
+            stim["grating_mask"] == bar_idx, 0, stripe_top["rectangle_mask"]
         )
 
     stripe_bottom = rectangle(
@@ -506,7 +512,7 @@ def anderson(
         if bar_idx < 0:
             bar_idx = int(stim["n_bars"]) + bar_idx
         stripe_bottom["rectangle_mask"] = np.where(
-            stim["grating_mask"] == bar_idx + 1, 0, stripe_bottom["rectangle_mask"]
+            stim["grating_mask"] == bar_idx, 0, stripe_bottom["rectangle_mask"]
         )
 
     try:
@@ -543,7 +549,7 @@ def howe(
     target_indices_bottom=None,
     target_center_offset=0,
     target_height=None,
-    intensity_stripes=(1.0, 0.0),
+    intensity_stripes=(0.0, 1.0),
     round_phase_width=True,
 ):
     """Howe variation of White's stimulus
@@ -637,7 +643,7 @@ def yazdanbakhsh(
     target_indices_bottom=None,
     target_center_offset=0,
     target_heights=None,
-    # intensity_stripes=(1.0, 0.0),
+    # intensity_stripes=(0.0, 1.0),
     gap_size=None,
     round_phase_width=True,
 ):
@@ -734,7 +740,7 @@ def yazdanbakhsh(
         # Reduce to just this bar: intersection between rect mask and bar mask
         if bar_idx < 0:
             bar_idx = int(stim["n_bars"]) + bar_idx
-        stim["gap_mask"] = (stim["grating_mask"] == bar_idx + 1) & rect["rectangle_mask"]
+        stim["gap_mask"] = (stim["grating_mask"] == bar_idx) & rect["rectangle_mask"]
 
         # Remove everywhere it intersects with target mask
         stim["gap_mask"] = np.where(stim["target_mask"] == t_idx + 1, 0, stim["gap_mask"])
@@ -758,7 +764,7 @@ def yazdanbakhsh(
         # Reduce to just this bar: intersection between rect mask and bar mask
         if bar_idx < 0:
             bar_idx = int(stim["n_bars"]) + bar_idx
-        stim["gap_mask"] = (stim["grating_mask"] == bar_idx + 1) & rect["rectangle_mask"]
+        stim["gap_mask"] = (stim["grating_mask"] == bar_idx) & rect["rectangle_mask"]
 
         # Remove everywhere it intersects with target mask
         stim["gap_mask"] = np.where(stim["target_mask"] == t_idx + 1, 0, stim["gap_mask"])
diff --git a/stimupy/utils/export.py b/stimupy/utils/export.py
index a7870aa6..0d83d794 100644
--- a/stimupy/utils/export.py
+++ b/stimupy/utils/export.py
@@ -2,6 +2,7 @@
 import json
 import pickle
 from hashlib import md5
+from pathlib import Path
 
 import numpy as np
 from PIL import Image
@@ -21,48 +22,43 @@
 
 
 def array_to_checksum(arr):
-    """Hash (md5) values, and save only the hex
+    """Hash (md5) array values, and return hex-checksum
 
     Parameters
     ----------
-    arr : np.ndarray
-        Array to be hashed.
+    arr : numpy.ndarray
+        array to be hashed.
 
     Returns
     ----------
-    hex
+    str
+        hex-string representation of hash (MD5) of given array
     """
     return md5(np.ascontiguousarray(arr.round(8))).hexdigest()
 
 
-def array_to_image(arr, filename, norm=True):
+def array_to_image(arr, filename, format=None, norm=True):
     """Save a 2D numpy array as a grayscale image file.
 
     Parameters
     ----------
-    arr : np.ndarray
-        Array to be saved. Values will be cropped to [0,255].
-    filename : str
-        full path to the file to be created.
+    arr : numpy.ndarray
+        array to be exported. Values will be cropped to [0,255].
+    filename : Path or str
+        (full) path to the file to be created.
     norm : bool
-        if True (default), multiply array by 255
+        multiply array by 255, by default True
     """
-    if filename[-4:] != ".png" and filename[-4:] != ".jpg":
-        filename += ".png"
+    try:
+        arr_to_write = np.ascontiguousarray(arr)
+    except TypeError as e:
+        raise e from ValueError("arr should be a numpy.ndarray(-like)")
 
-    if isinstance(arr, (np.ndarray, list)):
-        arr = np.array(arr)
+    if norm:
+        arr_to_write = arr_to_write * 255
 
-        if norm:
-            arr = arr * 255
-
-        if Image:
-            imsize = arr.shape
-            im = Image.new("L", (imsize[1], imsize[0]))
-            im.putdata(arr.flatten())
-            im.save(filename)
-    else:
-        raise ValueError("arr should be a np.ndarray")
+    im = Image.fromarray(arr_to_write.astype("uint8"), mode="L")
+    im.save(filename)
 
 
 def array_to_npy(arr, filename):
@@ -70,15 +66,19 @@ def array_to_npy(arr, filename):
 
     Parameters
     ----------
-    arr : np.ndarray
-        Array to be saved.
-    filename : str
-        full path to the file to be creaated.
+    arr : numpy.ndarray
+        array to be exported.
+    filename : Path or str
+        (full) path to the file to be created.
     """
-    if isinstance(arr, (np.ndarray, list)):
-        np.save(filename, arr)
-    else:
-        raise ValueError("arr should be a np.ndarray")
+    try:
+        arr_to_write = np.ascontiguousarray(arr)
+    except TypeError as e:
+        raise e from ValueError("arr should be a numpy.ndarray(-like)")
+
+    filepath = Path(filename).resolve().with_suffix(".npy")
+
+    np.save(filepath, arr_to_write)
 
 
 def array_to_mat(arr, filename):
@@ -86,18 +86,19 @@ def array_to_mat(arr, filename):
 
     Parameters
     ----------
-    arr : np.ndarray
-        Array to be saved.
-    filename : str
-        full path to the file to be creaated.
+    arr : numpy.ndarray
+        array to be exported.
+    filename : Path or str
+        (full) path to the file to be created.
     """
-    if filename[-4:] != ".mat":
-        filename += ".mat"
+    try:
+        arr_to_write = np.ascontiguousarray(arr)
+    except TypeError as e:
+        raise e from ValueError("arr should be a numpy.ndarray(-like)")
+
+    filepath = Path(filename).resolve().with_suffix(".mat")
 
-    if isinstance(arr, (np.ndarray, list)):
-        savemat(filename, {"arr": arr})
-    else:
-        raise ValueError("arr should be a np.ndarray")
+    savemat(filepath, {"arr": arr_to_write})
 
 
 def array_to_pickle(arr, filename):
@@ -105,19 +106,20 @@ def array_to_pickle(arr, filename):
 
     Parameters
     ----------
-    arr : np.ndarray
-        Array to be saved.
-    filename : str
-        full path to the file to be creaated.
+    arr : numpy.ndarray
+        array to be exported.
+    filename : Path or str
+        (full) path to the file to be created.
     """
-    if filename[-7:] != ".pickle":
-        filename += ".pickle"
+    try:
+        arr_to_write = np.ascontiguousarray(arr)
+    except TypeError as e:
+        raise e from ValueError("arr should be a numpy.ndarray(-like)")
 
-    if isinstance(arr, (np.ndarray, list)):
-        with open(filename, "wb") as handle:
-            pickle.dump({"arr": arr}, handle, protocol=pickle.HIGHEST_PROTOCOL)
-    else:
-        raise ValueError("arr should be a np.ndarray")
+    filepath = Path(filename).resolve().with_suffix(".pickle")
+
+    with filepath.open("wb") as file:
+        pickle.dump({"arr": arr_to_write}, file, protocol=pickle.HIGHEST_PROTOCOL)
 
 
 def array_to_json(arr, filename):
@@ -125,22 +127,20 @@ def array_to_json(arr, filename):
 
     Parameters
     ----------
-    arr : np.ndarray
-        Array to be saved.
-    filename : str
-        full path to the file to be creaated.
+    arr : numpy.ndarray
+        array to be exported.
+    filename : Path or str
+        (full) path to the file to be created.
     """
-    if filename[-5:] != ".json":
-        filename += ".json"
+    try:
+        arr_to_write = np.ascontiguousarray(arr)
+    except TypeError as e:
+        raise e from ValueError("arr should be a numpy.ndarray(-like)")
+
+    filepath = Path(filename).resolve().with_suffix(".json")
 
-    if isinstance(arr, np.ndarray):
-        with open(filename, "w", encoding="utf-8") as f:
-            json.dump(arr.tolist(), f, ensure_ascii=False, indent=4)
-    elif isinstance(arr, list):
-        with open(filename, "w", encoding="utf-8") as f:
-            json.dump(arr, f, ensure_ascii=False, indent=4)
-    else:
-        raise ValueError("arr should be a np.ndarray")
+    with filepath.open("w", encoding="utf-8") as file:
+        json.dump(arr_to_write.tolist(), file, ensure_ascii=False, indent=4)
 
 
 def arrays_to_checksum(stim, keys=["img", "mask"]):
@@ -149,7 +149,7 @@ def arrays_to_checksum(stim, keys=["img", "mask"]):
     Parameters
     ----------
     stim : dict
-        stimulus dictionary containing keys
+        stimulus dictionary to export.
     keys : str of list of str
         keys of dict for which the hashing should be performed
 
@@ -175,27 +175,24 @@ def to_json(stim, filename):
     Parameters
     ----------
     stim : dict
-        stimulus dictionary containing keys
-    filename : str
-        full path to the file to be creaated.
+        stimulus dictionary to export.
+    filename : Path or str
+        (full) path to the file to be created.
 
     """
-    if filename[-5:] != ".json":
-        filename += ".json"
+    if not isinstance(stim, dict):
+        raise ValueError("stim should be a dict")
 
-    # stimulus-dict(s) as (pretty) JSON
-    if isinstance(stim, dict):
-        stim2 = copy.deepcopy(stim)
+    filepath = Path(filename).resolve().with_suffix(".json")
 
-        for key in stim2.keys():
-            # np.ndarrays are not serializable; change to list
-            if isinstance(stim2[key], np.ndarray):
-                stim2[key] = stim2[key].tolist()
+    stim2 = copy.deepcopy(stim)
+    for key in stim2.keys():
+        # np.ndarrays are not serializable; change to list
+        if isinstance(stim2[key], np.ndarray):
+            stim2[key] = stim2[key].tolist()
 
-        with open(filename, "w", encoding="utf-8") as f:
-            json.dump(stim2, f, ensure_ascii=False, indent=4)
-    else:
-        raise ValueError("stim should be a dict")
+    with filepath.open("w", encoding="utf-8") as file:
+        json.dump(stim2, file, ensure_ascii=False, indent=4)
 
 
 def to_mat(stim, filename):
@@ -204,19 +201,18 @@ def to_mat(stim, filename):
     Parameters
     ----------
     stim : dict
-        stimulus dictionary containing keys
-    filename : str
-        full path to the file to be creaated.
+        stimulus dictionary to export.
+    filename : Path or str
+        (full) path to the file to be created.
 
     """
-    if filename[-4:] != ".mat":
-        filename += ".mat"
-
-    if isinstance(stim, dict):
-        savemat(filename, stim)
-    else:
+    if not isinstance(stim, dict):
         raise ValueError("stim should be a dict")
 
+    filepath = Path(filename).resolve().with_suffix(".mat")
+
+    savemat(filepath, stim)
+
 
 def to_pickle(stim, filename):
     """Save stimulus-dict(s) as pickle-file
@@ -224,23 +220,21 @@ def to_pickle(stim, filename):
     Parameters
     ----------
     stim : dict
-        stimulus dictionary containing keys
-    filename : str
-        full path to the file to be creaated.
+        stimulus dictionary to export.
+    filename : Path or str
+        (full) path to the file to be created.
 
     """
-    if filename[-7:] != ".pickle":
-        filename += ".pickle"
+    if not isinstance(stim, dict):
+        raise ValueError("stim should be a dict")
 
-    if isinstance(stim, dict):
-        stim2 = copy.deepcopy(stim)
+    filepath = Path(filename).resolve().with_suffix(".pickle")
 
-        for key in stim2.keys():
-            # certain classes can cause problems for pickles; change to list
-            if key in ["visual_size", "ppd", "shape"]:
-                stim2[key] = list(stim2[key])
+    stim2 = copy.deepcopy(stim)
+    for key in stim2.keys():
+        # certain classes can cause problems for pickles; change to list
+        if key in ["visual_size", "ppd", "shape"]:
+            stim2[key] = list(stim2[key])
 
-        with open(filename, "wb") as handle:
-            pickle.dump(stim2, handle, protocol=pickle.HIGHEST_PROTOCOL)
-    else:
-        raise ValueError("stim should be a dict")
+    with filepath.open("wb") as file:
+        pickle.dump(stim2, file, protocol=pickle.HIGHEST_PROTOCOL)
diff --git a/stimupy/utils/pad.py b/stimupy/utils/pad.py
index d462b73e..4d73418a 100644
--- a/stimupy/utils/pad.py
+++ b/stimupy/utils/pad.py
@@ -1,4 +1,5 @@
 import copy
+
 import numpy as np
 
 from .utils import resolution
@@ -268,7 +269,7 @@ def pad_dict_by_visual_size(dct, padding, ppd, pad_value=0.0, keys=("img", "*mas
                     img = pad_by_visual_size(img, padding, ppd, pad_value)
                 new_dict[key] = img
 
-    # Update visual_size and shape-keys
+    # Update resolution
     new_dict["visual_size"] = resolution.visual_size_from_shape_ppd(img.shape, ppd)
     new_dict["shape"] = resolution.validate_shape(img.shape)
     return new_dict
@@ -301,7 +302,13 @@ def pad_dict_to_visual_size(dct, visual_size, ppd, pad_value=0, keys=("img", "*m
     shape = resolution.shape_from_visual_size_ppd(visual_size=visual_size, ppd=ppd)
 
     # pad
-    return pad_dict_to_shape(dct=dct, shape=shape, pad_value=pad_value, keys=keys)
+    padded_dict = pad_dict_to_shape(dct=dct, shape=shape, pad_value=pad_value, keys=keys)
+
+    # update resolution
+    padded_dict["ppd"] = ppd
+    padded_dict["visual_size"] = visual_size
+
+    return padded_dict
 
 
 def pad_dict_by_shape(dct, padding, pad_value=0, keys=("img", "*mask")):
@@ -430,7 +437,7 @@ def pad_dict_to_shape(dct, shape, pad_value=0, keys=("img", "*mask")):
                     img = pad_by_shape(img, padding=padding, pad_value=pad_value)
                 new_dict[key] = img
 
-    # Update visual_size and shape-keys
+    # Update resolution
     new_dict["shape"] = resolution.validate_shape(shape)
     if "ppd" in dct.keys():
         new_dict["visual_size"] = resolution.visual_size_from_shape_ppd(shape, dct["ppd"])
diff --git a/stimupy/utils/plotting.py b/stimupy/utils/plotting.py
index 52feadac..a5e7df5a 100644
--- a/stimupy/utils/plotting.py
+++ b/stimupy/utils/plotting.py
@@ -9,27 +9,10 @@
 __all__ = [
     "plot_stim",
     "plot_stimuli",
-    "compare_plots",
+    "plot_comparison",
 ]
 
 
-def compare_plots(plots):
-    """Plot multiple plots in one plot for comparing.
-
-    Parameters
-    ----------
-    plots : list of plots
-        List containing plots which should be plotted
-
-    """
-    M = len(plots)
-    for i, (plot_name, plot) in enumerate(plots.items()):
-        plt.subplot(1, M, i + 1)
-        plt.title(plot_name)
-        plt.imshow(plot, cmap="gray")
-    plt.show()
-
-
 def plot_stim(
     stim,
     mask=False,
@@ -228,3 +211,43 @@ def plot_stimuli(
         plt.close()
     else:
         raise ValueError("save can be None or a filepath")
+
+
+def plot_comparison(original_img, new_img):
+    """Plots visual comparison of two image-arrays
+
+    Parameters
+    ----------
+    original_img : numpy.ndarray
+        original image-array
+    new_img : numpy.ndarray
+        new image-array
+
+    Returns
+    -------
+    matplotlib.Figure
+        Figure containing plots of images, and their comparison(s)
+    """
+    vmin, vmax = 0, 1
+
+    fig = plt.figure(figsize=(20, 6))
+    plt.subplot(1, 4, 1)
+    plt.imshow(original_img, cmap="gray", vmin=vmin, vmax=vmax)
+    plt.title("Original")
+
+    plt.subplot(1, 4, 2)
+    plt.imshow(new_img, cmap="gray", vmin=vmin, vmax=vmax)
+    plt.title("New")
+
+    plt.subplot(1, 4, 3)
+    plt.imshow(original_img - new_img, cmap="coolwarm", vmin=-vmax, vmax=vmax)
+    plt.colorbar()
+    num_pix_off = np.prod(original_img.shape) - np.sum(np.isclose(original_img, new_img))
+    plt.title(f"Difference: {num_pix_off} pix")
+
+    plt.subplot(1, 4, 4)
+    plt.plot(original_img[:, 128], label="Original")
+    plt.plot(new_img[:, 128], label="New")
+    plt.legend()
+
+    return fig
diff --git a/stimupy/utils/utils.py b/stimupy/utils/utils.py
index 0adb8423..77254bb8 100644
--- a/stimupy/utils/utils.py
+++ b/stimupy/utils/utils.py
@@ -1,4 +1,6 @@
+import collections
 import copy
+import functools
 import itertools
 
 import numpy as np
@@ -18,6 +20,7 @@
     "flip_dict",
     "roll_dict",
     "strip_dict",
+    "make_two_sided",
     "permutate_params",
     "create_stimspace_stimuli",
 ]
@@ -440,6 +443,103 @@ def strip_dict(
     return new_dict
 
 
+def make_two_sided(func, two_sided_params):
+    """Create two-sided version of a stimulus function
+
+    Where (some) parameters can be specified separately for each side.
+    These parameters should then be specified as a 2-Sequence (2-ple, list of len=2),
+    where entry [0] is the parameter value for left side, and [1] for right side.
+    This means that if the kwarg takes a Sequence itself, e.g., `intensities`,
+    then the two-sided specification must be, e.g.,
+    ((int_left_0, int_left_1), (int_right_0, int_right_1))
+
+    Will be left- and right-sided.
+
+
+    Parameters
+    ----------
+    func : function
+        stimulus function to double
+    two_sided_params : Sequence[str]
+        names of parameters (kwargs) of func
+        that can be specified separately for each side of the display.
+
+    Returns
+    -------
+    function
+        two-sided version of stimulus function
+    """
+
+    @functools.wraps(func)
+    def two_sided_func(**kwargs):
+        shape = kwargs.pop("shape", None)
+        visual_size = kwargs.pop("visual_size", None)
+        ppd = kwargs.pop("ppd", None)
+
+        # Resolve resolution
+        try:
+            shape, visual_size, ppd = resolution.resolve(
+                shape=shape, visual_size=visual_size, ppd=ppd
+            )
+        except resolution.TooManyUnknownsError:
+            shape = resolution.validate_shape(shape)
+            visual_size = resolution.validate_visual_size(visual_size)
+            ppd = resolution.validate_ppd(ppd)
+        if visual_size.width is not None:
+            _visual_size = (visual_size.height, visual_size.width / 2)
+        else:
+            _visual_size = None
+        if shape.width is not None:
+            _shape = (shape.height, shape.width / 2)
+        else:
+            _shape = None
+
+        # Left side
+        left_side_params = {}
+        for key, value in kwargs.items():
+            if key in two_sided_params and isinstance(value, collections.abc.Sequence):
+                left_side_params[key] = value[0]
+            else:
+                left_side_params[key] = value
+
+        left = func(
+            visual_size=_visual_size,
+            ppd=ppd,
+            shape=_shape,
+            **left_side_params,
+        )
+
+        # Right side
+        right_side_params = {}
+        for key, value in kwargs.items():
+            if key in two_sided_params and isinstance(value, collections.abc.Sequence):
+                right_side_params[key] = value[1]
+            else:
+                right_side_params[key] = value
+
+        right = func(
+            visual_size=_visual_size,
+            ppd=ppd,
+            shape=_shape,
+            **right_side_params,
+        )
+
+        stim = stack_dicts(left, right, direction="horizontal")
+        stim["shape"] = left["shape"].height, left["shape"].width + right["shape"].width
+        stim["visual_size"] = (
+            left["visual_size"].height,
+            left["visual_size"].width + right["visual_size"].width,
+        )
+        for key in two_sided_params:
+            if key in kwargs:
+                stim[key] = kwargs[key]
+            else:
+                stim[key] = (left[key], right[key])
+        return stim
+
+    return two_sided_func
+
+
 def permutate_params(params):
     if not isinstance(params, dict):
         raise ValueError("params needs to be a dict with all stimulus parameters")


Test patch:
diff --git a/tests/papers/RHS2007.json b/tests/papers/RHS2007.json
index eb4f8f51..3ec78f70 100644
--- a/tests/papers/RHS2007.json
+++ b/tests/papers/RHS2007.json
@@ -24,15 +24,15 @@
         "mask": "9ecd24302fd639a3b1d88bc8d2473384"
     },
     "WE_radial_thick_small": {
-        "img": "3358a7ee8baea343c183c672c91b80c0",
+        "img": "90e0cd70342442fd7ce956022f002f02",
         "mask": "9dfa712381dad3c534d9d122bceac93a"
     },
     "WE_radial_thick": {
-        "img": "42db484da7509847c1073e8581948165",
+        "img": "5ef9cda7c331a97d0dbc9071af9e2fae",
         "mask": "773b67ccbfbc66f3e980745554d090d9"
     },
     "WE_radial_thin_small": {
-        "img": "4e971a779c6f63b904891018c4418b64",
+        "img": "d963190f8bf290cfe2f53ede6b955ce8",
         "mask": "5194ec93724d6ce281c08a823e3faaa6"
     },
     "WE_radial_thin": {
@@ -100,12 +100,12 @@
         "mask": "6920625ccfcc07d7e904b4b49188d39f"
     },
     "todorovic_benary1_2": {
-        "img": "973b958be7468c55e2f86094d21690e8",
-        "mask": "c2b7f645d607b64f3d74b5d64156af4b"
+        "img": "4fcde7a00d26712efa8d9db9d117fe18",
+        "mask": "3796426b61d4a19ec3e08d36d1b39bf0"
     },
     "todorovic_benary3_4": {
-        "img": "02fb79c76a8ba63ef6d85a88037ecf8f",
-        "mask": "e93a68afad23c70d25d553c0b8884d63"
+        "img": "90cbebcb9dfac6397acda3f60bfe7bab",
+        "mask": "a3e65dc662927408caec80fa9ab9e5f4"
     },
     "todorovic_benary1_2_3_4": {
         "img": "89b7ff63872f65e426bf9632befcb26c",
diff --git a/tests/papers/gen_ground_truth.py b/tests/papers/gen_ground_truth.py
index 6d14e884..93f121aa 100644
--- a/tests/papers/gen_ground_truth.py
+++ b/tests/papers/gen_ground_truth.py
@@ -22,13 +22,18 @@
         mask_keys = [key for key in stim.keys() if key.endswith("mask")]
         if "target_mask" in mask_keys:
             mask = stim["target_mask"]
-        else:
+        elif mask_keys:
             mask = stim[mask_keys[0]]
+        else:
+            mask = None
 
         stim.clear()
         stim["img"] = img
-        stim["mask"] = mask
-        export.arrays_to_checksum(stim, ["img", "mask"])
+        keys = ["img"]
+        if mask is not None:
+            stim["mask"] = mask
+            keys += ["mask"]
+        export.arrays_to_checksum(stim, keys)
 
     # Save all as .JSON
     export.to_json(stims, f"{d}/{paper}.json")
diff --git a/tests/papers/modelfest.json b/tests/papers/modelfest.json
index d064029c..f9e0f85a 100644
--- a/tests/papers/modelfest.json
+++ b/tests/papers/modelfest.json
@@ -156,8 +156,8 @@
         "mask": "9fba384c3a786d6a314fa5c94033975f"
     },
     "Disk40": {
-        "img": "fc5672879d59b22076ff1a79887df0c5",
-        "mask": "d774d1e22f9c1e5606e86cfad05391ed"
+        "img": "660f810e3761c81fe9147c0c5ab4e70d",
+        "mask": "b1aee35abaab5aef69c4a79f66b24fa7"
     },
     "Bessel41": {
         "img": "ad3bfbbde8c3954e2975a0d9c1179be0",
diff --git a/tests/papers/modelfest_compare.py b/tests/papers/modelfest_compare.py
new file mode 100644
index 00000000..fc28bcb3
--- /dev/null
+++ b/tests/papers/modelfest_compare.py
@@ -0,0 +1,76 @@
+from pathlib import Path
+
+import matplotlib.pyplot as plt
+import numpy as np
+from PIL import Image
+
+from stimupy.papers import modelfest
+from stimupy.papers.modelfest import *  # noqa: F401,F403
+from stimupy.utils import plot_comparison
+
+
+def compare_by_name(stim_name, image_dir):
+    """Compare current stimulus implementation against original image from file
+
+    Parameters
+    ----------
+    stim_name : str
+        name of the stimulus, which is both the `modelfest.<func_name>` and `<filename>.tif`
+    image_dir : Path or str
+        where to find the `.tif`
+
+    Returns
+    -------
+    matplotlib.Figure
+        comparison plot(s)
+
+    See also
+    --------
+    stimupy.utils.plotting.plot_comparison
+    """
+
+    original_img_file = Path(image_dir) / f"{stim_name}.tif"
+    original_img = np.array(Image.open(original_img_file))  # .astype(int)
+
+    func = globals()[stim_name]
+    new_img = func()["img"]
+    new_img = (new_img * 255.0).round()
+
+    original_img = original_img.astype(float) / 255.0
+    new_img /= 255.0
+
+    return plot_comparison(original_img, new_img)
+
+
+def compare_all(image_dir):
+    """Create comparisons for all ModelFest stimuli
+
+    For each `modelfest.<stim>()`-function there has to be a corresponding `<stim>.tif`
+
+    Parameters
+    ----------
+    image_dir : Path or str
+        where to find the `.tif`s
+
+    Returns
+    -------
+    None
+        instead, saves all comparison plots in `<image_dir>/comparisons`
+    """
+
+    comparisons_dir = Path(image_dir) / "comparisons"
+
+    if not comparisons_dir.exists():
+        comparisons_dir.mkdir()
+
+    for stim_name in modelfest.__all__:
+        f = compare_by_name(stim_name, image_dir)
+
+        plt.savefig(comparisons_dir / f"{stim_name}.png")
+        plt.close(f)
+
+
+if __name__ == "__main__":
+    image_dir = Path(__file__).parent.resolve() / "modelfest"
+
+    compare_all(image_dir=image_dir)