% -*- mode: Noweb; noweb-code-mode: c-mode -*- \section{Geometric adaptive optics} \label{sec:open-loop-geometric} In this section, the performances of several AO systems are compared. All the systems employ geometric Shack--Hartmann wavefront sensors (SH--WFS). For each system, the on--axis wavefront is estimated several times with a different random draw of the phase screens. The components common to all the systems are defined first: \begin{itemize} \item the science source, <>= source src; src.setup("K",ARCSEC(0) , 0, INFINITY,(N_SIDE_LENSLET+1)*(N_SIDE_LENSLET+1), "SRC"); <>= src.cleanup(); @ \item the atmosphere, <>= atmosphere atm; //atm.setup(20e-2,30,10e3,10,0); atm.gmt_setup(15e-2,60); /* float altitude[] = {0, 10e3}, xi0[] = {0.5, 0.5}, wind_speed[] = {10, 10}, wind_direction[] = {0, 0}; atm.setup(20e-2,30,altitude, xi0, wind_speed, wind_direction); */ <>= atm.cleanup(); @ \item the diameter of the telescope, <>= float D = 25; // telescope diameter in meter @ leading to a lenslet size of: <>= float d = D/N_SIDE_LENSLET; @ \item the pupil mask. <>= mask pupil_mask; pupil_mask.setup( (N_SIDE_LENSLET+1)*(N_SIDE_LENSLET+1) ); <>= pupil_mask.cleanup(); @ \item the wavefront sensor centroid container, <>= centroiding cog; cog.setup(N_SIDE_LENSLET,1); <>= cog.cleanup(); @ \item the diffraction limited science imager, <>= imaging imager; imager.setup(N_SIDE_LENSLET+1,1,4,1,1); <>= imager.cleanup(); @ \item the turbulence limited science imager, <>= imaging imager_turb; imager_turb.setup(N_SIDE_LENSLET+1,1,4,1,1); <>= imager_turb.cleanup(); @ \item the science imager for NGSAO, <>= imaging imager_ngsao; imager_ngsao.setup(N_SIDE_LENSLET+1,1,4,1,1); <>= imager_ngsao.cleanup(); @ \item the science imager for LGSAO, <>= imaging imager_lgsao; imager_lgsao.setup(N_SIDE_LENSLET+1,1,4,1,1); <>= imager_lgsao.cleanup(); @ \item the science imager for LGSAO LMMSE, <>= imaging imager_lgsao_lmmse; imager_lgsao_lmmse.setup(N_SIDE_LENSLET+1,1,4,1,1); <>= imager_lgsao_lmmse.cleanup(); @ \item the science imager for LTAO, <>= imaging imager_ltao; imager_ltao.setup(N_SIDE_LENSLET+1,1,4,1,1); <>= imager_ltao.cleanup(); @ \item the statistical tool. <>= stats S; S.setup(); <>= S.cleanup(); @ \end{itemize} @ The Fried geometry for a circular pupil with the intensity [[threshold]] is enforced:: <>= float threshold = 0.5; cog.MASK_SET = fried_geometry_setup(cog.lenslet_mask, pupil_mask.m, N_SIDE_LENSLET, 16, threshold); @ The filtering properties associated with the pupil are set with: <>= pupil_mask.set_filter(); @ The science is propagated through the [[pupil_mask]] to the focal plane of the [[imager]]: <>= src.wavefront.M = &pupil_mask; src.wavefront.masked(); imager.propagate(&src); //imager.frame2file("refFrame.bin"); char plotly_name[64], plotly_dir[64]; sprintf(plotly_dir,"GEAOS/D=%.1fm, N=%d (%d samples)/",D,N_SIDE_LENSLET,N_SAMPLE); SET_PLOTLY_NAME(plotly_name,plotly_dir,"frames/diffraction limited"); imager.show_frame(plotly_name); @ A few useful variables are defined here: <>= int NP, NP2, k_SAMPLE; float tau=0.; NP = N_SIDE_LENSLET+1; NP2 = NP*NP; k_SAMPLE = 0; @ The science wavefront is propagated through the atmosphere from [[src]] to [[pupil_mask]]. <>= atm.get_phase_screen(&src,d,NP,d,NP,tau); src.wavefront.masked(); if (k_SAMPLE==(N_SAMPLE-1)) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"phases/science phase screen"); src.wavefront.show_phase(plotly_name); } <>= float science_wf_rms=0.; <>= science_wf_rms += S.var(src.wavefront.phase, &pupil_mask, NP2); @ and then propagated to the focal plane of the imager: <>= imager_turb.propagate(&src); if (k_SAMPLE==(N_SAMPLE-1)) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"frames/turbulence limited"); imager_turb.show_frame(plotly_name); } @ The turbulence wavefront is saved apart for later use: <>= complex_amplitude phase_screen; phase_screen.setup((N_SIDE_LENSLET+1)*(N_SIDE_LENSLET+1)); phase_screen.add_phase(1,src.wavefront.phase); @ All CEO programs must include the following headers which also contains the headers for all the CEO library modules. <
>= #ifndef __CEO_H__ #include "ceo.h" #endif @ The size of the lenslet array is defined in the header: <
>= #define N_SIDE_LENSLET 60 #define _N_LENSLET_ (N_SIDE_LENSLET*N_SIDE_LENSLET) @ The number of atmosphere sample <
>= #define N_SAMPLE 1000 #define PLOTLY_LIM (N_SAMPLE-1) // (N_SAMPLE-1) for plotting, higher for disabling plotting @ The main function is: <>= <
> void SET_PLOTLY_NAME(char *name,char *dir,char *path) { strcpy(name, dir); strcat(name,path); } int main(int argc,char *argv[]) { cudaSetDevice(1); <> fprintf(stderr," Samples: %d:\n",N_SAMPLE); for (k_SAMPLE=0;k_SAMPLE> <> <> <> } fprintf(stderr,"\n"); <> <> } @ \subsection{Natural guide star adaptive optics} \label{sec:ngs-ao} The natural guide star (NGS) is given by: <>= source ngs; ngs.setup("R",ARCSEC(0) , 0, INFINITY, "NGS"); <>= ngs.cleanup(); @ The NGS wavefronts gradients from the turbulence are generated from a dedicated slopes turbulence generator. <>= atm.get_phase_screen_gradient(&cog,N_SIDE_LENSLET,d,&ngs,tau); @ In the following, the wavefront is reconstructed from the NGS centroids: <>= LMMSE ngs_lmmse; ngs_lmmse.setup(&atm,&ngs,1,&src,1,d,N_SIDE_LENSLET,&pupil_mask,"MINRES"); int ngsao_cvgce_iteration=0; float ngsao_elapsed_time=0.; <>= ngs_lmmse.estimation(&cog); ngsao_elapsed_time += ngs_lmmse.elapsed_time; ngsao_cvgce_iteration += ngs_lmmse.iSolve.cvgce_iteration; @ The estimated wavefront is removed from the science wavefront <>= src.wavefront.add_phase(-1,ngs_lmmse.d__phase_est); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"phases/LMMSE NGS AO"); src.wavefront.show_phase(plotly_name); } <>= float ngsao_wfe_rms=0.; <>= ngsao_wfe_rms += S.var(src.wavefront.phase, &pupil_mask, NP2); <>= ngs_lmmse.cleanup(); @ and the residual wavefront corresponding image is computed. <>= imager_ngsao.propagate(&src); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"frames/LMMSE NGS AO"); imager_ngsao.show_frame(plotly_name,&imager); } @ \subsection{Laser guide star adaptive optics} \label{sec:lgs-ao} The laser guide star (LGS) is given by: <>= source lgs; lgs.setup("R",ARCSEC(0) , 0, 90e3, (N_SIDE_LENSLET+1)*(N_SIDE_LENSLET+1), "LGS"); <>= lgs.cleanup(); @ The LGS turbulence wavefront and wavefront gradient are computed with: <>= atm.get_phase_screen(&lgs,d,NP,d,NP,0); lgs.wavefront.M = &pupil_mask; lgs.wavefront.masked(); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"phases/LGS phase screen"); lgs.wavefront.show_phase(plotly_name); } <>= float lgs_wf_rms=0.; <>= lgs_wf_rms += S.var(lgs.wavefront.phase, &pupil_mask, NP2); @ and with: <>= atm.get_phase_screen_gradient(&cog,N_SIDE_LENSLET,d,&lgs,tau); @ The LGS wavefront is subtracted from the science wavefront: <>= src.wavefront.reset(phase_screen); src.wavefront.add_phase(-1,lgs.wavefront.phase); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"phases/LGS AO"); src.wavefront.show_phase(plotly_name); } <>= float lgs_wfe_rms=0.; <>= lgs_wfe_rms += S.var(src.wavefront.phase, &pupil_mask, NP2); @ and the residual wavefront corresponding image is computed. <>= imager_lgsao.propagate(&src); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"frames/LGS AO"); imager_lgsao.show_frame(plotly_name,&imager); } @ The wavefront is reconstructed from the LGS centroids: <>= LMMSE lgs_lmmse; lgs_lmmse.setup(&atm,&lgs,1,&src,1,d,N_SIDE_LENSLET,&pupil_mask,"MINRES"); int lgsao_cvgce_iteration=0; float lgsao_elapsed_time=0.; <>= lgs_lmmse.estimation(&cog); lgsao_elapsed_time += lgs_lmmse.elapsed_time; lgsao_cvgce_iteration += lgs_lmmse.iSolve.cvgce_iteration; @ The LMMMSE LGS wavefront estimate is subtracted from the science wavefront: <>= src.wavefront.reset(phase_screen); src.wavefront.add_phase(-1,lgs_lmmse.d__phase_est); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"phases/LMMSE LGS AO"); src.wavefront.show_phase(plotly_name); } <>= float lgs_wfe_rms_lmmse=0.; <>= lgs_wfe_rms_lmmse += S.var(src.wavefront.phase, &pupil_mask, NP2); //lgs_lmmse.toFile("phaseEstLgsItp.bin"); <>= lgs_lmmse.cleanup(); @ and the residual wavefront corresponding image is computed. <>= imager_lgsao_lmmse.propagate(&src); //imager.frame2file("lgsLmmseFrame.bin"); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"frames/LMMSE LGS AO"); imager_lgsao_lmmse.show_frame(plotly_name,&imager); } @ \subsection{Laser tomography adaptive optics} \label{sec:ltao} For the LTAO wavefront estimation, the LGS constellation is defined first. We will use 3 LGSs on a 30 arcsec radius ring. <>= int N_GS = 6; float gs_radius = 30; // 3 LGS on a ring /* float zenith[] = {ARCSEC(gs_radius),ARCSEC(gs_radius),ARCSEC(gs_radius)}, */ /* azimuth[] = {0,2.*PI/3.,4.*PI/3.}; */ // 6 LGS on a ring float zenith[] = {ARCSEC(30),ARCSEC(30),ARCSEC(30),ARCSEC(30),ARCSEC(30),ARCSEC(30)}, azimuth[] = {0,2.*PI/6.,4.*PI/6.,6.*PI/6.,8.*PI/6.,10.*PI/6.,12.*PI/6.}; source gs_ast; gs_ast.setup("R",zenith,azimuth,90e3,N_GS,NP2); <>= gs_ast.cleanup(); @ The wavefront sensor of the LGS asterism are setup next <>= centroiding gs_ast_cog; gs_ast_cog.setup(N_SIDE_LENSLET,N_GS); gs_ast_cog.MASK_SET = fried_geometry_setup(gs_ast_cog.lenslet_mask, pupil_mask.m, N_SIDE_LENSLET, 16, threshold); <>= gs_ast_cog.cleanup(); @ The 3 source are propagated through the atmosphere to the wavefront sensor. <>= atm.get_phase_screen(&gs_ast,N_GS,d,NP,d,NP,0); //gs_ast.phase2file("gsAstWavefronts.bin"); atm.get_phase_screen_gradient(&gs_ast_cog,N_SIDE_LENSLET,d,&gs_ast,N_GS,tau); @ The science wavefront is reconstructed with the tomographic estimator: <>= LMMSE gs_ast_lmmse; gs_ast_lmmse.setup(&atm,&gs_ast,N_GS,&src,1,d,N_SIDE_LENSLET,&pupil_mask,"MINRES"); int ltao_cvgce_iteration=0; float ltao_elapsed_time=0.; <>= gs_ast_lmmse.estimation(&gs_ast_cog); ltao_elapsed_time += gs_ast_lmmse.elapsed_time; ltao_cvgce_iteration += gs_ast_lmmse.iSolve.cvgce_iteration; @ The tomographic wavefront is subtracted from the science wavefront: <>= src.wavefront.reset(phase_screen); src.wavefront.add_phase(-1,gs_ast_lmmse.d__phase_est); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"phases/LMMSE LTAO"); src.wavefront.show_phase(plotly_name); } <>= float ltao_wfe_rms=0.; <>= ltao_wfe_rms += S.var(src.wavefront.phase, &pupil_mask, NP2); <>= gs_ast_lmmse.cleanup(); @ and the residual wavefront corresponding image is computed. <>= imager_ltao.propagate(&src); if (k_SAMPLE==PLOTLY_LIM) { SET_PLOTLY_NAME(plotly_name,plotly_dir,"frames/LMMSE LTAO"); imager_ltao.show_frame(plotly_name,&imager); } @ \subsection{Results} \label{sec:results} <>= science_wf_rms = 1E9*sqrtf(science_wf_rms/N_SAMPLE); ngsao_wfe_rms = 1E9*sqrtf(ngsao_wfe_rms/N_SAMPLE); lgs_wf_rms = 1E9*sqrtf(lgs_wf_rms/N_SAMPLE); lgs_wfe_rms = 1E9*sqrtf(lgs_wfe_rms/N_SAMPLE); lgs_wfe_rms_lmmse = 1E9*sqrtf(lgs_wfe_rms_lmmse/N_SAMPLE); ltao_wfe_rms = 1E9*sqrtf(ltao_wfe_rms/N_SAMPLE); ngsao_elapsed_time /= N_SAMPLE; ngsao_cvgce_iteration /= N_SAMPLE; lgsao_elapsed_time /= N_SAMPLE; lgsao_cvgce_iteration /= N_SAMPLE; ltao_elapsed_time /= N_SAMPLE; ltao_cvgce_iteration /= N_SAMPLE; printf("------------------------------\n"); printf("___ TURBULENCE WAVEFRONT ___\n"); printf("\n NGS WF RMS: %7.2fnm\n", science_wf_rms); printf("\n LGS WF RMS: %7.2fnm\n", lgs_wf_rms); printf("\n___ ON-AXIS WAVEFRONT ESTIMATE FROM NGS ___\n"); printf("\n NGSAO WFE RMS: %8.3fnm in %.2fms with %d iterations\n", ngsao_wfe_rms,ngsao_elapsed_time,ngsao_cvgce_iteration); printf("\n___ ON-AXIS WAVEFRONT ESTIMATE FROM LGS ___\n"); printf("\n NGS-LGS WFE RMS : %8.3fnm\n",lgs_wfe_rms); printf("\n LGS MMSE WFE RMS: %8.3fnm in %.2fms with %d iterations\n", lgs_wfe_rms_lmmse,lgsao_elapsed_time,lgsao_cvgce_iteration); printf("\n___ ON-AXIS WAVEFRONT ESTIMATE FROM LGS ASTERISM ___\n"); printf("\n WFE RMS: %8.3fnm in %.2fms with %d iterations\n", ltao_wfe_rms,ltao_elapsed_time,ltao_cvgce_iteration); printf("------------------------------\n");