Help with Gradient Descent Implementation MATLAB

Question

I am trying to generate a window for a time single with a specific (asymmetric shape); my ideal window would have a frequency response as shown below.

Step function frequency response

Ideally whenever a target is present the amplitude of the signal will increase by an amount proportional to the target magnitude. What I am struggling to do now is realize this filter using gradient descent. I want to initialize a window with random values and then use gradient descent to achieve the optimum filter in a least squares sense

The loss function I am using is

$ L(f) = \sum(|I(f)|-|X(f)|)^2 $ where $ I(f)$ is the ideal frequency response and $X(f)$ is the frequency response of the current iteration. Then I can calculate the update to my current window values $x[n]$ as

$\frac{\partial L(f)}{\partial x[n]} =\frac{\partial L(f)}{\partial |X(f)|}\frac{\partial |X(f)|}{\partial X(f)}\frac{\partial X(f)}{\partial x[n]} $

where

$ \frac{\partial L(f)}{\partial |X(f)|} = -2(|I(f)|-|X(f)|)$

$\frac{\partial |X(f)|}{\partial X(f)}= 1/2(X(f)X^*(f))^{-\frac{1}{2}}(X^*(f))$

$\frac{\partial X(f)}{\partial x[n]} = e^{-i2\pi fTn}$

where $f$ is the frequency I of interest and $T$ is the sampling period. Once I calculate $\frac{\partial L(f)}{\partial x[n]}$ for all frequencies and n then I sum across frequencies and update $x[n]$ based on the stepsize and the calculated derivative. For some reason though I do not seem to be converging to anything that remotely resmbles my ideal function. Below is the matlab script I am using to calculate my gradient descent.

My question is what is going on, why am I not converging to something that at least appears similar to the expected function? I have tried a variety of step sizes/iterations and either the loss function grows or the resulting filter does not look as expected. Should there be some kind of normalization which I am missing or is there something fundamentally wrong with how I am approaching the problem.

Loss vs Iteration Fourier Transform of Optimized Window

clear
clc

% Radar parameters setup
ts = 0; %Time start
te = 25.8e-6; %Time end
fADC = 5000e3; %Samples per second
TADC = 1/fADC; %Sampling period
t = ts:TADC:te; %time 
fs = 60e9; % starting frequency
slope = 144.984e12; %Radar slope

% Physical constants setup
c = 3e8;

%% Define Ideal Window Frequency Response
% Assume that the data is collected 
N = 256; % The number of frequencies to calculate the discrete fourier transform at 
f = linspace(-fADC/2,fADC/2,N); % The frequencies
I = [zeros(1,N/2),ones(1,N/2)];
figure(2)
plot(f,I,'*r')
xlabel('frequency (Hz)')
ylabel('FFT magnitude')
set(gca,'fontsize',18)


%% Initialize Window
xn = ones(size(t)); %the window intially

figure(3)
plot(t,real(xn),'b',t,imag(xn),'r','linewidth',2)
set(gca,'fontsize',18);grid on

xnFFT = fftshift(fft(xn,N));
figure(4)
plot(abs(xnFFT))

%% Optimize Window
% Set descent parameters
stepSize = 1e-6; % The step size to take
iter = 500000;     % The number of iterations
L = zeros(iter,1); % Holds the loss for each iteration

% Define the fourier transform matrix to multiply the signal by %xn is
% 1x130. So 
DFTmat = exp(-1i*2*pi*repmat(f',1,length(t)).*repmat(t,N,1)); %t increase with column, %f increases by row
dXf_xn = DFTmat;

for i = 1:iter
    % Step 1 calculate the DFT
    Xf = (DFTmat*xn.').';
%     figure(5)
%     plot(abs(Xf))
    
    % Step 2 calculate the loss
    L(i) = sum((abs(I)-abs(Xf)).^2);

    % Step 3 calculate partial derivatives
    dL_dmXf = -2*(abs(I)-abs(Xf)); % Partial derivative of loss with magnitude |Xf|
    dmXf_dXf = 1/2*(Xf.*conj(Xf)).^(-1/2).*conj(Xf); %d|Xf|/dXf
    dL_dXf = dL_dmXf.*dmXf_dXf; %dL/dXf
    
    % Step 4 calculate the update
    xnUp = sum(repmat(transpose(dL_dXf),1,length(t)).*dXf_xn,1)/length(t);
    
    % Step 5 apply the update
    xn = xn-stepSize*xnUp;
end

%%
figure(6)
plot(L)
xlabel('Iteration')
ylabel('Loss')

xnFFT = fftshift(fft(xn));
figure(7)
plot(abs(xnFFT))

figure(8)
plot(t,real(xn),'r',t,imag(xn),'b')```