beam_2d=prim_beam(ad,diameter,apfunc,freq,boost)
Calculate a stack of primary beams using the given
aperture field function at the set of specified frequencies.
diameter in meters
eg: ad=calc_ad(20,256);
beam_2d=prim_beam(ad,3.5,'sza_apfield',30e9);
boost is optional arg to reduce the effect of aliasing at the expense
of extra computing time - defeult is 2 - 4 is better.0001 function beam_2d=prim_beam(ad,diameter,apfunc,freq,boost) 0002 % beam_2d=prim_beam(ad,diameter,apfunc,freq,boost) 0003 % 0004 % Calculate a stack of primary beams using the given 0005 % aperture field function at the set of specified frequencies. 0006 % 0007 % diameter in meters 0008 % 0009 % eg: ad=calc_ad(20,256); 0010 % beam_2d=prim_beam(ad,3.5,'sza_apfield',30e9); 0011 % 0012 % boost is optional arg to reduce the effect of aliasing at the expense 0013 % of extra computing time - defeult is 2 - 4 is better. 0014 0015 if(~exist('boost')) 0016 boost=2; 0017 end 0018 0019 for i=1:length(freq) 0020 beam_2d(:,:,i)=calc_beam(ad,diameter,apfunc,freq(i),boost); 0021 end 0022 0023 0024 return 0025 0026 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0027 function beam_2d=calc_beam(ad,diameter,apfunc,freq,boost) 0028 0029 even = ~(mod(ad.N_pix,2)); 0030 0031 center = floor(ad.N_pix/2) + 1; 0032 0033 % Calc wavelength for this band 0034 lambda=3e8/freq; 0035 0036 % Convert ad.del_u which is in wavelengths to delta_r for 0037 % a given frequency. 0038 delta_r=ad.del_u*lambda; 0039 0040 % Boost the uv space resolution (and hence the image space span) 0041 delta_r=delta_r/boost; 0042 0043 % We want an odd number 0044 N_r_val=2*floor((diameter(1)/2)/delta_r)+1; 0045 0046 r_val=-delta_r*(N_r_val/2-0.5):delta_r:delta_r*(N_r_val/2-0.5); 0047 0048 % Fill grid with r values 0049 [x,y]=meshgrid(r_val,r_val); 0050 r=sqrt(x.^2+y.^2); 0051 0052 % Get aperture field values 0053 % Pass antenna diameter as argument in case want to use it in the function 0054 a_val_2d=feval(apfunc,r,diameter); 0055 0056 % Place this array in the middle of a big array for FFT 0057 uv=zeros(ad.N_pix*boost); 0058 i=(ad.N_pix*boost/2+1)-(N_r_val/2-0.5); 0059 j=(ad.N_pix*boost/2+1)+(N_r_val/2-0.5); 0060 uv(i:j,i:j)=a_val_2d; 0061 0062 % Invert to get beam pattern 0063 im=f2i(ad,uv); 0064 0065 % Pull out central region where aliasing is hopefully negligible 0066 i=ceil((ad.N_pix/2)*(boost-1)+1); 0067 j=ceil((ad.N_pix/2)*(boost+1)); 0068 im=im(i:j,i:j); 0069 0070 % Beam is the square of the transform 0071 beam_2d=im.^2; 0072 0073 % Normalize zero angle to be unit efficiency 0074 beam_2d=beam_2d/beam_2d(center,center); 0075 return