# import necessary packages
%matplotlib inline
%load_ext autoreload
%autoreload 2
import sys
sys.path.append('..')
import sfp 
import pandas as pd
import torch
import numpy as np
import matplotlib.pyplot as plt
from tqdm.auto import tqdm
import seaborn as sns
import pyrtools as pt

"""
Explanation: This paper presents a novel set of stimuli, which we refer to as "logpolar gratings" or "scaled gratings". These stimuli are sinusoidal gratings whose spatial frequency decreases as the reciprocal of the eccentricity. To be specific local spatial frequency is 
$$\omega_l(r,\theta)=\frac{\sqrt{\omega_r^2+\omega_a^2}}{r}$$
where coordinates $(r,\theta)$ specify the eccentricity and polar angle, relative to the center of the image. The angular frequency $\omega_a$ is an integer specifying the number of grating cycles per revolution around the image, while the radial frequency $\omega_r$ specifies the number of radians per unit increase in $\ln(r)$.
See the paper for more details about the stimuli and their specific use in this experiment. This notebook shows how you can use the software included in this repo to create these stimuli and presents linear approximations of them (but see spatial-frequency-stimuli for a stand-alone repo with a command-line interface for creating sets of these stimuli).
End of explanation
"""


stim = sfp.stimuli.log_polar_grating(200, w_a=3)
pt.imshow(stim);

"""
Explanation: Creation of stimuli
Our stimuli can be created by calling sfp.stimuli.log_polar_grating, and specifying the size and several other parameters (note that we only support square gratings):
End of explanation
"""


stim = []
for i in [3, 6]:
    stim.append(sfp.stimuli.log_polar_grating(200, w_r=i))
pt.imshow(stim);

"""
Explanation: In the visualization above, with w_r=0,w_a=3, we can see clearly that w_a specifies the number of grating cycles per revolution across the image.
End of explanation
"""


w_r = 6
w_a = 11
phase = 0
stim = sfp.stimuli.log_polar_grating(500, w_r=w_r, w_a=w_a, phi=phase)
dx, dy, _, _ = sfp.stimuli.create_sf_maps_cpp(stim.shape[0], w_r=w_r, w_a=w_a)

pt.imshow([stim], zoom=.5);

"""
Explanation: w_r is a little bit harder to grok, but in the above we can see that increasing it increases the number of times the grating cycles in the radial direction.
In our experiment, we created an anti-aliasing mask to put at the center of each stimulus image, as well as rescaled the stimuli to account for the projector's modulation transfer function (MTF).
Linear approximations of stimuli
The following shows some linear approximations of our stimuli. First, we generate the stimuli, and then we calculate the analytic spatial frequency with respect to x and y everywhere in the image.
End of explanation
"""


masked_grating, masked_approx = sfp.plotting.plot_grating_approximation(stim, dx, dy, num_windows=20, w_r=w_r, w_a=w_a, phase=phase)

"""
Explanation: We then use these spatial frequency maps to create a linear approximation: the top shows a windowed view of our actual stimuli, while the bottom shows a windowed view of the linear approximation. The approximation is created by creating a standard sinusoid (that is, one with a constant spatial frequency and orientation) with the analytic dx and dy calculated above at the center of each of the little windows. For relatively small values of w_r,w_a and relatively large numbers of windows, the approximation willb e pretty good, though you can see it start to fail near the center of the image, where the actual stimulus's spatial frequency is changing rapidly, and so the linear approximation (which has a constant spatial frequency in each image) is notably different.
End of explanation
"""


masked_grating, masked_approx = sfp.plotting.plot_grating_approximation(stim, dx, dy, num_windows=10, w_r=w_r, w_a=w_a, phase=phase)

"""
Explanation: If we decrease the number of windows, we can see this more clearly:
End of explanation
"""


