## 2 Materials and Methods
### The Hierarchical Bayesian MAP
Next, we solve \(\widetilde{\mathbf{h}}\) in Fourier domain,  
\[p\Big{(}\overline{\mathbf{P}}\big{|}\widetilde{\mathbf{f}},\Sigma_{n}, \widetilde{\mathbf{h}}\Big{)}\ =\ \frac{1}{(2\pi)^{\frac{n}{2}}|\Sigma_{n}|^{\frac{1}{2}}}\exp \biggl{(}-\frac{1}{2}\Big{(}\overline{\mathbf{P}}-\Sigma_{Af}\widetilde{ \mathbf{h}}\Big{)}^{T}\Sigma_{n}^{-1}\Big{(}\overline{\mathbf{P}}-\Sigma_{Af} \widetilde{\mathbf{h}}\Big{)}\biggr{)} \tag{16}\]  
where \(\mathbf{y}\ =\ \mathbf{A}\mathbf{f}\in\mathbb{R}^{m\times 1}\), \(\overline{\mathbf{y}}\) is the Fourier transform of \(\mathbf{y}\) with \(\Sigma_{Af}\)\(=\) diag\((\overline{\mathbf{y}})\). By calculating by (15) and (16), the conditional distribution \(\widetilde{\mathbf{h}}\) is followed as  
\[p\Big{(}\widetilde{\mathbf{h}}\big{|}\widetilde{\mathbf{f}},\overline{ \mathbf{P}},\Sigma_{n}\Big{)}\propto\exp\biggl{(}-\frac{1}{2}\Big{(} \overline{\mathbf{P}}-\Sigma_{Af}\widetilde{\mathbf{h}}\Big{)}^{T}\Sigma_{n} ^{-1}\Big{(}\overline{\mathbf{P}}-\Sigma_{Af}\widetilde{\mathbf{h}}\Big{)} \biggr{)}\exp\biggl{(}-\frac{1}{2}\Big{(}\widetilde{\mathbf{h}}-\widetilde{ \mathbf{h}}_{0}\Big{)}^{T}\Sigma_{h}^{-1}\Big{(}\widetilde{\mathbf{h}}- \widetilde{\mathbf{h}}_{0}\Big{)}\biggr{)} \tag{17}\]  
According to (17), the conditional distribution \(\widetilde{\mathbf{h}}\) can be described by a multivariate Gaussian distribution. For simplicity, we assume that \(\sum_{n}\ =\ \sigma_{n}^{2}\mathbf{I}\). Since \(\mathbf{H}\) is circular, it holds \(\Sigma_{h}\ =\ \sigma_{h}^{2}\mathbf{I}\). The mean vector and covariance matrix of the multivariate Gaussian distribution is given respectively,  
\[\widetilde{\mathbf{m}}\ =\ \Bigg{(}\frac{\mathbf{I}}{\sigma_{h}^{2}}+\frac{ \left|\Sigma_{Af}\right|^{2}}{\sigma_{n}^{2}}\Bigg{)}^{-1}\Bigg{(}\frac{ \widetilde{\mathbf{h}}_{0}}{\sigma_{h}^{2}}+\frac{\Sigma_{Af}^{T}\overline{ \mathbf{P}}}{\sigma_{n}^{2}}\Bigg{)} \tag{18}\]  
\[\mathbf{M}(\mathbf{x})\ =\ \mathbf{A}\mathbf{H}\psi\mathbf{G}(\mathbf{x})\sum \limits_{v}\mathbf{G}(\mathbf{x})\psi^{\mathrm{T}}\mathbf{H}^{\mathrm{T}} \mathbf{A}^{\mathrm{T}}+\sum\
olimits_{m^{\prime}}\quad\mathbf{M}(\mathbf{x}) \in\mathbb{R}^{m\times m} \tag{10}\]0  
The multivariate Gaussian distribution is easy to sample, and the conditional distribution \(\widetilde{\mathbf{h}}\) can be obtained. The inverse Fourier transform of the conditional distribution of \(\widetilde{\mathbf{h}}\) was implemented to obtain the point spread function \(\mathbf{h}\). Furthermore, the cyclic shift \(\mathbf{h}\) was implemented to obtain the degradation matrix \(\mathbf{H}\). The proposed algorithm is described in the \[\mathbf{M}(\mathbf{x})\ =\ \mathbf{A}\mathbf{H}\psi\mathbf{G}(\mathbf{x})\sum \limits_{v}\mathbf{G}(\mathbf{x})\psi^{\mathrm{T}}\mathbf{H}^{\mathrm{T}} \mathbf{A}^{\mathrm{T}}+\sum\
olimits_{m^{\prime}}\quad\mathbf{M}(\mathbf{x}) \in\mathbb{R}^{m\times m} \tag{10}\]1  
To reconstruct CT image with far fewer projection data, a compressed sensing image reconstruction algorithm including a blind image restoration process was proposed. Image reconstruction is achieved by estimating a sparse coefficient vector. We utilized the hierarchical Bayesian model and transformed the image estimation into two independent variables estimation. The hierarchical Bayesian MAP method solved the problem of parameter estimation. The introduction of blind image restoration mechanism improves the quality of the output image in each iteration and improves the quality of the final reconstructed image. Therefore, under the same image quality requirements, the proposed method can effectively reduce the radiation dose.