detection/subset_domain.m

   1 function red=subset_domain(w,varargin)
   2 % red=subset_domain(w)
   3 % find rectangular domain around fire with some user guidance
   4 % and convert fire arrival time to datenum
   5 %
   6 % input: w structure with fields
   7 %    fxlat, fxlong - latitude, longitude
   8 %    tign_g        - fire arrival time
   9 %    nfuel_cat     - fuel categories
  10 %    times         - time at simulation end as string
  11 %
  12 % ouput: red structure with fields
  13 %    ispan,jspan   - list i,j indices selected
  14 %    fxlat,fxlong  - coordinates on submesh
  15 %    tign_g        - fire arrival time on submesh
  16 %    min_lon,max_lon,min_lat,max_lat - selected box
  17 %    time          - times as datenum
  18 %    max_tign_g    - max fire arrival time in w, corresponds to times
  19 %    tign          - tign_g as datenum
  20 %    min_tign      - min
  21 %    max_tign      - max
  22 %    base_time     - base time for displays
  23 %
  24 % Changes made %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  25 % original had red=subset_domain(w,red)
  26 % changed to red=subset_domain(varargin)
  27 %
  28 % comment out margin
  29 % comment out bounds, take default_bounds{2}
  30 %
  31
  32 % arguments
  33 if nargin >=2
  34     force = varargin{1},
  35 else
  36     force = 0;
  37 end
  38
  39 %for lots of evaluations of a single fire:
  40 %load red.mat
  41
  42 sim.min_lat = min(w.fxlat(:));
  43 sim.max_lat = max(w.fxlat(:));
  44 sim.min_lon = min(w.fxlong(:));
  45 sim.max_lon = max(w.fxlong(:));
  46 sim.min_tign= min(w.tign_g(:));
  47 sim.max_tign= max(w.tign_g(:));
  48 act.x=find(w.tign_g(:)<sim.max_tign);
  49 act.min_lat = min(w.fxlat(act.x));
  50 act.max_lat = max(w.fxlat(act.x));
  51 act.min_lon = min(w.fxlong(act.x));
  52 act.max_lon = max(w.fxlong(act.x));
  53 margin=input_num('relative margin around the fire',0.5,force);
  54 min_lon=max(sim.min_lon,act.min_lon-margin*(act.max_lon-act.min_lon));
  55 min_lat=max(sim.min_lat,act.min_lat-margin*(act.max_lat-act.min_lat));
  56 max_lon=min(sim.max_lon,act.max_lon+margin*(act.max_lon-act.min_lon));
  57 max_lat=min(sim.max_lat,act.max_lat+margin*(act.max_lat-act.min_lat));
  58
  59 default_bounds{1}=[min_lon,max_lon,min_lat,max_lat];
  60 default_bounds{2}=[sim.min_lon,sim.max_lon,sim.min_lat,sim.max_lat];
  61 default_bounds{3}=[-112.75414 -112.43779 40.27268 40.50655]
  62 for i=1:length(default_bounds),fprintf('default bounds %i: %8.5f %8.5f %8.5f %8.5f\n',i,default_bounds{i});end
  63
  64 bounds=input_num('enter bounds [min_lon,max_lon,min_lat,max_lat] or number of bounds above (1)> ',1,force);
  65 if isempty(bounds),bounds=1;end
  66 if length(bounds)==1,
  67    bounds=default_bounds{bounds};
  68 end
  69 %take the second default_bounds
  70 %bounds = default_bounds{2};
  71 [ii,jj]=find(w.fxlong>=bounds(1) & w.fxlong<=bounds(2) & w.fxlat >=bounds(3) & w.fxlat <=bounds(4));
  72 ispan=min(ii):max(ii);
  73 jspan=min(jj):max(jj);
  74 if isempty(ispan) | isempty(jspan), error('selection empty'),end
  75
  76 % restrict simulation
  77
  78 red.max_tign=sim.max_tign;
  79 red.ispan=ispan;
  80 red.jspan=jspan;
  81 red.fxlat=w.fxlat(ispan,jspan);
  82 red.fxlong=w.fxlong(ispan,jspan);
  83 red.tign_g=w.tign_g(ispan,jspan);
  84 red.lfn=w.lfn(ispan,jspan);
  85 %not all w.mat files will have fmc_g, nfuel_cat
  86 if isfield(w,'nfuel_cat')
  87     red.nfuel_cat=w.nfuel_cat(ispan,jspan);
  88 end
  89 if isfield(w,'fmc_g')
  90     red.fmc_g = w.fmc_g(ispan,jspan);
  91 end
  92 if isfield(w,'fhgt')
  93     red.fhgt = w.fhgt(ispan,jspan);
  94 end
  95 if isfield(w,'ros')
  96       if isempty(w.ros)
  97           w.ros = ones(size(w.tign_g));
  98       end
  99     red.ros = w.ros(ispan,jspan);
 100 end
 101 red.min_lat = min(red.fxlat(:));
 102 red.max_lat = max(red.fxlat(:));
 103 red.min_lon = min(red.fxlong(:));
 104 red.max_lon = max(red.fxlong(:));
 105
 106 % convert tign_g to datenum
 107
 108
 109 red.end_datenum=datenum(char(w.times(:))'); % this time step end
 110 %red.end_time=w.dt*w.itimestep; % time from simulation start in seconds
 111 red.end_time = [];
 112 red.start_time=0;
 113 % for reading wrfinput file
 114 if isempty(red.end_time)
 115     red.end_time = max(red.tign_g(:));
 116 end
 117 red.start_datenum=red.end_datenum-red.end_time/(24*3600);
 118
 119 fprintf('simulation start seems to be %s\n',datestr(red.start_datenum,'dd-mmm-yyyy HH:MM:SS'));
 120
 121 red.max_tign_g=max(w.tign_g(:));
 122 %red.tign=(red.tign_g - red.max_tign_g)/(24*60*60) + red.time;
 123 red.tign=time2datenum(red.tign_g,red);  % the tign array, in datenum
 124 red.min_tign=min(red.tign(:));
 125 red.max_tign=max(red.tign(:));
 126 end