Hi Lindsey,
Sorry for the delay in responding. It took some time to convince myself I understood what was happening, and I was working on putting together some examples to demonstrate things. But here’s the short version:
The issue is not with the mask or how the alpha channel is rendered. The problem is the base grating.
Your genSinusoid function generates a normal sinusoid ranging from -1 to 1. This is fine if you’re working in linear color (which is what MWorks does internally). However, you’re putting these values (after converting them to integers in the range [0,255]) in to a PNG file with no colorspace information. The standard way to interpret a PNG (or any other image) with no colorspace info is to assume that the colors are in the sRGB colorspce, which uses a nonlinear gamma.
This means that when color management is enabled in MWorks (which it is by default), MWorks applies a sRGB-to-linear conversion of the PNG’s colors before rendering it. Since the image colors aren’t actually in sRGB, they end up looking quite different than you’d expect.
The simplest way to solve this is to convert the computed grating values from linear gamma to sRGB gamma before storing them in the PNG. Here’s a modified version of your image generation code that does this:
for omega_x, cpd_x in zip(omega_list, cpd_str):
omega = [omega_x, 0.0]
rho = 0
sinusoidParam = {'A':1, 'omega':omega, 'rho':rho, 'sz':(stim_diameter, stim_diameter)}
sin = genSinusoid(**sinusoidParam)
# Normalize range to [0, 1]
sin = (sin + 1.0) / 2.0
# Convert from linear color to sRGB
sin = np.piecewise(sin,
[sin < 0.0031308],
[(lambda x: 12.92 * x),
(lambda x: 1.055 * np.power(x, 1.0/2.4) - 0.055)])
# Convert from float64 [0, 1] to uint8 [0, 255]
sin = np.round_(255.0 * sin).astype(np.uint8)
no_mask_file_name = f"/grat_SF_{cpd_x}_no_mask.png"
im_sin = Image.fromarray(sin, mode='L').convert('RGB')
#
# Include the PNG sRGB chunk in the output image. While not strictly
# necessary, this makes it clear that we really mean for the image colors
# to be in the sRGB color space.
# https://www.w3.org/TR/png/#srgb-standard-colour-space
#
pnginfo = PngInfo()
pnginfo.add(b'sRGB', b'\x00')
im_sin.save(save_path + no_mask_file_name, pnginfo=pnginfo)
rgba = np.array(im_sin.convert('RGBA'))
mask_0_base = (mask_flat_peak - max(mask_flat_peak[0,:])) * 255 # convert mask from 0-1 to 0-255
mask_0_base[mask_0_base<0] = 0 # force edge (outside of incircle) = 0
mask_0_base = mask_0_base.astype(np.uint8)
rgba[:, :, -1] = mask_0_base # assign mask to the last channel (alpha)
img_rgba = Image.fromarray(rgba)
# Again, include the PNG sRGB chunk in the output image
img_rgba.save(save_path + f"/grat_alpha_{cpd_x}.png", pnginfo=pnginfo)
FYI, I found that the matplotlib function you were using to generate the mask-free image alters the colors in ways I don’t understand. I recommend using PIL to create both images, which is what my modified code does.
If you think more detailed examples would be helpful, I can finish putting those together and send them along. But I’m pretty sure that switching to the above code will eliminate the discrepancy between your generated gratings and MWorks’ gratings.
Cheers,
Chris