function plot_RAIM( filePrefix, d35_19 )

% 
% plot_cycle_slips(filePrefix, d35_19)
%
% Input:
%  filePrefix - String that will be at the beginning of each plot filename
%  d35_19     - viewdat output of record 39:19
%
% viewdat -d35_19 -h --mat=Obs.mat Input.T04 
%  
% In Matlab:-
% Obs = load_viewdat('Obs.mat');
% plot_RAIM('MyString',Obs);
%
% Useful utility to plot RAIM flags for GPS, GLONASS ... observables.
%
%
% Inspired by plot_cycle_slips.m developed by Andrew Cartmell
% Copyright Trimble Inc
%
% Tested on Matlab: 9.2.0.556344 (R2017a)
%

colorArray(1,1:3)  = [0.0   0.5   0.0  ];
colorArray(2,1:3)  = [1.0   0.0   0.0  ];
colorArray(3,1:3)  = [1.0   0.75  0.0  ];
colorArray(4,1:3)  = [0.0   1.0   0.0  ];
colorArray(5,1:3)  = [0.0   1.0   1.0  ];
colorArray(6,1:3)  = [1.0   0.0   1.0  ];
colorArray(7,1:3)  = [0.0   0.0   0.0  ];
colorArray(8,1:3)  = [0.75  0.75  0.75 ];
colorArray(9,1:3)  = [1.0   0.0   0.6  ];

% Record 35:19 fields
key.WN       = 1;  % GPS Week Number
key.TIME     = 2;  % GPS Time [Sec]
key.CLKBIAS  = 3;  % Clock Offset [ms]
key.RAIM     = 4;  % RAIM INFO
key.SV       = 5;  % SV ID
key.FDMA     = 6;  % FDMA Channel (GLONASS)
key.channel  = 7;  % Receiver Channel
key.SAT_TYPE = 8;  % Satellite Type
key.SV_FLAGS = 9;  % Satellite Flags
key.AZ       = 10; % Azimuth [deg]
key.EL       = 11; % Elevation [deg]
key.FREQ     = 12; % blockType = Frequency Band
key.TRACK    = 13; % trackType = Track Type
key.PHASE    = 14; % Phase [Cycles]
key.RANGE    = 15; % Pseudorange [m]
key.DOPP     = 16; % Doppler [Hz]
key.CNO      = 17; % C/No [dBHz]
key.SLIP_CNT = 18; % Cycle Slip Counter
key.MEAS     = 19; % Measurement Flags
key.ANTENNA  = 20; % Antenna #
key.CSTATE   = 21; % CSTATE in rec35:19 - not st_iface*.h value

% RAIM flags within SV_FLAGS
flag.dSV_FLAGS_FAULT_MM_RAIM  = 64;
flag.dSV_FLAGS_FAULT_WLS_RAIM = 128;
flag.dSV_FLAGS_FAULT_KF_RAIM  = 256;
flag.dSV_FLAGS_FAULT_SBAS_MSG = 512;
flag.dSV_FLAGS_FAULT_CNAV     = 1024;
flag.dSV_FLAGS_FAULT_FFT_CORR = 2048;
flag.dSV_FLAGS_FAULT_RTX_POS  = 4096;
flag.dSV_FLAGS_FAULT_SYS_RAIM = 8192;

%
% A few options (could make parameters on the input)
%

% Plot the X-axis in hours (0) or seconds (1)
plot_sec_axis  = 1;
% Only valid when the x-axis is in hours, force the X-axis to 24 hours
force_to_one_day = 0;
% When set align the time to be seconds within the current day
align_to_day = 0;

% Take a guess at the data rate
drate = 1.0;
stime = d35_19(floor(length(d35_19)/2),key.TIME);
for tmpx=floor(length(d35_19)/2):floor(length(d35_19)/2)+100
  if ( stime ~= d35_19(tmpx,key.TIME) )
    drate = d35_19(tmpx,key.TIME) - stime;
    break;
  end
end

fprintf('Guessed data is rate %.2f\n',drate);

% Align in days and convert to seconds
day= floor(d35_19(1,key.TIME)/86400);

if(plot_sec_axis)
  Scale = 1;
else
  Scale = 3600;
end

if(align_to_day == 1)
  d35_19(:,key.TIME) = (d35_19(:,key.TIME) - day*86400); 
end

taxis = round( (d35_19(1,key.TIME):drate:d35_19(end,key.TIME)) * 1000 * Scale);

isOctave = exist('OCTAVE_VERSION', 'builtin') ~= 0;
for sy = [0, 1, 2, 3, 4, 10, 9] % RT_SatType_e
  if(isOctave)
    figure("visible","off");
  else
    figure;
  end

  clear figLabel;
  clear pltHandle
  clear plotType;
  clear plotStr;
  figIndex = 1;
  legendIndex = 1;

  clf; hold on; grid on;
  yind = 1;
  ytlbl = char(1,4);

  if ( sy == 1 )
    sv_range=120:158;
    SVString = 'SBAS';
    sigString = 'L1 C/A';
  elseif(sy  == 0)
    sv_range=1:32;
    SVString = 'GPS';
    sigString = 'L1 C/A';
  elseif(sy  == 2)
    sv_range=1:24;
    SVString = 'GLN';
    sigString = 'L1 C/A';
  elseif(sy  == 3)
    sv_range=1:36;
    SVString = 'GAL';
    sigString = 'E1';
  elseif(sy  == 4)
    sv_range=193:197;
    SVString = 'QZSS';
    sigString = 'L1 C/A';
  elseif(sy  == 10)
    sv_range=1:30;
    SVString = 'BDS';
    sigString = 'B1';
  elseif(sy  == 9)
    sv_range=1:7;
    SVString = 'IRNSS';
    sigString = 'L5 C/A';
  else
    sv_range=1:32;
    SVString = 'Unknown';
    sigString = 'Unknown';
  end

  for sv = sv_range
    ind = find((d35_19(:,key.SV)==sv)&(d35_19(:,key.SAT_TYPE)==sy));
    if ( ~isempty(ind) )
      iL1CA=[]; iL1C=[]; iL1E=[]; iL2E=[]; iL2C=[]; iL5=[];

      if ( sy == 0 )
        [iL1CA, iL1C, iL1E, iL2E, iL2C, iL5] = find_gps_ind(d35_19(ind,:));
        sychar = 'G';
      elseif ( sy == 1 )
        [iL1CA, iL5] = find_sbas_ind(d35_19(ind,:));
        sychar = 'S';
      elseif ( sy == 2 )
        [iL1CA, iL1E, iL2C, iL2E, iL5] = find_gln_ind(d35_19(ind,:));
        sychar = 'R';
      elseif ( sy == 3 )
        [iL1CA, iL1E, iL2E, iL5, ~] = find_gal_ind(d35_19(ind,:));
        sychar = 'E';
      elseif ( sy == 4 )
        [iL1CA, iL1C, iL1E, iL2C, iL5, ~] = find_qzss_ind(d35_19(ind,:));
        sychar = 'Q';
      elseif ( sy == 5 || sy == 7 || sy == 10 )
        [iL1CA, iL5, iL2E] = find_cmp_ind(d35_19(ind,:));
        sychar = 'B';
      elseif ( sy == 9 )
        [iL1CA, iL1C] = find_irnss_ind(d35_19(ind,:));
        sychar = 'I';
      end

      yoff = .0;
      if (~isempty(iL1CA))

        d35_19ms=round(d35_19(ind(iL1CA),key.TIME)*1000*Scale);
        [~,A,B]=intersect(taxis,d35_19ms);
        tmpx = NaN*taxis; tmpx(A)=d35_19(ind(iL1CA(B)),key.TIME) / Scale;
        if(exist('OCTAVE_VERSION', 'builtin') ~= 0)
          pltHandle{legendIndex} = plot(tmpx,ones(length(tmpx),1)*yind + yoff + 0.05,'b','linewidth',10);
        else
          pltHandle{legendIndex} = plot(tmpx,ones(length(tmpx),1)*yind + yoff,'b');
        end
        
        plotType(legendIndex) = 0;
        plotStr{legendIndex} = sigString;
        legendIndex = legendIndex + 1;
        
        flagBit = flag.dSV_FLAGS_FAULT_MM_RAIM;
        for i=1:8
          yoff = yoff + 0.1;

          % Get the SV flags
          cind = find(bitand(d35_19(ind(iL1CA),key.SV_FLAGS),flagBit));
          if(length(cind) > 0)
            pltHandle{legendIndex} = plot(d35_19(ind(iL1CA(cind)),key.TIME) / Scale,ones(length(cind),1)*yind + yoff,'o','color',colorArray(i,:));
            plotType(legendIndex) = yoff;
          
            if(flagBit == flag.dSV_FLAGS_FAULT_MM_RAIM)
              plotStr{legendIndex} = 'MM RAIM';
            elseif(flagBit == flag.dSV_FLAGS_FAULT_WLS_RAIM)
              plotStr{legendIndex} = 'WLS RAIM';
            elseif(flagBit == flag.dSV_FLAGS_FAULT_KF_RAIM)
              plotStr{legendIndex} = 'KF RAIM';
            elseif(flagBit == flag.dSV_FLAGS_FAULT_SBAS_MSG)
              plotStr{legendIndex} = 'SBAS MSG';
            elseif(flagBit == flag.dSV_FLAGS_FAULT_CNAV)
              plotStr{legendIndex} = 'CNAV';
            elseif(flagBit == flag.dSV_FLAGS_FAULT_FFT_CORR)
              plotStr{legendIndex} = 'FFT CORR';
            elseif(flagBit == flag.dSV_FLAGS_FAULT_RTX_POS)
              plotStr{legendIndex} = 'RTX POS';
            elseif(flagBit == flag.dSV_FLAGS_FAULT_SYS_RAIM)
              plotStr{legendIndex} = 'SYS RAIM';
            end
        
            fprintf('%s%3d %s Epochs = %d | Percent = %.4f\n', ...
                     SVString, ...
                     sv, ...
                     plotStr{legendIndex}, ...
                     length(cind), ...
                     100*length(cind)/length(iL1CA));

            legendIndex = legendIndex + 1;
          end

          flagBit = flagBit * 2;
        end

        ytlbl(yind,1:4) = sprintf('%c%3d',sychar,sv);
        yind = yind + 1;
      end
    end
  end
  set(gca(),'YTick',1:yind-1);
  set(gca(),'YTickLabel',ytlbl);

  if(plot_sec_axis)
    xlabel('Time [GPS Secs]');
  else
    if(force_to_one_day)
      a = axis;
      a(1) = 0;
      a(2) = 24;
      a(3) = 0;
      a(4) = yind;
      axis(a);
    end
    if(align_to_day == 1)
      xlabel('Time [Hour into GPS day]');
    else
      xlabel('Time [Hour into GPS Week]');
    end
  end
 
  if(exist('pltHandle'))
    % Find the unique data types and use to form the legend
    [C,ia,ic] = unique(plotType, 'first');
    if(exist('OCTAVE_VERSION', 'builtin') ~= 0)
      % in V4 octave the legend is really large and transparent.
      % Set outside to make the plot more readable
      legend([pltHandle{ia}], plotStr{ia}, 'location', 'eastoutside')
    else
      legend([pltHandle{ia}], plotStr{ia}, 'location', 'best')
    end
  end

  title([SVString,' Tracking']);

  if(exist('OCTAVE_VERSION', 'builtin') ~= 0)
    print('-dpng',[filePrefix,'-',SVString,'-RAIM.png'],'-S960,720');
  else
    print('-dpng',[filePrefix,'-',SVString,'-RAIM.png']);
  end
end



