I see you have commented out the definition of "deltat" and have an incorrect computation for "sigmap" that follows.
Assuming you define "deltat" the corrected portion of your code will be :
% Assuming deltat is defined, I have taken an arbitrary value for mock computation:
deltat = 50e-12; % pulse width in seconds
% Correct calculation of sigmaP for a Gaussian pulse in the frequency domain
sigmaP = 1 / (deltat * sqrt(2 * log(2)));
% Corrected alpha calculation using the Gaussian envelope formula
alpha = exp(-((wP - (wS + wI))^2) / (2 * sigmaP^2));
% Continue with phi and f as before
phi = (sinc(Deltak * L / 2)) .* exp(1i * Deltak * L / 2); % phase matching function
f = alpha * phi; % joint spectral amplitude
joint_spectral_intensity = abs((f).^2); % joint spectral intensity
Hope this Helps!
Soumnath
