fold_rfi

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  function r = fold_rfi(startMjd,endMjd,folddir)

   25-Jan-2012 (MAS): Created.


   Applies the folding procedure to a length of data, 
   eventually calculating the four histograms used for 
   choosing parameters.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

0001 function r = fold_rfi(startMjd,endMjd,basedir,folddir)
0002 
0003 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0004 %  function r = fold_rfi(startMjd,endMjd,folddir)
0005 %
0006 %   25-Jan-2012 (MAS): Created.
0007 %
0008 %
0009 %   Applies the folding procedure to a length of data,
0010 %   eventually calculating the four histograms used for
0011 %   choosing parameters.
0012 %
0013 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0014 
0015 % We need to use a reduction script for pipeline steps prior to RFI
0016 % flagging.
0017 % [home installeddir] = where_am_i();
0018 % redscript_paths = cell(3,1);
0019 % redscript_paths{1} = [home '/' installeddir '/rfi_tuning/redScript_preRFI_0.red'];
0020 % redscript_paths{2} = [home '/' installeddir '/rfi_tuning/redScript_preRFI_1.red'];
0021 % redscript_paths{3} = [home '/' installeddir '/rfi_tuning/redScript_preRFI_2.red'];
0022 
0023 
0024 % How many time sections will we use?
0025 nTimes = round(10 * (endMjd - startMjd));
0026 
0027 
0028 % Prep the cell array for output file names:
0029 outFileC = cell(nTimes,2);
0030 
0031 
0032 % Now loop over the time sections:
0033 %=================================
0034 for k=1:nTimes
0035 
0036     startSec = mjd2string(startMjd + 0.1 * (k-1));
0037     endSec = mjd2string(startMjd + 0.1 * k);
0038     
0039     disp(['Data from ' startSec ' to ' endSec]);
0040 
0041     % Generate the output file names.
0042     outFile_f = [folddir startSec '_filtered.mat'];
0043     outFile_p = [folddir startSec '_polonly.mat'];
0044 
0045     outFileC{k,1} = outFile_f;
0046     outFileC{k,2} = outFile_p;
0047     
0048     
0049     % Check if an output file for this time period already exists...
0050     %----------------------------------
0051     if exist(outFile_p,'file')
0052         disp('Time period already ran... on to next one.');
0053         continue;
0054     end
0055     
0056     disp(nTimes);
0057     disp(k);
0058     
0059     % Read in the requested data.
0060     %----------------------------------
0061     d = read_arc(startSec,endSec); 
0062 
0063     
0064     % Now we're going to apply alpha and load corrections for the
0065     % 'FILTERED' and 'POLONLY' modes.
0066     %----------------------------------
0067     disp('Applying alpha database and load correction to the data...');
0068     d1 = assembleAlphaStreams(d,'FILTERED');
0069     d1 = applyAlpha(d1,'FILTERED');
0070     d1 = load_template_corr(d1,0,0);
0071     
0072     d2 = assembleAlphaStreams(d,'POLONLY'); % 8
0073     d2 = applyAlpha(d2,'POLONLY');
0074     d2 = load_template_corr(d2,2,0);
0075     
0076     % And construct a new data vector containing I1, Q3, U3, and I2 from
0077     % each of the two modes:
0078     d.antenna0.receiver.data = ...
0079         [d1.antenna0.receiver.data(:,[1 6 7 8]) d2.antenna0.receiver.data(:,[1 6 7 8])];
0080     
0081     
0082     % Use the folding technique to look for RFI in the filtered data.
0083     %----------------------------------
0084     disp('Folding the data...');
0085     foldResult = fold_eval(d1);
0086 
0087     
0088     % Save some memory.
0089     clear d1;
0090     clear d2;
0091 
0092 
0093     % Get rid of noise diode events.
0094     %----------------------------------
0095     d = noiseRemove(d);
0096     
0097     
0098     % Do this for FILTERED and then for POLONLY data...
0099     disp('Comparing to statistical flagging...');
0100     for l=1:2
0101         switch l
0102             case 1
0103                 disp('FILTERED data...');
0104             case 2
0105                 disp('POLONLY data...');
0106         end
0107         
0108         
0109         % OK.  Now let us compare with RFI flagging stats ...
0110         %----------------------------------
0111         disp('Running statistical flagging...');
0112         [dataUtc, shortI_A, shortI_B, shortL_A, shortL_B, ...
0113             ~, longI_B, longI_B_min, longI_B_max, longL_A, longL_B] = ...
0114             stat_eval(d.antenna0.receiver.utc,d.antenna0.receiver.data(:,[1 2 3 4] + 4*(l-1)));    
0115 
0116 
0117         % Now produce the colormaps and histograms for the plots:
0118         %----------------------------------
0119         disp('Constructing the Histograms...');
0120         [s_A s_B l_A l_B] = stat_hist(foldResult, ...
0121             shortI_A, shortI_B, shortL_A, shortL_B, ...
0122         longI_B, longI_B_min, longI_B_max, longL_A, longL_B);
0123 
0124 
0125 
0126         % Now write the relevant info to output file.
0127         %----------------------------------
0128         disp('Writing data for this time series...');
0129 
0130         % Prepare the output...
0131 
0132         dataUtc = dataUtc(foldResult > 1);    
0133         falseNum = sum(foldResult == 1);
0134 
0135         switch l
0136             case 1
0137                 outFile = outFile_f;
0138             case 2
0139                 outFile = outFile_p;
0140             otherwise
0141                 % Shouldn't happen.
0142                 outFile = [folddir startSec '.mat'];
0143         end
0144         
0145         
0146         save(outFile, 's_A', 's_B', 'l_A', 'l_B', 'falseNum', 'dataUtc');
0147     end
0148     
0149 
0150     % Clean up and move along...
0151     clear d;
0152 end
0153 
0154 
0155 % Now compile the outputs into a three files...
0156 %----------------------------------
0157 disp('Compiling the folding outputs into three files...')
0158 fold_compile(basedir,outFileC(:,1),'filtered');
0159 % fold_compile(basedir,outFileC{:,k},'classic');
0160 fold_compile(basedir,outFileC(:,2),'polonly');
0161 
0162 r = 1;
0163 
0164 end