First things first. A couple of notes about this problem.
- The "new" transformation is nothing but the chirp $z$-transform introduced by Rabiner et. al. in "The chirp z-transform algorithm".
- The efficient algorithm to compute CZT is based on the Bluestein's algorithm.
Below I attempt to explain it.
Now let's go back to the problem. Let's consider the $z$-transform of a sequence $x[n]$
$$
X(z_{n}) = \sum_{k=0}^{N-1}x[k]z_{n}^{-k},~n = 0,\ldots, K-1
$$
where $z_{n} \triangleq A\hat{W}^{n} = e^{i\theta_{0}}e^{-in\Delta\theta} $. Substituting $z_{n}$ yields the transformation
\begin{equation}
X^{c}[n] = \sum_{k=0}^{N-1}x[k]A^{-k}\hat{W}^{-kn}
\tag{1}
\label{eq:1}
\end{equation}
Note that if $A = 1, \hat{W} = e^{2\pi i\frac{1}{K}}$, then $X^{c}[n]$ is DFT.
Let's use the relationship
$$
kn = \frac{k^{2}+n^{2}-(k-n)^{2}}{2}
$$
in \eqref{eq:1} which yields
$$
\begin{eqnarray}
X^{c}[n] & = & \sum_{k=0}^{N-1}x[k]A^{-k}\hat{W}^{-\frac{n^{2}+k^{2}-(n-k)^{2}}{2}} \\
& = & \hat{W}^{-\frac{n^{2}}{2}}\sum_{k=0}^{N-1}\underbrace{x[k]A^{-k}\hat{W}^{-\frac{k^{2}}{2}}}_{g[k]}\underbrace{\hat{W}^{\frac{(n-k)^{2}}{2}}}_{h[n-k]} \\
& = & \hat{W}^{-\frac{n^{2}}{2}}\underbrace{\sum_{k=0}^{N-1}g[k]h[n-k]}_{(g * h)[n]\triangleq v[n]},~n = 0,\ldots, K-1 \tag{2}\label{eq:2}
\end{eqnarray}
$$
The block diagram below shows how $X^{c}[n]$ is obtained.
This is what you have already shown. To efficiently compute this transform using radix-2 FFT, let's focus on the sum in \eqref{eq:2}.
The signals involved in this sum (convolution) are
- $g[n] = x[n]A^{-n}\hat{W}^{-\frac{n^{2}}{2}},~n = 0,\ldots,N-1$
- $h[n] = \hat{W}^{\frac{n^{2}}{2}},~n = -N+1,-N+2,\ldots,0,1,K-1$
Recall that linear convolution of two finite sequences of length $P$ and $Q$ can be done using DFT of length $R > P+Q-1$ which is needed to avoid time-aliasing when periodizing the sequences to perform circular convolution.
In the problem at hand, the filter $h[n]$, which is a chirp sequence yielding the name of chirp transform, is not a finite sequence. But for a certain pair $(N,K)$, $h[n]$ can be truncated as below
$$
h[n] = \left\{\begin{array}{ll}\hat{W}^{\frac{n^{2}}{2}} & n = -N+1,\ldots, K-1 \\ 0 & \text{otherwise}\end{array}\right .
$$
The figure below pictorially shows $g[n]$ and $h[n]$ involved in the convolution for the pair $(N,K) = (4,3)$
Let $R > N+K-1 = 6$ be the next integer of power of 2, thus $R = 8$. Hence, to perform the linear convolution in \eqref{eq:2} in an efficient way, zero-pad $g[n]$ up to length $R = 8$ and periodize it with period $R$
$$
\tilde{g}[n] = \sum_{r=-\infty}^{+\infty}g[n-rR]
$$
which is pictorially shown below
Note that if the Fourier transform $G(e^{i\omega})$ of a sequence $g[n]$ with length $R$ is sampled at frequencies $\omega_{m} = 2\pi m/R$, then the resulting sequence corresponds to the discrete Fourier series coefficients of the periodic sequence $\tilde{g}[n]$. From the definition of the discrete Fourier transform (DFT), it follows that the finite-length sequence
$$
G[m] = \left\{\begin{array}{ll}G(e^{i\omega_{m}}) & 0\leq m\leq R-1 \\
0 & \text{otherwise}\end{array}\right .
$$
is the DFT of one period of $\tilde{g}[n]$, $g_{p}[n]$, given as
$$
g_{p}[n] = \left \{\begin{array}{ll}\tilde{g}[n] & 0\leq n \leq R-1 \\
0 & \text{otherwise}\end{array}\right .
$$
Perform the same operation on $h[n]$, i.e., zero-pad $h[n]$ to a length $R = 8$ and periodize.
$$
\tilde{h}[n] = \sum_{r=-\infty}^{+\infty}h[n-rR]
$$
For the example $(N,K) = (4,3)$, the zero-padded and periodized $h[n]$ is shown below
Note that one period of $\tilde{h}[n]$ with the support $0\leq n \leq R-1$ is given by
\begin{eqnarray}
h_{p}[n] & = & \left\{\begin{array}{ll}\tilde{h}[n] & 0\leq n \leq R-1 \\
0 & \text{otherwise}\end{array} \right . \\
& = & \left\{\begin{array}{ll}\hat{W}^{\frac{n^{2}}{2}} & 0\leq n \leq K-1 \\
\hat{W}^{\frac{(n-R)^{2}}{2}} & R-N+1\leq n \leq R-1 \\
0 & \text{otherwise}\end{array}\right . \\
& = & \left \{\begin{array}{ll}h[n] & 0\leq n \leq K-1 \\
h[n-R] & R-N+1\leq n \leq R-1 \\
0 & \text{otherwise}\end{array}\right .
\tag{3}
\label{eq:3}
\end{eqnarray}
Similarly, the finite-length sequence
$$
H[m] = \left\{\begin{array}{ll}H(e^{i\omega_{m}}) & 0\leq m\leq R-1 \\
0 & \text{otherwise}\end{array}\right .
$$
is the DFT of $h_{p}[n]$.
From circular convolution theorem, $H[m]G[m]$ is the DFT of the one-period sequence $v_{p}[n]$
$$
v_{p}[n] = \left\{\begin{array}{ll}\tilde{v}[n] & 0\leq n\leq R-1 \\
0 & \text{otherwise}
\end{array}\right .
$$
where $\tilde{v}[n] = \sum_{k=0}^{R-1}\tilde{g}[k]\tilde{h}[n-k]$. Note that $v_{p}[n]$ is the circular convolution of $g_{p}$ and $h_{p}$, i.e.,
$$
v_{p}[n] = (g_{p}\circledast_{R} h_{p})[n] = \sum_{k=0}^{R-1}g_{p}[k]h_{p}[(n-k)_{R}],~0\leq n\leq R-1
$$
where $\circledast_{R}$ denotes $R$-point circular convolution and $(n-k)_{R}$ is $\text{mod}(n-k,R)$.
Hence,
$$
v[n] = \underline{\cal{F}}_{R}^{-1}\left\{\underline{\cal{F}}_{R}\{g_{p}\}\underline{\cal{F}}_{R}\{h_{p}\}\right\}[n]
$$
Here is a snippet MATLAB code
% Parameters
N = 4;
K = 3;
R = 2^nextpow2(N+K-1);
deltaTheta = 2*pi*rand;
theta0 = rand;
A = exp(1i*theta0);
What = exp(-1i*deltaTheta);
% Signal
x = randn(N,1);
% Setup g and h
NN = (0:N-1)';
NK = (-N+1:K-1)';
g = x.*( ...
(A.^(-NN)).* ...
(What.^(-(NN.^2/2))) ...
); % x*A^(-n)*What^(-n^2/2)
h = What.^((NK.^2/2));
% Zero-pad and periodize
gp = [g;zeros(R-N,1)]; % note that MATLAB fft does zero-pad. I did it here for illustration.
hp = circshift([h;zeros(R-(N+K-1),1)],-(N-1)); % Eq. (3)
% Compute convolution
vp = ifft(fft(gp,R).*fft(hp,R),R);
% Transform at n = 0,1,...,K-1
KK = (0:K-1)';
Xc = What.^(-KK.^2/2).*vp(1:K);
% In MATLAB, there is a CZT function
tXc = czt(x,K,1/What,A);
