cycling/perim2gran.m

   1 function p_struct = perim2gran(n,perim)
   2 %reads in a directory of perimieters and makes them into "satellite granules"
   3 %n - number of points from each perimeter to use
   4 %perim - path to the perimeter data
   5
   6
   7 %%% reading the perimeters
   8 if perim(end) == 'l'
   9     fprintf('Using a kml file for perimeters \n')
  10     %%%%%% convert kml file to a struct
  11     temp_struct = kml2struct(perim);
  12
  13
  14 elseif perim(end) == '/'
  15
  16     perim_dat = ['kml';'shp'];
  17     p_type = input_num('Type of perimeter file to use? (1) kml (2) shp', 2,1);
  18     fprintf('Reading %s files in  directory %s \n',perim_dat(p_type,:),perim);
  19     if p_type == 1
  20         d = dir([perim,'*.kml']);
  21     else
  22         a = shape2struct(perim);
  23     end
  24 end % if perim_path...
  25
  26 %count perims, gather information
  27 %find perimeters from perim file
  28 p_count = 0;
  29
  30 %n gives  number of points in the perimeter to use
  31 %n = 100;
  32 for i = 1:length(a)
  33     %i, a(i)
  34     if strcmp(a(i).Geometry,'Polygon')
  35         p_count = p_count + 1;
  36
  37         % get perimeter time
  38         a(i).p_string = a(i).Name(end-12:end);
  39         formatIn = 'yyyymmdd HHMM';
  40         %perim times are local, need to convert to UTC
  41         %patch fire
  42         zone_shift = 6;
  43         % cougar creek , camp fire
  44         if strcmp(a(i).Name(1:2),'ca') | strcmp(a(i).Name(1:2),'wa')
  45             %a(i).p_string = a(i).Name(end-12:end);
  46             zone_shift = 8;
  47             %formatIn = 'yyyymmdd HHMM';
  48             %a(i).p_string = a(i).p_string(end-12:end)
  49         end
  50         %datenum format as used by TIGN
  51         a(i).p_time = datenum(a(i).p_string,formatIn)+zone_shift/24;
  52
  53         %set decimate to an  postive integer to use just a subset of points
  54         %  in perimeter
  55         fprintf('Perimeter %s: %d points in the perimeter \n',a(i).p_string,length(a(i).Lon));
  56         %filter out the NaN
  57         latnan = find(isnan(a(i).Lat));
  58         a(i).Lat(latnan) = [];
  59         a(i).Lon(latnan) = [];
  60         lonnan = find(isnan(a(i).Lon));
  61         a(i).Lat(lonnan) = [];
  62         a(i).Lon(lonnan) = [];
  63         if length(a(i).Lon) ~= length(a(i).Lat)
  64             fprintf('Size mismatch between Lon and Lat\n')
  65         end
  66         if length(a(i).Lat) > n
  67             decimate = round(length(a(i).Lat)/n);
  68             lats = a(i).Lat(1:decimate:end);
  69             lons = a(i).Lon(1:decimate:end);
  70         else
  71             lats = a(i).Lat;
  72             lons = a(i).Lon;
  73         end
  74
  75         %find data inside of perimeter
  76 %         [x,y] = meshgrid(xa,ya);
  77 %         x = x(:);
  78 %         y = y(:);
  79 %         [in,on] = inpolygon(x,y,lons,lats);
  80 %         fires = logical(in+on);
  81 %         data = reshape(fires,n,n);
  82         %make all high confidence fires
  83 %         data = uint8(9.0*data);
  84 %         a(i).data = data;
  85 %         geotransform = [ a(i).BoundingBox(1,1) dx 0  a(i).BoundingBox(2,2) 0 dy];
  86 %         a(i).geotransform = geotransform;
  87 %         %save the file for use in data assimilation
  88 %         %save(a(i).TIF_name,'data','geotransform');
  89 %         %plot results
  90
  91         %store perimter structure
  92         nums = length(lats);
  93         p_struct(p_count).det = [0 0 1 1 1];
  94         p_struct(p_count).power = ones(1,nums)*50;
  95         p_struct(p_count).data = ones(1,nums)*9;
  96         p_struct(p_count).conf = ones(1,nums)*95;
  97         p_struct(p_count).time = a(i).p_time;
  98         p_struct(p_count).lat = lats;
  99         p_struct(p_count).lon = lons;
 100         p_struct(p_count).file = replace(a(i).Name,' ','_');
 101         %p_struct(p_count).Name = a(i).Name;
 102         p_struct(p_count).xlon = [];
 103         p_struct(p_count).xlat = [];
 104         p_struct(p_count).fxdata = [];
 105         p_struct(p_count).axis = [];
 106
 107     end
 108 end %for
 109 end % function