makeMapArrays

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

  function makeMapArrays(mapFile, varFile, healpixCoordFile, outFileName)

  function takes the skymap and corresponding healpixCoordFile to split it
  up longitude to make for faster searching when decomposing the sky into
  the timestream.

   sjcm
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

0001 function makeMapArrays(mapFile, varFile, healpixCoordFile, outFileName)
0002 
0003 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0004 %
0005 %  function makeMapArrays(mapFile, varFile, healpixCoordFile, outFileName)
0006 %
0007 %  function takes the skymap and corresponding healpixCoordFile to split it
0008 %  up longitude to make for faster searching when decomposing the sky into
0009 %  the timestream.
0010 %
0011 %   sjcm
0012 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0013 
0014 % healpixCoordFiles can be obtained from:
0015 % http://lambda.gsfc.nasa.gov/toolbox/tb_pixelcoords.cfm
0016 
0017 d2 = fitsread(healpixCoordFile, 'BinTable');
0018 data = fitsread(mapFile, 'BinTable');
0019 vardata = fitsread(varFile, 'BinTable');
0020 
0021 % data file
0022 Imap = (data{1}');
0023 Imap = Imap(:);
0024 Qmap = data{2}';
0025 Umap = data{2}';
0026 Qmap = Qmap(:);
0027 Umap = Umap(:);
0028 Pmap = pyth(Qmap, Umap);
0029 
0030 % var data file -- this might not be right.
0031 Ivar = (vardata{1}' + vardata{2}' + vardata{3}');
0032 Ivar = Ivar(:);
0033 Qvar = (vardata{2}' + vardata{4}' + vardata{5}');
0034 Uvar = (vardata{3}' + vardata{5}' + vardata{6}');
0035 Qvar = Qvar(:);
0036 Uvar = Uvar(:);
0037 Pvar = Qvar + Uvar;
0038 
0039 % pixel correspondance
0040 latmap = d2{2};
0041 longmap = d2{1};
0042 
0043 % remove bad data
0044 indbad = Imap < -1e10 | Qmap < -1e10 | Umap < -1e10;
0045 Imap(indbad) = [];
0046 Umap(indbad) = [];
0047 Qmap(indbad) = [];
0048 Pmap(indbad) = [];
0049 Ivar(indbad) = [];
0050 Uvar(indbad) = [];
0051 Qvar(indbad) = [];
0052 Pvar(indbad) = [];
0053 latmap(indbad) = [];
0054 longmap(indbad) = [];
0055 [uniquelat, I, J] = unique(latmap);
0056 display('makeMapArrays:: splitting up map');
0057 % split it up first
0058 tic
0059 for m=1:length(uniquelat)
0060   longSplit{m} = single(longmap(J==m));
0061   Isplit{m}    = single(Imap(J==m));
0062   Psplit{m}    = single(Pmap(J==m));
0063   Ivarsplit{m} = single(Ivar(J==m));
0064   Pvarsplit{m} = single(Pvar(J==m));
0065 end
0066 toc
0067       
0068 eval(sprintf('save %s Isplit Psplit Ivarsplit Pvarsplit longSplit uniquelat', outFileName));
0069 
0070 % SJCM:  I checked whether or not I could save this in text format for ease
0071 % of CVS update checks, but it creates a file that's 776MB big instead of
0072 % 129MB.  So I'm just going to add the .mat file to the repository, and
0073 % change the names when I need to updateit.
0074 
0075 
0076 return;
0077 
0078 
0079       
0080