Home > cbassSouthFunctions > JLFunctions > jl_gridder.m

jl_gridder

PURPOSE ^

SYNOPSIS ^

function jl_gridder(object_string,start_date,stop_date)

DESCRIPTION ^

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function jl_gridder(object_string,start_date,stop_date)  
0002 
0003 save_dir='comm_rasters';
0004 save_path = [save_dir,'/',start_date,'_',object_string]
0005 which_roach=3;
0006 cell_size= 0.05;
0007 n=0;
0008 plot_contours=1;
0009 savefileon=1;
0010 map_size=16;
0011 
0012 disp(['Reducing a CBASS-S observation'])
0013 d = read_arcSouth(start_date,stop_date);
0014 d = determineIndicesSouth(d);
0015 
0016 d.index.radio_point_scan.fast = d.index.beammap.fast;
0017 d.index.radio_point_scan.slow = d.index.beammap.slow;
0018 
0019 if(which_roach==1)
0020   disp('reducing roach1')
0021 LL = [d.antenna0.roach1.LL]; 
0022 RR = [d.antenna0.roach1.RR]; 
0023 QQ = [d.antenna0.roach1.Q];
0024 UU = [d.antenna0.roach1.U];
0025 d.antenna0.receiver.utc = d.antenna0.roach1.utc;
0026 elseif (which_roach==2)
0027   disp('reducing roach2')
0028 LL = [d.antenna0.roach2.LL]; 
0029 RR = [d.antenna0.roach2.RR]; 
0030 QQ = [d.antenna0.roach2.Q];
0031 UU = [d.antenna0.roach2.U];
0032 d.antenna0.receiver.utc = d.antenna0.roach2.utc;
0033 elseif (which_roach==3)
0034 LL = [d.antenna0.roach2.LL,d.antenna0.roach1.LL]; 
0035 RR = [d.antenna0.roach2.RR,d.antenna0.roach1.RR]; 
0036 QQ = [d.antenna0.roach2.Q,d.antenna0.roach1.Q];
0037 UU = [d.antenna0.roach2.U,d.antenna0.roach1.U];
0038 d.antenna0.receiver.utc = d.antenna0.roach1.utc;
0039 end
0040 
0041 d.antenna0.receiver.data(:,1)=  mean(LL'); %LL
0042 d.antenna0.receiver.data(:,2)=  mean(RR'); %QQ
0043 d.antenna0.receiver.data(:,3)=  mean(QQ'); %UU
0044 d.antenna0.receiver.data(:,4)=  mean(UU'); %RR
0045 
0046 d.antenna0.servo.apparent(:,1)=d.antenna0.servo.az; 
0047 d.antenna0.servo.apparent(:,2)=d.antenna0.servo.el;
0048 
0049 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0050 % The basic deglitch routine
0051 %
0052 no_stds_to_reject=700;
0053 d.antenna0.receiver.data(:,3)=-1*d.antenna0.receiver.data(:,3);
0054 d.antenna0.receiver.data(:,4)=-1*d.antenna0.receiver.data(:,4);
0055 for chan=1:4;
0056 d_mean=mean(d.antenna0.receiver.data(:,chan));
0057 d_rms=rms(d.antenna0.receiver.data(:,chan));
0058 flagarray =( (d.antenna0.receiver.data(:,chan)>(d_mean-no_stds_to_reject*(d_rms-d_mean))) & (d.antenna0.receiver.data(:,chan)<(d_mean+no_stds_to_reject*(d_rms-d_mean))));
0059  d.antenna0.receiver.data(~flagarray,chan)=nan; 
0060 end
0061 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0062 
0063 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0064 % Plot the data we have
0065 hFig = figure(1);
0066 set(hFig, 'Visible', 'off');
0067 %figure(1)
0068 for i=1:4
0069 subplot(5,1,i)
0070 plot(d.antenna0.receiver.utc,d.antenna0.receiver.data(:,i))
0071 if (i==1) 
0072 title('Mean LL');
0073 end
0074 if (i==2)
0075   title('Mean RR');
0076 end
0077 
0078 if (i==3)
0079   title('Mean QQ');
0080 end
0081 if (i==4)
0082   title('Mean UU');
0083 end
0084 
0085 end
0086 % Plot wether it thinks it is a raster or not
0087 subplot(5,1,5)
0088 plot(d.antenna0.receiver.utc,d.index.radio_point_scan.fast);
0089 title('d.index.radio\_point\_scan.fast');
0090 
0091 gtitle('Raw data')
0092 saveas(hFig,[save_path,'_raw.png'])
0093 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0094 % Plot the bigger version of data(1)
0095 
0096 hFig = figure();
0097 set(hFig, 'Visible', 'off')
0098 plot(d.antenna0.receiver.utc,d.antenna0.receiver.data(:,1));
0099 title('Mean LL');
0100 %points = (dgood.index.radio_point_scan.fast==1);
0101 
0102 
0103 [s e] = findStartStop(d.index.radio_point_scan.fast);
0104 
0105 points = (d.index.radio_point_scan.fast==1);
0106 contour_max(1) = max(d.antenna0.receiver.data(points,1));
0107 contour_max(2) = max(d.antenna0.receiver.data(points,2));
0108 contour_max(3) = max(d.antenna0.receiver.data(points,3));
0109 contour_max(4) = max(d.antenna0.receiver.data(points,4));
0110 saveas(hFig,[save_path,'_raw_big.png'])
0111 %keyboard;
0112 %disp('Length of S is');
0113 %length(s)
0114 
0115 for m= 1:length(s)
0116 %   m
0117     n = n+1;
0118     ind = zeros(size(d.index.radio_point_scan.fast));
0119     ind(s(m):e(m)) = 1;
0120     ind = logical(ind);
0121 
0122     dcut  = framecut(d, ind);
0123 
0124     azApp = interp1(dcut.antenna0.tracker.utc, ...
0125         dcut.antenna0.tracker.horiz_topo(:,1), dcut.antenna0.receiver.utc);
0126    % keyboard;
0127     azOffSave = wrap180(azApp) - dcut.antenna0.servo.apparent(:,1);
0128     azOffSave = -azOffSave;
0129     % Note that dcut.antenna0.servo.el = dcut.antenna0.servo.apparent(:,2)
0130     azOffSave = azOffSave.*cos(dcut.antenna0.servo.el*pi/180);% Scale azimuth by cos(elevation)???CHECK
0131     
0132     %keyboard;
0133 
0134     elApp = interp1(dcut.antenna0.tracker.utc, ...
0135         dcut.antenna0.tracker.horiz_topo(:,2), dcut.antenna0.receiver.utc);
0136     elOffSave = elApp - dcut.antenna0.servo.apparent(:,2);
0137     elOffSave = -elOffSave;
0138 
0139     angle = sqrt(elOffSave.^2 + azOffSave.^2);
0140  
0141 
0142        
0143     % Make sure we ignore NaNs
0144     a = (isnan(azOffSave) | isnan(elOffSave)| isnan(dcut.antenna0.receiver.data(:,1))  );
0145     
0146     
0147     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0148     % Calculate the mean parallactic angle (NOT USED)
0149     AZ = deg2rad(mean(azApp(~a)));
0150     EL = deg2rad(mean(elApp(~a)));
0151     
0152     LAT_D = 37.233889; % OVRO Latitude in Degrees
0153     LON_D = 118.28222; % OVRO Longitude in Degrees
0154     
0155     LAT = deg2rad(LAT_D);
0156     LON = deg2rad(LON_D);
0157     
0158     DEC = deg2rad(mean((dcut.antenna0.tracker.equat_geoc(:,2))));
0159     sinPA = ((sin(AZ).*cos(LAT))./cos(DEC)); % radians
0160     sinPA(abs(sinPA)>=1)=1;
0161     PA = asin(sinPA);
0162     PA_D(m)= rad2deg(PA); % deg
0163     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0164 
0165     
0166     %Define a grid in degrees
0167     ti=-map_size/2:cell_size:map_size/2;
0168     
0169     % Make sure we ignore NaNs
0170     a = (isnan(azOffSave) | isnan(elOffSave) | isnan(dcut.antenna0.receiver.data(:,1)) );
0171 
0172 %   keyboard;
0173     disp('About to grid data');
0174          
0175     
0176     for (chan=1:1) % only plot the I channel data for speed
0177         if(size(dcut.antenna0.receiver.data(~a,chan),1)>0)
0178             % Grid the data
0179         disp('Gridding data');
0180             [XI,YI]=meshgrid(ti,ti);
0181             ZI = griddata(azOffSave(~a),elOffSave(~a),dcut.antenna0.receiver.data(~a,chan),XI,YI);
0182 
0183         disp ('about to plot contours');
0184             hFig = figure();
0185             set(hFig, 'Visible', 'off');
0186 
0187             
0188             %subplot(2,1,1)
0189             colormap(hot(256))
0190             %imagesc(ti,ti,10*log10(ZI));
0191             imagesc(ti,ti,ZI);
0192             xlabel('Angular offset (degrees) in Az direction');
0193             ylabel('Angular offset (degrees) in El.');
0194             title('Mean RR map');             
0195             saveas(hFig,[save_path,'_',num2str(m),'_map.png'])
0196             %figure
0197             %subplot(2,1,2)
0198             %contourf(XI,YI,ZI);
0199             hFig = figure();
0200             set(hFig, 'Visible', 'off'); 
0201             
0202 %            plot(azOffSave(~a),10*log10(dcut.antenna0.receiver.data(~a,chan))  )
0203              title('Mean RR Az scans');
0204              plot(azOffSave(~a),dcut.antenna0.receiver.data(~a,chan))
0205             saveas(hFig,[save_path,'_',num2str(m),'_az_scans.png'])
0206 %            if(plot_contours)
0207 %          disp('plot contours');
0208 %                % Plot mini contour maps on 3x3 grid on each page
0209 %                fig_no = ceil(m/16);
0210 %                subplot_no = m-(fig_no-1)*16;
0211 %
0212 %                %subplot(nfigs_x,nfigs_y,n);
0213 %                %Make contours based on the max_value in the data set for I,Q,U,V
0214 %                VI = [1:-0.1:0]*contour_max(chan);
0215 %
0216 %                figure(fig_no+2);
0217 %                subplot(4,4,subplot_no)
0218 %                contour(XI,YI,ZI,VI);
0219 %                axis square;
0220 %                grid on;
0221 %                xlabel('Az Offset,deg');
0222 %                ylabel('El Offset,deg');
0223 %                title(num2str(m));
0224 %            end
0225 
0226 
0227         end
0228     end
0229 end
0230 disp('Finished');

Generated on Sun 14-Jun-2015 17:12:45 by m2html © 2005