Commit 5980d194 authored by Andrey Filippov's avatar Andrey Filippov

Before starting elevations fitting

parent 440719d9
......@@ -727,9 +727,13 @@ min_str_neib_fpn 0.35
public boolean terr_um_en = true;
public double terr_um_sigma = 1.0;
public double terr_um_weight = 0.8;
public double terr_nan_tolerance = 0.001; // set undefined terrain (out of capture area) to nan, it is almost zero after UM
public int terr_nan_grow = 20; // grow undefined terrain after detection
public double terr_nan_tolerance = 0.001; // set undefined terrain (out of capture area) to nan, it is almost zero after UM. if !terr_um_en will be forced ==0.0
public int terr_nan_grow = 20; // grow undefined scenes terrain after detection
public int terr_shrink_veget = 20; // shrink accumulated vegetation for filling terrain
public int terr_shrink_terrain= 20; // shrink accumulated terrain after cutting shrunk vegetation (added to terr_shrink_veget)
public double terr_vegetation_over = 35; // initial vegetation over (hotter) filled terrain
public int terr_filter_veget = 10; // shrink+grow shrunk vegetation to remove small clusters
public boolean terr_tile_woi = true; // if false - use woi50;
public Rectangle terr_woi_enclos = new Rectangle(80, 210, 210, 230); // will be tiled, using width/height from woi_step;
public Rectangle terr_woi_step = new Rectangle(10,10,20,20);
......@@ -737,6 +741,12 @@ min_str_neib_fpn 0.35
public boolean terr_skip_exist= true; // skip existing woi/parameters, already saved.
public boolean terr_continue= false; // start from last/single
public boolean terr_fit_terr = true; // adjust terrain pixels
public boolean terr_fit_veget = true; // adjust vegetation pixels
public boolean terr_fit_alpha = true; // adjust vegetation alpha pixels
public boolean terr_fit_scenes = true; // adjust scene offsets (start from 0 always?)
public boolean terr_fit_elevations = false; // adjust elevation pixels (not yet implemented)
public int terr_min_scenes = 1; // minimal number of scenes (inside woi) vegetation pixel must influence
public int terr_min_samples_scene = 10;
public int terr_min_total_scenes = 2;
......@@ -1970,11 +1980,19 @@ min_str_neib_fpn 0.35
gd.addStringField ("Model state file", terr_model_state, 50, "Model vegetation source data (w/o extension).");
gd.addStringField ("Segments subdir", terr_segments_dir, 50, "Model vegetation source data (w/o extension).");
gd.addMessage ("Unsharp mask filter of the input data (currently not used as both DC and HF fitting are implemented simultaneously).");
gd.addCheckbox ("Enable UM", terr_um_en, "Enable unsharp mask filter.");
gd.addNumericField("UM sigma", terr_um_sigma, 5,7, "pix", "Unsharp mask sigma.");
gd.addNumericField("UM weight", terr_um_weight, 5,7, "", "Unsharp mask weight.");
gd.addNumericField("NaN tolerance", terr_nan_tolerance, 5,7,"", "Replace continuous almost zeros with NaN-s.");
gd.addNumericField("NaN grow", terr_nan_grow, 0,3,"pix", "Grow \"almost\" zero areas after detection (+2 - 8 directions, +1 - only ortho).");
gd.addMessage ("Preparaion of the initial LMA approximation");
gd.addNumericField("NaN tolerance", terr_nan_tolerance, 5,7,"", "Replace continuous almost zeros with NaN-s (ignored and forced to 0 if no UM filter).");
gd.addNumericField("NaN grow", terr_nan_grow, 0,3,"pix", "Grow zero/\"almost\" zero rendered scene areas (+2 - 8 directions, +1 - only ortho), replace with NaNs.");
gd.addNumericField("Shrink vegetation", terr_shrink_veget, 0,3,"pix", "Shrink detected vegetation areas for terrain separation.");
gd.addNumericField("Shrink terrain", terr_shrink_terrain, 0,3,"pix", "Shrink remaining detected terrain areas before filling with laplacian (terr_shrink_veget is added to start from the original size).");
gd.addNumericField("Vegetation hotter", terr_vegetation_over, 5,7,"", "Initially consider to be vegetation if it is this hotter than filled in terrain holes.");
gd.addNumericField("Vegetation shrink+grow", terr_filter_veget, 0,3,"pix", "Shrink + grow hot enough vegetation to filter from small clusters - needed for initial alpha.");
gd.addMessage ("Scan WOI parameters");
gd.addCheckbox ("Tiled WOI", terr_tile_woi, "Process tiled WOI (False - a single one).");
......@@ -1985,6 +2003,13 @@ min_str_neib_fpn 0.35
gd.addCheckbox ("Skip existing WOI", terr_skip_exist, "Skip already esisting WOIs during scanning.");
gd.addMessage ("LMA losses and goals for alpha, terrain and vegetation pixels");
gd.addCheckbox ("Adjust terrain", terr_fit_terr, "Adjust terrain pixels.");
gd.addCheckbox ("Adjust vegetation", terr_fit_veget, "Adjust vegetation pixels.");
gd.addCheckbox ("Adjust alpha", terr_fit_alpha, "Adjust vegetation alpha pixels.");
gd.addCheckbox ("Adjust scene offsets", terr_fit_scenes, "Adjust scene offsets (start from 0 always?).");
gd.addCheckbox ("Adjust elevation", terr_fit_elevations,"Adjust elevation pixels (not yet implemented).");
gd.addNumericField("Min influenced scenes",terr_min_scenes, 0,3, "", "Minimal number of scenes (inside woi) vegetation pixel must influence.");
gd.addNumericField("Minimal samples/scene",terr_min_samples_scene, 0,3,"","Minimal samples per scene used for fitting (skip scene if less).");
gd.addNumericField("Minimum scenes used", terr_min_total_scenes, 0,3,"", "Minimal total number of scenes used (skip segment if less).");
......@@ -2668,6 +2693,10 @@ min_str_neib_fpn 0.35
terr_nan_tolerance = gd.getNextNumber();// double
terr_nan_grow = (int) gd.getNextNumber();// int
terr_shrink_veget = (int) gd.getNextNumber();// int
terr_shrink_terrain =(int) gd.getNextNumber();// int
terr_vegetation_over = gd.getNextNumber();// double
terr_filter_veget = (int) gd.getNextNumber();// int
terr_tile_woi = gd.getNextBoolean();// boolean
terr_continue = gd.getNextBoolean();// boolean
......@@ -2675,6 +2704,12 @@ min_str_neib_fpn 0.35
terr_woi_step = stringToRectangle(gd.getNextString());// Rectangle
terr_woi_last = stringToRectangle(gd.getNextString());// Rectangle
terr_skip_exist = gd.getNextBoolean();// boolean
terr_fit_terr = gd.getNextBoolean();// boolean
terr_fit_veget = gd.getNextBoolean();// boolean
terr_fit_alpha = gd.getNextBoolean();// boolean
terr_fit_scenes = gd.getNextBoolean();// boolean
terr_fit_elevations = gd.getNextBoolean();// boolean
terr_min_scenes = (int) gd.getNextNumber();// int
terr_min_samples_scene=(int) gd.getNextNumber();// int
......@@ -3327,6 +3362,11 @@ min_str_neib_fpn 0.35
properties.setProperty(prefix+"terr_um_weight", terr_um_weight+""); // double
properties.setProperty(prefix+"terr_nan_tolerance", terr_nan_tolerance+""); // double
properties.setProperty(prefix+"terr_nan_grow", terr_nan_grow+""); // int
properties.setProperty(prefix+"terr_nan_grow", terr_nan_grow+""); // int
properties.setProperty(prefix+"terr_nan_grow", terr_nan_grow+""); // int
properties.setProperty(prefix+"terr_nan_tolerance", terr_nan_tolerance+""); // double
properties.setProperty(prefix+"terr_nan_grow", terr_nan_grow+""); // int
properties.setProperty(prefix+"terr_tile_woi", terr_tile_woi+""); // boolean
properties.setProperty(prefix+"terr_continue", terr_continue+""); // boolean
......@@ -3334,6 +3374,12 @@ min_str_neib_fpn 0.35
properties.setProperty(prefix+"terr_woi_step", rectangleToString(terr_woi_step)+""); // Rectangle
properties.setProperty(prefix+"terr_woi_last", rectangleToString(terr_woi_last)+""); // Rectangle
properties.setProperty(prefix+"terr_skip_exist", terr_skip_exist+""); // boolean
properties.setProperty(prefix+"terr_fit_terr", terr_fit_terr+""); // boolean
properties.setProperty(prefix+"terr_fit_veget", terr_fit_veget+""); // boolean
properties.setProperty(prefix+"terr_fit_alpha", terr_fit_alpha+""); // boolean
properties.setProperty(prefix+"terr_fit_scenes", terr_fit_scenes+""); // boolean
properties.setProperty(prefix+"terr_fit_elevations", terr_fit_elevations+""); // boolean
properties.setProperty(prefix+"terr_min_scenes", terr_min_scenes+""); // int
properties.setProperty(prefix+"terr_min_samples_scene", terr_min_samples_scene+""); // int
......@@ -4008,6 +4054,10 @@ min_str_neib_fpn 0.35
if (properties.getProperty(prefix+"terr_um_weight")!= null) terr_um_weight=Double.parseDouble(properties.getProperty(prefix+"terr_um_weight"));
if (properties.getProperty(prefix+"terr_nan_tolerance")!= null) terr_nan_tolerance=Double.parseDouble(properties.getProperty(prefix+"terr_nan_tolerance"));
if (properties.getProperty(prefix+"terr_nan_grow")!= null) terr_nan_grow=Integer.parseInt(properties.getProperty(prefix+"terr_nan_grow"));
if (properties.getProperty(prefix+"terr_shrink_veget")!= null) terr_shrink_veget=Integer.parseInt(properties.getProperty(prefix+"terr_shrink_veget"));
if (properties.getProperty(prefix+"terr_shrink_terrain")!= null) terr_shrink_terrain=Integer.parseInt(properties.getProperty(prefix+"terr_shrink_terrain"));
if (properties.getProperty(prefix+"terr_vegetation_over")!= null) terr_vegetation_over=Double.parseDouble(properties.getProperty(prefix+"terr_vegetation_over"));
if (properties.getProperty(prefix+"terr_filter_veget")!= null) terr_filter_veget=Integer.parseInt(properties.getProperty(prefix+"terr_filter_veget"));
if (properties.getProperty(prefix+"terr_tile_woi")!= null) terr_tile_woi=Boolean.parseBoolean(properties.getProperty(prefix+"terr_tile_woi"));
if (properties.getProperty(prefix+"terr_continue")!= null) terr_continue=Boolean.parseBoolean(properties.getProperty(prefix+"terr_continue"));
......@@ -4015,7 +4065,12 @@ min_str_neib_fpn 0.35
if (properties.getProperty(prefix+"terr_woi_step")!= null) terr_woi_step=stringToRectangle((String) properties.getProperty(prefix+"terr_woi_step"));
if (properties.getProperty(prefix+"terr_woi_last")!= null) terr_woi_last=stringToRectangle((String) properties.getProperty(prefix+"terr_woi_last"));
if (properties.getProperty(prefix+"terr_skip_exist")!= null) terr_skip_exist=Boolean.parseBoolean(properties.getProperty(prefix+"terr_skip_exist"));
if (properties.getProperty(prefix+"terr_fit_terr")!= null) terr_fit_terr=Boolean.parseBoolean(properties.getProperty(prefix+"terr_fit_terr"));
if (properties.getProperty(prefix+"terr_fit_veget")!= null) terr_fit_veget=Boolean.parseBoolean(properties.getProperty(prefix+"terr_fit_veget"));
if (properties.getProperty(prefix+"terr_fit_alpha")!= null) terr_fit_alpha=Boolean.parseBoolean(properties.getProperty(prefix+"terr_fit_alpha"));
if (properties.getProperty(prefix+"terr_fit_scenes")!= null) terr_fit_scenes=Boolean.parseBoolean(properties.getProperty(prefix+"terr_fit_scenes"));
if (properties.getProperty(prefix+"terr_fit_elevations")!= null) terr_fit_elevations=Boolean.parseBoolean(properties.getProperty(prefix+"terr_fit_elevations"));
if (properties.getProperty(prefix+"terr_min_scenes")!= null) terr_min_scenes=Integer.parseInt(properties.getProperty(prefix+"terr_min_scenes"));
if (properties.getProperty(prefix+"terr_min_samples_scene")!= null) terr_min_samples_scene=Integer.parseInt(properties.getProperty(prefix+"terr_min_samples_scene"));
if (properties.getProperty(prefix+"terr_min_total_scenes")!= null) terr_min_total_scenes=Integer.parseInt(properties.getProperty(prefix+"terr_min_total_scenes"));
......@@ -4658,6 +4713,10 @@ min_str_neib_fpn 0.35
imp.terr_um_weight = this.terr_um_weight;
imp.terr_nan_tolerance = this.terr_nan_tolerance;
imp.terr_nan_grow = this.terr_nan_grow;
imp.terr_shrink_veget = this.terr_shrink_veget;
imp.terr_shrink_terrain = this.terr_shrink_terrain;
imp.terr_vegetation_over = this.terr_vegetation_over;
imp.terr_filter_veget = this.terr_filter_veget;
imp.terr_tile_woi = this.terr_tile_woi;
imp.terr_continue = this.terr_continue;
......@@ -4666,6 +4725,12 @@ min_str_neib_fpn 0.35
imp.terr_woi_last = new Rectangle(this.terr_woi_last);
imp.terr_skip_exist = this.terr_skip_exist;
imp.terr_fit_terr = this.terr_fit_terr;
imp.terr_fit_veget = this.terr_fit_veget;
imp.terr_fit_alpha = this.terr_fit_alpha;
imp.terr_fit_scenes = this.terr_fit_scenes;
imp.terr_fit_elevations = this.terr_fit_elevations;
imp.terr_min_scenes = this.terr_min_scenes;
imp.terr_min_samples_scene = this.terr_min_samples_scene;
imp.terr_min_total_scenes = this.terr_min_total_scenes;
......
......@@ -24,8 +24,12 @@ import ij.io.FileSaver;
public class VegetationLMA {
public static final int TVAO_TERRAIN = 0;
public static final int TVAO_VEGETATION = 1;
public static final int TVAO_VEGETATION_ALPHA = 2;
public static final int TVAO_SCENE_OFFSET = 3;
public static final int TVAO_ALPHA = 2;
public static final int TVAO_ELEVATION = 3;
public static final int TVAO_SCENE_OFFSET = 4;
public static final int TVAO_TYPES = 5;
public static final String PAR_EXT = ".par-tiff";
public static final int DATA_SOURCE_HEAD = 0; // header row with { scene, full_indx, terr_indx, scale_indx}, always present
public static final int DATA_SOURCE_CORN_VEGET =1; // [4]/null Z-shape 4 corners vegetation, either >= as parameter index or -1 - (x + image_width*y) of the unmodified full index
......@@ -37,12 +41,17 @@ public class VegetationLMA {
public static final int DATA_SOURCE_HEAD_SINDEX =3; // header row frame offset subindex
public static final int DATA_SOURCE_ITEMS = DATA_SOURCE_NEIB + 1; // == 4 number or rows in data source element
public static final int DATA_SOURCE_HEAD_SIZE = DATA_SOURCE_HEAD_SINDEX+1;
public static final int [] EMPTY4= {-1,-1,-1,-1};
public static final int [] SAVE_TYPES = {TVAO_TERRAIN, TVAO_VEGETATION, TVAO_ALPHA, TVAO_SCENE_OFFSET};
public final Rectangle full;
public final int image_length;
public final int num_scenes;
public Rectangle woi = null;
public double [] vegetation_average; // full image average vegetation (with nulls)
public double [] vegetation_filtered; // same as vegetation_average, but has more NaN to be used to select initial terrain
public double [] vegetation_pull; // have nan far from vegetation, use to pull LMA to this
public double [] terrain_average; // full image average terrain (no nulls)
public double [][] terrain_rendered; // terrain rendered for scenes (has nulls)
public double [][] terrain_rendered_hf; // terrain rendered for scenes (has nulls)
......@@ -71,12 +80,15 @@ public class VegetationLMA {
public int [][] par_rindex; // parameter -> pair of type (0..3) and full window pixel index
private int num_pars_terrain;
private int num_pars_vegetation;
private int num_pars_vegetation_alpha; //== num_pars_vegetation
private int num_pars_alpha; //== num_pars_vegetation
private int num_pars_elevations; //== num_pars_vegetation
private int num_pars_scenes;
private int ind_pars_terrain = 0; // index of the first terrain parameter
private int ind_pars_vegetation; // index of the first vegetation parameter
private int ind_pars_vegetation_alpha; // index of the first alpha parameter
private int ind_pars_alpha; // index of the first alpha parameter
private int ind_pars_elevations; // index of the first elevation parameter
private int ind_pars_scenes; // index of the first scene parameter
private int num_pars_scenes;
private int [][][] data_source; // [index][windx][DATA_SOURCE_ITEMS]:
// 0:{ scene, full_indx, terr_indx, scale_indx} - always present as parameters
// 1: [4]/null Z-shape 4 corners vegetation, either >= as parameter index or -1 - (x + image_width*y) of the unmodified full index
......@@ -95,28 +107,37 @@ public class VegetationLMA {
private double weight_pure = 0;
private double [] weights; // normalized so sum is 1.0 for all - samples and extra regularization
public boolean start_warm_veget = true; // start with vegetation warmer than terrain
public boolean start_warm_veget = true; // start with vegetation warmer than terrain
public double terrain_warmest = 90; // pull vegetations to warm, terrain to cold
public double initial_split = 0.1; // pull vegetations to warm, terrain to cold
public double min_split_frac = 0.15; // minimal modality fraction to use split by temperature
public double terr_difference = 100; // pull vegetation to be this warmer
public double terr_pull_cold = 0; // pull vegetations to warm, terrain to cold
public double initial_split = 0.1; // pull vegetations to warm, terrain to cold
public double min_split_frac = 0.15; // minimal modality fraction to use split by temperature
public double terr_difference = 100; // pull vegetation to be this warmer
public double terr_pull_cold = 0; // pull vegetations to warm, terrain to cold
// next two just for saving?
public double hifreq_weight; // 22.5 0 - do not use high-freq. Relative weight of laplacian components
public double reg_weights; // fraction of the total weight used for regularization
public double alpha_loss = 0; // not used with cosine alpha
public double alpha_offset = 0; // if >0, start losses above 0.0 and below 1.0;
public double alpha_lpf = 0;
public boolean [] fits;
/*
public boolean fit_terr = true;
public boolean fit_veget = true;
public boolean fit_alpha = true;
public boolean fit_scenes = true;
public boolean fit_elevations = false;
*/
public double alpha_loss = 0; // not used with cosine alpha
public double alpha_offset = 0; // if >0, start losses above 0.0 and below 1.0;
public double alpha_lpf = 0;
public boolean alpha_piece_linear = true;
public double alpha_scale_avg = 1.0; // 1.1; // scale average alpha (around 0.5) when pulling to it
public double alpha_push = 12.0; // push from alpha==0.5
public double alpha_scale_avg = 1.0; // 1.1; // scale average alpha (around 0.5) when pulling to it
public double alpha_push = 12.0; // push from alpha==0.5
public double alpha_push_neutral = 0.8; // alpha point from which push (closer to opaque)
public double alpha_push_center = 1.5; // weight of center alpha pixel relative to each of the 4 ortho ones
public double alpha_mm_hole = 0.1; // NaN to disable. Local "almost minimum" (lower than this fraction between min and max neighbor) is not subject to alpha_lpf
public boolean alpha_en_holes = true; // Search for small semi-transparent holes, disable diffusion of local alpha minimums
public double alpha_diff_hole = 0.01; // minimal alpha difference between min and max neighbor to be considered a hole
public double alpha_push_center = 1.5; // weight of center alpha pixel relative to each of the 4 ortho ones
public double alpha_mm_hole = 0.1; // NaN to disable. Local "almost minimum" (lower than this fraction between min and max neighbor) is not subject to alpha_lpf
public boolean alpha_en_holes = true; // Search for small semi-transparent holes, disable diffusion of local alpha minimums
public double alpha_diff_hole = 0.01; // minimal alpha difference between min and max neighbor to be considered a hole
public boolean from_file = false;
......@@ -126,8 +147,8 @@ public class VegetationLMA {
public double veget_lpf = 0;
// when unsharp mask is applied , pulling to 0 (when alpha is 0 (for vegetation) or 1.0 (for terrain) makes sense
public double terr_pull0 = 0;
public double veget_pull0 = 0;
public double terr_pull0 = 0; // now - pull to filled terrain - terrain_average
public double veget_pull0 = 0; // now - pull to vegetation_pull (extended vegetation)
public double boost_parallax = 1;
......@@ -174,14 +195,16 @@ public class VegetationLMA {
diff_offsets = diff_mode;
num_scenes = vegetation_offsets.length;
this.vegetation_average = vegetation_average;
this.vegetation_filtered= vegetation_average; // obsolete
this.terrain_average = terrain_average;
this.terrain_rendered = terrain_rendered;
this.vegetation_offsets = vegetation_offsets;
tvao = new double[4][];
tvao = new double[TVAO_TYPES][];
tvao[TVAO_TERRAIN] = this.terrain_average.clone();
tvao[TVAO_VEGETATION] = this.vegetation_average.clone();
tvao[TVAO_VEGETATION_ALPHA] = new double [image_length]; // 0 - use terrain
tvao[TVAO_ALPHA] = new double [image_length]; // 0 - use terrain
tvao[TVAO_SCENE_OFFSET] = new double [num_scenes];
setupLaplacians(0.0); // double weight_diag)
}
......@@ -191,14 +214,19 @@ public class VegetationLMA {
num_scenes = vegetationModel.terrain_scenes_render.length;
vegetation_average =vegetationModel.vegetation_full; // vegetation_average_render;
vegetation_filtered=vegetationModel.vegetation_filtered; // same with more NaNs for initial selection
terrain_average = vegetationModel.terrain_filled; // terrain_average_render;
terrain_rendered = vegetationModel.terrain_scenes_render;
vegetation_pull = vegetationModel.vegetation_pull; // pull LMA to this level
terrain_average = vegetationModel.terrain_filled; // terrain filled under vegetation, use to pull terrain level to it
terrain_rendered = vegetationModel.terrain_scenes_render;
vegetation_offsets =vegetationModel.vegetation_warp;
tvao = new double[4][];
tvao = new double[TVAO_TYPES][];
tvao[TVAO_TERRAIN] = this.terrain_average.clone();
tvao[TVAO_VEGETATION] = this.vegetation_average.clone();
tvao[TVAO_VEGETATION_ALPHA] = new double [image_length]; // 0 - use terrain
// tvao[TVAO_ALPHA] = new double [image_length]; // 0 - use terrain
tvao[TVAO_SCENE_OFFSET] = new double [num_scenes];
tvao[TVAO_ALPHA] =getInitialAlpha(
vegetation_filtered, // double [] vegetation_filtered,
vegetationModel.initial_transparent, // double transparent,
vegetationModel.initial_opaque); // double opaque)
setupLaplacians(0.0); // double weight_diag)
diff_offsets = vegetationModel.diff_mode;
this.vegetationModel = vegetationModel; // to access scene names, directories, reference index
......@@ -216,6 +244,28 @@ public class VegetationLMA {
}
}
private static double [] getInitialAlpha(
double [] vegetation_filtered,
double transparent,
double opaque) {
final Thread[] threads = ImageDtt.newThreadArray();
final AtomicInteger ai = new AtomicInteger(0);
final double [] alpha = new double [vegetation_filtered.length];
// vegetation_filtered
// tvao[TVAO_ALPHA]
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int nPix = ai.getAndIncrement(); nPix < vegetation_filtered.length; nPix = ai.getAndIncrement()) {
alpha[nPix] = Double.isNaN( vegetation_filtered[nPix])? transparent: opaque;
}
}
};
}
ImageDtt.startAndJoin(threads);
return alpha;
}
public int prepareLMA(
......@@ -225,17 +275,21 @@ public class VegetationLMA {
final int min_total_scenes,
final int min_samples_scene, // 10
final int min_valid_pixels,
final boolean start_warm_veget, // start with vegetation warmer than terrain -> USED now by vegetation_filtered
final double terrain_warmest, // warmest terrain (above is initially vegetation) NOT USED
final double initial_split, // initial alpha: terrain 0.0+, vegetation 1.0-. USED
final double min_split_frac, // minimal modality fraction to use split by temperature NOT USED?
final double terr_difference, // pull vegetation to be this warmer
final double terr_pull_cold, // pull vegetation to warm, terrain to cold
final double default_alpha,
final double hifreq_weight, // 22.5 0 - do not use high-freq. Relative weight of laplacian components
final boolean fit_terr,
final boolean fit_veget,
final boolean fit_alpha,
final boolean fit_scenes,
final boolean fit_elevations,
final double reg_weights, // fraction of the total weight used for regularization
final double alpha_loss, // quadratic loss when alpha reaches -1.0 or 2.0
final double alpha_offset, // quadratic loss when alpha reaches -1.0 or 2.0
......@@ -258,15 +312,22 @@ public class VegetationLMA {
String parameters_read_path,
final int debugLevel) {
this.woi = woi;
double alpha_contrast = 1.0;
this.start_warm_veget = start_warm_veget;
this.terrain_warmest = terrain_warmest;
this.initial_split = initial_split;
this.min_split_frac = min_split_frac;
this.terr_difference = terr_difference;
this.terr_pull_cold = terr_pull_cold;
this.hifreq_weight = hifreq_weight; // 22.5 0 - do not use high-freq. Relative weight of laplacian components
fits = new boolean[TVAO_TYPES];
fits[TVAO_TERRAIN] = fit_terr;
fits[TVAO_VEGETATION] = fit_veget;
fits[TVAO_ALPHA] = fit_alpha;
fits[TVAO_SCENE_OFFSET] =fit_scenes;
fits[TVAO_ELEVATION] = fit_elevations;
this.reg_weights = reg_weights; // fraction of the total weight used for regularization
this.alpha_loss = alpha_loss;
this.alpha_offset = alpha_offset;
......@@ -319,9 +380,9 @@ public class VegetationLMA {
setupParametersIndices(valid_woi); // needs to know number of used scenes // all but scenes
alpha_neibs = getNeighbors(
TVAO_VEGETATION_ALPHA, // final int tvao, // TVAO_VEGETATION_ALPHA
ind_pars_vegetation_alpha, // final int ind_samples, // ind_pars_vegetation_alpha
num_pars_vegetation_alpha); // final int num_samples); // num_pars_vegetation_alpha
TVAO_ALPHA, // final int tvao, // TVAO_VEGETATION_ALPHA
ind_pars_alpha, // final int ind_samples, // ind_pars_vegetation_alpha
num_pars_alpha); // final int num_samples); // num_pars_vegetation_alpha
terr_neibs = getNeighbors(
TVAO_TERRAIN, // final int tvao, // TVAO_VEGETATION_ALPHA
ind_pars_terrain, // final int ind_samples, // ind_pars_vegetation_alpha
......@@ -343,12 +404,20 @@ public class VegetationLMA {
}
if (!keep_parameters) {
/*
setupParametersVector (
start_warm_veget,// final boolean start_warm_veget,// start with vegetation warmer than terrain
terrain_warmest, // final double terrain_warmest, // warmest terrain (above is initially vegetation)
initial_split, // final double initial_split, // initial alpha: terrain 0.0+, vegetation 1.0-.
min_split_frac, // final double min_frac, // minimal modality fraction to use split by temperature
default_alpha); // final double default_alpha // use in areas where both terrain and vegetation are available
*/
// new version - uses preset alpha for the full image
if (alpha_contrast == 0) {
setupParametersVector ();
} else {
setupParametersVector(alpha_contrast); // double alpha_contrast);
}
}
from_file = false;
if (parameters_read_path != null) {
......@@ -1052,7 +1121,7 @@ public class VegetationLMA {
if (cw != null) {
for (int i = 0; i < 4; i++ ) {
int ia=indx_alpha[i];
alpha[i] = cw[i] * ((ia >= 0) ? vector[ia]: tvao[TVAO_VEGETATION_ALPHA][-1-ia]);
alpha[i] = cw[i] * ((ia >= 0) ? vector[ia]: tvao[TVAO_ALPHA][-1-ia]);
sum_a += alpha[i];
}
double tw = Math.max(0, 1.0 - sum_a);
......@@ -1108,7 +1177,10 @@ public class VegetationLMA {
jt[i] = new double [weights.length]; // weights.length];
}
}
final int dbg_n = -69992;
int dbg_scene = -20;
final int dbg_n = -7402; // -69992;
final int fdbg_scene = dbg_scene;
final int dbg_cols=2,dbg_rows=2;
final Thread[] threads = ImageDtt.newThreadArray();
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
......@@ -1116,27 +1188,50 @@ public class VegetationLMA {
public void run() {
double [] vegetation = new double [4];
double [] alpha = new double [4];
// for (int n = ai.getAndIncrement(); n < y_vector.length; n = ai.getAndIncrement()) {
for (int n = ai.getAndIncrement(); n < data_source.length; n = ai.getAndIncrement()) { // low-frequency
if (n == dbg_n) {
System.out.println("getFxDerivs(): n="+n);
}
// int nscene = data_source[n][0][0];
// int indx = data_source[n][0][1];
double terrain = vector[data_source[n][0][2]];
int indx = data_source[n][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_FINDEX];
int tindx = data_source[n][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_TINDEX];
int scene = data_source[n][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SCENE];
if (scene == fdbg_scene) {
int drow = indx / full.width -woi.y;
int dcol = indx % full.width -woi.x;
if ((dcol < dbg_cols) && (drow < dbg_rows)) {
System.out.println("getFxDerivs(): scene="+scene+", drow="+drow+" dcol="+dcol);
}
}
double terrain = (tindx >= 0) ? vector[tindx] : tvao[TVAO_TERRAIN][indx];
double [] cw = corners_weights[n];
double d;
int [] indx_vegetation = data_source[n][1];
int [] indx_alpha = data_source[n][2];
double sum_v =0, sum_a =0;
int [] indx_vegetation = data_source[n][DATA_SOURCE_CORN_VEGET];
int [] indx_alpha = data_source[n][DATA_SOURCE_CORN_ALPHA];
double sum_v =0, sum_a =0, sum_vw =0, sum_aw =0;
if (cw != null) {
for (int i = 0; i < 4; i++ ) {
int iv = indx_vegetation[i], ia=indx_alpha[i];
/*
double v = (iv >= 0) ? vector[iv]: tvao[TVAO_VEGETATION][-1-iv];
double a = (ia >= 0) ? vector[ia]: tvao[TVAO_ALPHA][-1-ia];
if (!Double.isNaN(v)) {
sum_vw += cw[i];
sum_v +=cw[i] * v;
}
if (!Double.isNaN(a)) {
sum_aw += cw[i];
sum_a +=cw[i] * a;
}
*/
vegetation[i] = cw[i] * ((iv >= 0) ? vector[iv]: tvao[TVAO_VEGETATION][-1-iv]);
alpha[i] = cw[i] * ((ia >= 0) ? vector[ia]: tvao[TVAO_VEGETATION_ALPHA][-1-ia]);
alpha[i] = cw[i] * ((ia >= 0) ? vector[ia]: tvao[TVAO_ALPHA][-1-ia]);
sum_v += vegetation[i];
sum_a += alpha[i];
}
}
if ((cw != null) && !Double.isNaN(sum_v) && !Double.isNaN(sum_a)){
double k;
if (alpha_piece_linear) {
k = (sum_a < 0)? 0: ((sum_a > 1) ? 1.0: sum_a);
......@@ -1144,42 +1239,43 @@ public class VegetationLMA {
k = (sum_a < 0)? 0: ((sum_a > 1) ? 1.0: 0.5 * (1.0 - Math.cos(sum_a*Math.PI)));
}
// d = terrain * (1.0 - sum_a) + sum_v * sum_a;
d = terrain * (1.0 - k) + sum_v * k;
if (jt != null) {
jt[data_source[n][0][2]][n] = 1 - k;// sum_a; // d/dterrain
if (tindx >= 0) {
jt[tindx][n] = 1 - k;// sum_a; // d/dterrain
}
for (int i = 0; i < 4; i++ ) {
if (indx_vegetation[i] >= 0) {
jt[data_source[n][1][i]][n] = cw[i] * k; // sum_a; // d/dvegetation[i]
if ((indx_vegetation != null) && (indx_vegetation[i] >= 0)) {
jt[data_source[n][DATA_SOURCE_CORN_VEGET][i]][n] = cw[i] * k; // sum_a; // d/dvegetation[i]
}
if ((indx_alpha[i] >= 0) && (sum_a > 0) && (sum_a < 1.0)) {
if ((indx_alpha != null) && (indx_alpha[i] >= 0) && (sum_a > 0) && (sum_a < 1.0)) {
if (alpha_piece_linear) {
jt[data_source[n][2][i]][n] = cw[i] * (sum_v - terrain); // d/dalpha[i]
jt[data_source[n][DATA_SOURCE_CORN_ALPHA][i]][n] = cw[i] * (sum_v - terrain); // d/dalpha[i]
} else {
jt[data_source[n][2][i]][n] = cw[i] * (sum_v - terrain) *0.5*Math.PI *Math.sin(sum_a*Math.PI); // d/dalpha[i]
jt[data_source[n][DATA_SOURCE_CORN_ALPHA][i]][n] = cw[i] * (sum_v - terrain) *0.5*Math.PI *Math.sin(sum_a*Math.PI); // d/dalpha[i]
}
}
}
}
} else {
d = terrain;
if (jt != null) {
jt[data_source[n][0][2]][n] = 1; // d/dterrain
d = terrain;
if ((jt != null) &&(tindx >= 0)) {
jt[tindx][n] = 1;// sum_a; // d/dterrain
}
}
if (data_source[n][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX] >= 0) {
if (data_source[n][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX] >= 0) { // can only be reference scene? use offset = 0 here.
double scene_offs = vector[data_source[n][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX]];
fX[n] = d + scene_offs;
if (Double.isNaN(fX[n])) {
System.out.println("fX["+n+"]= NaN");
}
d += scene_offs;
if (jt != null) {
jt[data_source[n][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX]][n] = 1;
}
}
fX[n] = d;
if (Double.isNaN(fX[n])) {
System.out.println("fX["+n+"]= NaN");
}
}
}
};
......@@ -1192,8 +1288,6 @@ public class VegetationLMA {
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
// double [] vegetation = new double [4];
// double [] alpha = new double [4];
for (int n_lf = ai.getAndIncrement(); n_lf < data_source.length; n_lf = ai.getAndIncrement()) {
if (n_lf == dbg_n) {
System.out.println("getFxDerivs(): n_lf="+n_lf);
......@@ -1202,33 +1296,24 @@ public class VegetationLMA {
int [] neibs = data_source[n_lf][DATA_SOURCE_NEIB];
if (neibs != null) { // nothing to do if null
double d_lf = fX[n_lf]; // center value
// int nscene = data_source[n_lf][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SCENE];
double swd = 0, sw=0;
for (int dir = 0; dir < 4; dir++) {
int nindx = neibs[dir];
if (nindx != -1) {
double d_hf;
if (nindx >= 0) {
d_hf = fX[neibs[dir]]; // neighbor value
} else {
d_hf = tvao[TVAO_TERRAIN][-2-nindx];
// add same scene offset as is added to all fX[] values so it will not influence differences
if (data_source[n_lf][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX] >= 0) {
double scene_offs = vector[data_source[n_lf][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX]];
d_hf += scene_offs;
}
}
if (nindx >=0) {
//see, maybe discard all border tiles
double d_hf = fX[neibs[dir]]; // neighbor value
swd += d_hf;
sw += 1;
}
}
if (sw > 0) { // nothing to do if 0
// if (sw > 0) { // nothing to do if 0
if (sw >=4) { // try discarding all with partial neighbors
swd /= sw;
fX[n_hf] = d_lf - swd;
if (jt != null) {
// center derivative same as _lf with respect to all it depends on except scene offset
int pindx =data_source[n_lf][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_TINDEX];
if (pindx >= 0) jt [pindx][n_hf] = jt[pindx][n_lf]; // d/dterrain always?
if (pindx >= 0) jt [pindx][n_hf] = jt[pindx][n_lf]; // d/dterrain always? not anymore
for (int i = 0; i < 4; i++ ){
if ( data_source[n_lf][DATA_SOURCE_CORN_VEGET] != null) {
pindx = data_source[n_lf][DATA_SOURCE_CORN_VEGET][i];
......@@ -1272,16 +1357,16 @@ public class VegetationLMA {
// splitting alpha_lpf from alpha_loss+alpha_push
int ind_next = y_vector.length;
if ((alpha_loss > 0) || (alpha_push > 0)) {
if (fits[TVAO_ALPHA] && ((alpha_loss > 0) || (alpha_push > 0))) {
int dbg_nx = -76340;
final int ind_y_alpha_loss = ind_next;
ind_next += num_pars_vegetation_alpha;
ind_next += num_pars_alpha;
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int n = ai.getAndIncrement(); n < num_pars_vegetation_alpha; n = ai.getAndIncrement()) {
int np = ind_pars_vegetation_alpha + n; // index of the alpha parameter
for (int n = ai.getAndIncrement(); n < num_pars_alpha; n = ai.getAndIncrement()) {
int np = ind_pars_alpha + n; // index of the alpha parameter
int nx = n + ind_y_alpha_loss; // y_vector.length; // x - index
if (nx == dbg_nx) {
System.out.println("getFxDerivs(): n="+n+", nx="+nx);
......@@ -1315,7 +1400,7 @@ public class VegetationLMA {
if (d > neib_max) neib_max = d;
nn++;
} else if (di < -1) {
d = tvao[TVAO_VEGETATION_ALPHA][-di - 2];
d = tvao[TVAO_ALPHA][-di - 2];
avg+=d;
if (d < neib_min) neib_min = d;
if (d > neib_max) neib_max = d;
......@@ -1375,16 +1460,16 @@ public class VegetationLMA {
if (alpha_lpf >= 0) {
if (fits[TVAO_ALPHA] && (alpha_lpf >= 0)) {
int dbg_nx = -76340;
final int ind_y_alpha_lpf = ind_next;
ind_next += num_pars_vegetation_alpha;
ind_next += num_pars_alpha;
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int n = ai.getAndIncrement(); n < num_pars_vegetation_alpha; n = ai.getAndIncrement()) {
int np = ind_pars_vegetation_alpha + n; // index of the alpha parameter
for (int n = ai.getAndIncrement(); n < num_pars_alpha; n = ai.getAndIncrement()) {
int np = ind_pars_alpha + n; // index of the alpha parameter
int nx = n + ind_y_alpha_lpf; // y_vector.length; // x - index
if (nx == dbg_nx) {
System.out.println("getFxDerivs(): n="+n+", nx="+nx);
......@@ -1404,7 +1489,7 @@ public class VegetationLMA {
if (d > neib_max) neib_max = d;
nn++;
} else if (di < -1) {
d = tvao[TVAO_VEGETATION_ALPHA][-di - 2];
d = tvao[TVAO_ALPHA][-di - 2];
avg+=d;
if (d < neib_min) neib_min = d;
if (d > neib_max) neib_max = d;
......@@ -1465,7 +1550,7 @@ public class VegetationLMA {
ImageDtt.startAndJoin(threads);
} // if (alpha_lpf >= 0) {
if (terr_lpf >= 0) {
if (fits[TVAO_TERRAIN] && (terr_lpf >= 0)) {
final int ind_y_terr = ind_next;
ind_next += num_pars_terrain;
ai.set(0);
......@@ -1473,7 +1558,7 @@ public class VegetationLMA {
threads[ithread] = new Thread() {
public void run() {
for (int n = ai.getAndIncrement(); n < num_pars_terrain; n = ai.getAndIncrement()) {
int np = ind_pars_terrain + n; // index of the alpha parameter
int np = ind_pars_terrain + n; // index of the terrain parameter
int nx = n + ind_y_terr; // y_vector.length; // x - index
double d = 0;
if (terr_lpf > 0) {
......@@ -1493,7 +1578,16 @@ public class VegetationLMA {
}
}
avg /= nn; // average
fX[nx] += terr_lpf * (vector[np] - avg) + terr_pull0 * vector[np];
fX[nx] += terr_lpf * (vector[np] - avg); // + terr_pull0 * (vector[np]);
if (terr_pull0 > 0) {
int findx = data_source[np][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_FINDEX];
double terr_pull = terrain_average[findx]; // maybe use tvao[TVAO_TERRAIN][findx] ?
if (Double.isNaN(terr_pull)) {
terr_pull= 0; // should not happen
}
fX[nx] += terr_pull0 * (vector[np] - terr_pull);
}
if (jt != null) {
jt[np][nx] += terr_lpf + terr_pull0;
for (int i = 0; i < terr_neibs[n].length; i++) { // now 4, may be increased
......@@ -1510,7 +1604,7 @@ public class VegetationLMA {
}
ImageDtt.startAndJoin(threads);
}
if (veget_lpf >= 0) { // should be positive for pull0 and terr_pull_cold (difference between vegetation and terrain)
if (fits[TVAO_VEGETATION] && (veget_lpf >= 0)) { // should be positive for pull0 and terr_pull_cold (difference between vegetation and terrain)
final int ind_y_veget = ind_next;
ind_next += num_pars_vegetation;
ai.set(0);
......@@ -1538,7 +1632,16 @@ public class VegetationLMA {
}
}
avg /= nn; // average
fX[nx] += veget_lpf * (vector[np] - avg) + veget_pull0 * vector[np];
fX[nx] += veget_lpf * (vector[np] - avg); // + veget_pull0 * vector[np];
if (veget_pull0 > 0) {
int findx = data_source[np][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_FINDEX];
double veget_pull = vegetation_pull[findx]; // maybe use tvao[TVAO_TERRAIN][findx] ?
if (Double.isNaN(veget_pull)) {
veget_pull= 0; // should not happen unless too far from the vegetation
}
fX[nx] += veget_pull0 * (vector[np] - veget_pull);
}
if (jt != null) {
jt[np][nx] += veget_lpf + veget_pull0;
for (int i = 0; i < veget_neibs[n].length; i++) { // now 4, may be increased
......@@ -1600,23 +1703,41 @@ public class VegetationLMA {
return jt;
}
private int [] getParamDebugIndices(
private int [][] getParamDebugIndices(
int gap) {
int woi_width3=woi.width*3 + 2 * gap;
// int num_scene_rows = (int) Math.ceil(1.0*par_index[3].length/woi_width3);
int num_scene_rows = (int) Math.ceil(1.0*num_scenes/woi_width3); // some scenes have NaN
int [] indices = new int [woi_width3 * (woi.height + gap+ num_scene_rows)];
Arrays.fill(indices, -1);
int [][] indices = new int [woi_width3 * (woi.height + gap+ num_scene_rows)][2];
for (int i = 0; i < indices.length; i++) {
indices[i][0] = -1; // unused
}
// Arrays.fill(indices, new int [] {});
// int [] save_types = {TVAO_TERRAIN, TVAO_VEGETATION, TVAO_ALPHA, TVAO_SCENE_OFFSET};
for (int h = 0; h < woi.height; h++) {
for (int n = 0; n < 3; n++) {
for (int w = 0; w < woi.width; w++) {
int indx = (h + woi.y) * full.width + (w + woi.x);
indices [(h * woi_width3) + (woi.width + gap) * n + w] = par_index[n][indx]; // h * woi.width + w];
int pindx = par_index[SAVE_TYPES[n]][indx];
if (pindx >= 0) {
indices [(h * woi_width3) + (woi.width + gap) * n + w][0] = 0;
indices [(h * woi_width3) + (woi.width + gap) * n + w][1] = pindx;
} else {
indices [(h * woi_width3) + (woi.width + gap) * n + w][0] = SAVE_TYPES[n] + 1; // type + 1
indices [(h * woi_width3) + (woi.width + gap) * n + w][1] = indx;
}
}
}
}
for (int i = 0; i < par_index[3].length; i++) {
indices [((woi.height + gap) * woi_width3) + i] = par_index[3][i]; // some will be -1
for (int i = 0; i < par_index[TVAO_SCENE_OFFSET].length; i++) {
int pindx = par_index[TVAO_SCENE_OFFSET][i]; // some will be -1
if (pindx >= 0) {
indices [((woi.height + gap) * woi_width3) + i][0] = 0;
indices [((woi.height + gap) * woi_width3) + i][1] = pindx;
} else {
indices [((woi.height + gap) * woi_width3) + i][0] = TVAO_SCENE_OFFSET + 1; // type + 1
indices [((woi.height + gap) * woi_width3) + i][1] = i;
}
}
return indices;
}
......@@ -1685,14 +1806,20 @@ public class VegetationLMA {
int gap,
int [] wh,
String title) {
int [] indices = getParamDebugIndices (gap);
int [][] indices = getParamDebugIndices (gap);
int woi_width3=woi.width*3 + 2 * gap;
double [] dbg_img = new double [indices.length];
for (int i = 0; i < indices.length; i++) {
if (indices[i] < 0) {
if (indices[i][0] < 0) { // unused, no mapping to parameter of fixed value
dbg_img[i] = Double.NaN;
} else if (indices[i][0] > 0) {
dbg_img[i] = Double.NaN;
int type = indices[i][0] -1;
dbg_img[i] = tvao[type][indices[i][1]];
} else {
dbg_img[i]= vector[indices[i]];
dbg_img[i]= vector[indices[i][1]];
}
if (!Double.isNaN(dbg_img[i])) {
// scale alpha for the same image
int w = i % woi_width3;
int h = i / woi_width3;
......@@ -1744,7 +1871,7 @@ public class VegetationLMA {
if (vector == null) {
vector = parameters_vector;
}
int height = woi.height;
if (num_scenes > woi_length) { // very unlikely but still may happen
......@@ -1752,11 +1879,12 @@ public class VegetationLMA {
ext_length = woi.width * height;
}
// assuming
String [] titles = {"terrain","vegetation","alpha","offset","terrain-index","vegetation-index","alpha-index","offset-index"};
// int [] save_types = {TVAO_TERRAIN, TVAO_VEGETATION, TVAO_ALPHA, TVAO_SCENE_OFFSET};
String [] titles = {"terrain","vegetation","alpha","elevation", "offset","terrain-index","vegetation-index","alpha-index","elevation-index", "offset-index"};
double [][] data = new double [titles.length][ext_length];
for (int i = 0; i < 4; i++) Arrays.fill(data[i],Double.NaN);
int [][] indices = new int [4][ext_length];
for (int n = 0; n < 4; n++) {
for (int i = 0; i < TVAO_TYPES; i++) Arrays.fill(data[i],Double.NaN);
int [][] indices = new int [TVAO_TYPES][ext_length];
for (int n = 0; n < TVAO_TYPES; n++) {
for (int i = 0; (i < indices[n].length) && (i < par_index[n].length); i++) {
int indx = i;
if (n < TVAO_SCENE_OFFSET) {
......@@ -1767,15 +1895,29 @@ public class VegetationLMA {
int pindx = par_index[n][indx];
indices[n][i] = pindx;
if (pindx < 0) {
if (n == TVAO_SCENE_OFFSET) {
// if (n == TVAO_SCENE_OFFSET) {
// data[n][i] = Double.NaN;
// } else {
if (tvao[n] == null) {
data[n][i] = Double.NaN;
} else {
data[n][i] = tvao[n][indx];
}
/*
if (data[n] == null) {
System.out.println("data["+n+"] == null");
}
if (tvao[n] == null) {
System.out.println("tvao["+n+"] == null");
}
data[n][i] = tvao[n][indx];
// }
*/
} else {
data[n][i] = vector[pindx];
}
data[n+4][i] = indices[n][i];
data[n+TVAO_TYPES][i] = indices[n][i];
}
}
ImagePlus imp = ShowDoubleFloatArrays.makeArrays(
......@@ -1795,11 +1937,11 @@ public class VegetationLMA {
imp.setProperty("DIFF_OFFSETS", ""+diff_offsets);
imp.setProperty("NUM_PARS_TERRAIN", ""+num_pars_terrain);
imp.setProperty("NUM_PARS_VEGETATION", ""+num_pars_vegetation);
imp.setProperty("NUM_PARS_VEGETATION_ALPHA", ""+num_pars_vegetation_alpha);
imp.setProperty("NUM_PARS_VEGETATION_ALPHA", ""+num_pars_alpha);
imp.setProperty("NUM_PARS_SCENES", ""+num_pars_scenes);
imp.setProperty("IND_PARS_TERRAIN", ""+ind_pars_terrain);
imp.setProperty("IND_PARS_VEGETATION", ""+ind_pars_vegetation);
imp.setProperty("IND_PARS_VEGETATION_ALPHA", ""+ind_pars_vegetation_alpha);
imp.setProperty("IND_PARS_VEGETATION_ALPHA", ""+ind_pars_alpha);
imp.setProperty("IND_PARS_SCENES", ""+ind_pars_scenes);
......@@ -1833,6 +1975,15 @@ public class VegetationLMA {
return;
}
public int getSaveIndex(int tvao_type) {
for (int i = 0; i < SAVE_TYPES.length; i++) {
if (SAVE_TYPES[i] == tvao_type) return i;
}
return -1;
}
public double [][] restoreParametersFile( //FIXME: Not finished for real import !
String path,
Rectangle file_woi) {// if not null, will return woi data and not input parameters to this instance
......@@ -1841,25 +1992,36 @@ public class VegetationLMA {
throw new IllegalArgumentException("Could not read "+path);
}
int num_slices = imp_pars.getStack().getSize();
if (num_slices != 8) {
throw new IllegalArgumentException("Expecting an 8-slice file, got "+num_slices);
boolean no_elevation = num_slices == 8;
boolean has_elevation = num_slices == 10;
if (!no_elevation && !has_elevation) {
throw new IllegalArgumentException("Expecting an 10- or 8-slice file, got "+num_slices);
}
int num_types = has_elevation? 5 : 4;
int width = imp_pars.getWidth();
int height = imp_pars.getHeight();
double [][] data = new double[4][width*height];
int [][] indices = new int[4][width*height];
for (int n = 0; n < 4; n++) {
double [][] data = new double[TVAO_TYPES][width*height];
int [][] indices = new int[TVAO_TYPES][width*height];
for (int n = 0; n < num_types; n++) {
int nt = has_elevation ? n : SAVE_TYPES[n];
float [] pixels = (float[]) imp_pars.getStack().getPixels(n+1);
for (int i = 0; i < data[n].length; i++) {
data[n][i] = pixels[i];
data[nt][i] = pixels[i];
}
}
for (int n = 0; n < 4; n++) {
float [] pixels = (float[]) imp_pars.getStack().getPixels(n+5);
for (int n = 0; n < num_types; n++) {
float [] pixels = (float[]) imp_pars.getStack().getPixels(n+1+ num_types);
int nt = has_elevation ? n : SAVE_TYPES[n];
for (int i = 0; i < data[n].length; i++) {
indices[n][i] = (int) pixels[i];
indices[nt][i] = (int) pixels[i];
}
}
if (!has_elevation) {
Arrays.fill(indices[TVAO_ELEVATION],-1);
Arrays.fill(data[TVAO_ELEVATION],Double.NaN);
}
JP46_Reader_camera.decodeProperiesFromInfo(imp_pars);
Rectangle full = new Rectangle (
0,
......@@ -1877,8 +2039,6 @@ public class VegetationLMA {
+num_scenes+".");
}
Rectangle woi = new Rectangle (
Integer.parseInt((String) imp_pars.getProperty("WOI_X")),
Integer.parseInt((String) imp_pars.getProperty("WOI_Y")),
......@@ -1891,7 +2051,7 @@ public class VegetationLMA {
file_woi.width = woi.width;
file_woi.height = woi.height;
// Remove undefined values that were not adjusted for this woi
for (int n = 0; n < 4; n++) {
for (int n = 0; n < TVAO_TYPES; n++) {
for (int i = 0; i < data[n].length; i++) {
if (indices[n][i] < 0) {
data[n][i] = Double.NaN;
......@@ -1902,7 +2062,7 @@ public class VegetationLMA {
// int num_pars_scenes = Integer.parseInt((String) imp_pars.getProperty("NUM_PARS_SCENES"));
// double [] scene_offsets = new double [num_pars_scenes];
double [] scene_offsets = new double [num_scenes];
System.arraycopy(data[TVAO_SCENE_OFFSET], 0, scene_offsets, 0, scene_offsets.length);
System.arraycopy(data[getSaveIndex(TVAO_SCENE_OFFSET)], 0, scene_offsets, 0, scene_offsets.length);
// temporary patch, later they will be NaNs already - Already done !
for (int i = 0; i < scene_offsets.length; i++) {
if (indices[TVAO_SCENE_OFFSET][i] < 0) {
......@@ -1910,7 +2070,7 @@ public class VegetationLMA {
}
}
for (int n = 0; n < 3; n++) {
for (int n = 0; n < SAVE_TYPES.length; n++) if (SAVE_TYPES[n] != TVAO_SCENE_OFFSET){
if (woi.width * woi.height < data[n].length ) {
double [] d = new double [woi.width * woi.height];
System.arraycopy(data[n], 0, d, 0, d.length);
......@@ -1924,11 +2084,11 @@ public class VegetationLMA {
alpha_piece_linear = Boolean.parseBoolean((String) imp_pars.getProperty("ALPHA_PIECE_LINEAR"));
num_pars_terrain = Integer.parseInt((String) imp_pars.getProperty("NUM_PARS_TERRAIN"));
num_pars_vegetation = Integer.parseInt((String) imp_pars.getProperty("NUM_PARS_VEGETATION"));
num_pars_vegetation_alpha = Integer.parseInt((String) imp_pars.getProperty("NUM_PARS_VEGETATION_ALPHA"));
num_pars_alpha = Integer.parseInt((String) imp_pars.getProperty("NUM_PARS_VEGETATION_ALPHA"));
num_pars_scenes = Integer.parseInt((String) imp_pars.getProperty("NUM_PARS_SCENES"));
ind_pars_terrain = Integer.parseInt((String) imp_pars.getProperty("IND_PARS_TERRAIN"));
ind_pars_vegetation = Integer.parseInt((String) imp_pars.getProperty("IND_PARS_VEGETATION"));
ind_pars_vegetation_alpha = Integer.parseInt((String) imp_pars.getProperty("IND_PARS_VEGETATION_ALPHA"));
ind_pars_alpha = Integer.parseInt((String) imp_pars.getProperty("IND_PARS_VEGETATION_ALPHA"));
ind_pars_scenes = Integer.parseInt((String) imp_pars.getProperty("IND_PARS_SCENES"));
if (imp_pars.getProperty("HIGHFREQ_WEIGHT") != null) hifreq_weight = Double.parseDouble((String) imp_pars.getProperty("HIGHFREQ_WEIGHT"));
......@@ -2480,7 +2640,7 @@ public class VegetationLMA {
String path,
double [] vector,
int gap) {
int [] indices = getParamDebugIndices (gap);
int [][] indices = getParamDebugIndices (gap);
ImagePlus imp_pars = new ImagePlus (path);
if (imp_pars.getWidth()==0) {
throw new IllegalArgumentException("Could not read "+path);
......@@ -2497,7 +2657,7 @@ public class VegetationLMA {
float [] pixels = (float[]) imp_pars.getStack().getPixels(num_slices); // last slice
for (int i = 0; i < indices.length; i++) {
if (indices[i] >= 0) {
if (indices[i][0] == 0) { // parameter, > 0 nonmodified type+1, then [1] is full index.
double d = pixels[i];
// scale alpha for the same image
int w = i % woi_width3;
......@@ -2505,7 +2665,7 @@ public class VegetationLMA {
if ((w > 2 *( woi.width+ gap)) && (h < woi.height)) {
d = (d- debug_alpha_scales[0])/(debug_alpha_scales[1] - debug_alpha_scales[0]);
}
vector[indices[i]] = d;
vector[indices[i][1]] = d;
}
}
from_file = true;
......@@ -2564,8 +2724,8 @@ public class VegetationLMA {
// using >=0 no use 0 as NOP but reserve space, <0 - do not reserve space
//(alpha_loss > 0)
if ((alpha_loss > 0) || (alpha_push > 0)) extra_samples+= num_pars_vegetation_alpha; // need to split loss (always positive) from alpha_lpf
if (alpha_lpf >= 0) extra_samples+= num_pars_vegetation_alpha;
if ((alpha_loss > 0) || (alpha_push > 0)) extra_samples+= num_pars_alpha; // need to split loss (always positive) from alpha_lpf
if (alpha_lpf >= 0) extra_samples+= num_pars_alpha;
if (terr_lpf >= 0) extra_samples+= num_pars_terrain;
if (veget_lpf >= 0) extra_samples+= num_pars_vegetation;
......@@ -2670,7 +2830,7 @@ public class VegetationLMA {
@SuppressWarnings("unused")
private int [][] getAlphaNeighbors(){
final int [][] alpha_neibs = new int[num_pars_vegetation_alpha][4];
final int [][] alpha_neibs = new int[num_pars_alpha][4];
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
......@@ -2678,16 +2838,16 @@ public class VegetationLMA {
threads[ithread] = new Thread() {
TileNeibs tn = new TileNeibs(full.width, full.height);
public void run() {
for (int i = ai.getAndIncrement(); i < num_pars_vegetation_alpha; i = ai.getAndIncrement()) {
int pindx = ind_pars_vegetation_alpha + i; // full parameter index
for (int i = ai.getAndIncrement(); i < num_pars_alpha; i = ai.getAndIncrement()) {
int pindx = ind_pars_alpha + i; // full parameter index
int indx = par_rindex[pindx][1]; // full pixel index
for (int dir4 = 0; dir4 < 4; dir4++) { // ortho only
int dir = 2 * dir4;
int nindx = tn.getNeibIndex(indx, dir);
if (nindx < 0) {
alpha_neibs[i][dir4] = -1;
} else if (par_index[TVAO_VEGETATION_ALPHA][nindx] >= 0) {
alpha_neibs[i][dir4] = par_index[TVAO_VEGETATION_ALPHA][nindx];
} else if (par_index[TVAO_ALPHA][nindx] >= 0) {
alpha_neibs[i][dir4] = par_index[TVAO_ALPHA][nindx];
} else { // pull to fixed alpha value
alpha_neibs[i][dir4] = -2 -indx; // OFFSET by 2 !
}
......@@ -2760,23 +2920,22 @@ public class VegetationLMA {
int col = woi.x + dcol;
int windx =dcol + drow * woi.width;
int indx = row*full.width + col;
if (!Double.isNaN(terrain_rendered[nScene][indx]) && (par_index[TVAO_TERRAIN][indx] >= 0)) { //
scene_samples[nScene]++;
data_src[nScene][windx] = new int[DATA_SOURCE_ITEMS][];
// maybe use windx instead of indx below?
// data_src[nScene][windx][0] = new int [] {nScene, indx, par_index[TVAO_TERRAIN][indx],nScene};
data_src[nScene][windx][DATA_SOURCE_HEAD] = new int [DATA_SOURCE_HEAD_SIZE];
data_src[nScene][windx][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SCENE] = nScene;
data_src[nScene][windx][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_FINDEX] = indx;
data_src[nScene][windx][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_TINDEX] = par_index[TVAO_TERRAIN][indx];
data_src[nScene][windx][DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX] = nScene;
// next 4 not needed - they are already nulls
// data_src[nScene][windx][1] = null;
// data_src[nScene][windx][2] = null;
data_src[nScene][windx][DATA_SOURCE_CORN_VEGET] = null;
data_src[nScene][windx][DATA_SOURCE_CORN_ALPHA] = null;
data_src[nScene][windx][DATA_SOURCE_NEIB] = null;
corn_w[nScene][windx] = null;
// now par_index[TVAO_TERRAIN][indx] may be -1
// if (!Double.isNaN(terrain_rendered[nScene][indx]) && (par_index[TVAO_TERRAIN][indx] >= 0)) { //
if (!Double.isNaN(terrain_rendered[nScene][indx])) { //
boolean has_data = false;
int [][] data_record = new int[DATA_SOURCE_ITEMS][];
data_record[DATA_SOURCE_HEAD] = new int [DATA_SOURCE_HEAD_SIZE];
data_record[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_TINDEX] = par_index[TVAO_TERRAIN][indx]; // may be -1
if (par_index[TVAO_TERRAIN][indx] >= 0) {
has_data=true;
} else if (fits[TVAO_TERRAIN]) {
continue; // needs terrain data
}
data_record[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SCENE] = nScene;
data_record[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_FINDEX] = indx;
data_record[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX] = nScene;
// try building record, if no data - do not add it
double [] veg_xy = vegetation_offsets[nScene][indx];
if (veg_xy != null) {
veg_xy = veg_xy.clone();
......@@ -2794,26 +2953,32 @@ public class VegetationLMA {
double fx = veg_xy[0] - x0;
double fy = veg_xy[1] - y0;
corn_w[nScene][windx] = new double []{(1-fx)*(1-fy), fx*(1-fy), (1-fx)*fy, fx*fy};
// data_src[nScene][windx][1] = new int [4];
// data_src[nScene][windx][2] = new int [4];
data_src[nScene][windx][DATA_SOURCE_CORN_VEGET] = new int [4];
data_src[nScene][windx][DATA_SOURCE_CORN_ALPHA] = new int [4];
data_record[DATA_SOURCE_CORN_VEGET] = new int [4];
data_record[DATA_SOURCE_CORN_ALPHA] = new int [4];
int veg_indx = x0+full.width*y0;
int [] veg_indx4 = {veg_indx, veg_indx+1, veg_indx+full.width, veg_indx+full.width+1};
for (int i = 0; i < 4; i++ ) {
// data_src[nScene][windx][1][i] = (par_index[TVAO_VEGETATION][veg_indx4[i]] >= 0)?
// par_index[TVAO_VEGETATION][veg_indx4[i]] : (-1-veg_indx4[i]);
// data_src[nScene][windx][2][i] = (par_index[TVAO_VEGETATION_ALPHA][veg_indx4[i]] >= 0)?
// par_index[TVAO_VEGETATION_ALPHA][veg_indx4[i]] : (-1-veg_indx4[i]);
data_src[nScene][windx][DATA_SOURCE_CORN_VEGET][i] = (par_index[TVAO_VEGETATION][veg_indx4[i]] >= 0)?
par_index[TVAO_VEGETATION][veg_indx4[i]] : (-1-veg_indx4[i]);
data_src[nScene][windx][DATA_SOURCE_CORN_ALPHA][i] = (par_index[TVAO_VEGETATION_ALPHA][veg_indx4[i]] >= 0)?
par_index[TVAO_VEGETATION_ALPHA][veg_indx4[i]] : (-1-veg_indx4[i]);
int pi = par_index[TVAO_VEGETATION][veg_indx4[i]];
if (pi >= 0) {
data_record[DATA_SOURCE_CORN_VEGET][i] = pi;
has_data = true;
} else {
data_record[DATA_SOURCE_CORN_VEGET][i] = -1-veg_indx4[i];
}
pi = par_index[TVAO_ALPHA][veg_indx4[i]];
if (pi >= 0) {
data_record[DATA_SOURCE_CORN_ALPHA][i] = pi;
has_data = true;
} else {
data_record[DATA_SOURCE_CORN_ALPHA][i] = -1-veg_indx4[i];
}
}
}
}
if (has_data) {
scene_samples[nScene]++;
data_src[nScene][windx] = data_record;
}
}
}
}
......@@ -2823,21 +2988,37 @@ public class VegetationLMA {
}
ImageDtt.startAndJoin(threads);
// now condense removing nulls in data_src (cornw may still have nulls)
int max_samples = 0;
for (int nscene = 0; nscene < num_scenes; nscene++) {
max_samples = Math.max(max_samples, scene_samples[nscene]);
}
System.out.println("Maximal samples in a scene = "+max_samples);
int allowed_less = 0;
int need_samples = max_samples-allowed_less;
int num_samples = 0;
// int dbg_early_sample = 0; // 10;
used_scenes = new boolean [num_scenes];
final int [] used_scene_indices = new int [num_scenes];
int num_used_scenes = 0;
for (int nscene = 0; nscene < num_scenes; nscene++) {
int num_prev = num_samples;
int defined = scene_samples[nscene];
boolean enough = defined >= need_samples;
if (!enough) {
System.out.println ("Not enough samples in scene "+nscene+" ("+defined+" < "+need_samples+"), removing it.");
}
used_scene_indices[nscene] = num_used_scenes;
if (defined >= min_samples_scene) {
// if ((defined >= min_samples_scene) && (nscene >= dbg_early_sample)) {
if (enough) {
used_scenes[nscene] = true;
num_samples+= scene_samples[nscene];
scene_samples[nscene] = num_prev; // start index
num_used_scenes++;
if (nscene != vegetationModel.reference_index) {
num_used_scenes++;
} else {
used_scene_indices[nscene] = -1;
}
}
}
data_source = new int [num_samples][][];
corners_weights = new double[num_samples][];
......@@ -2851,8 +3032,11 @@ public class VegetationLMA {
if (data_src[nScene][windx] != null) {
corners_weights[out_indx] = corn_w[nScene][windx];
int [][] dsrc = data_src[nScene][windx];
// dsrc[0][3] = ind_pars_scenes + used_scene_indices[dsrc[0][3]]; // replace scene number with the corresponding parameter index
dsrc[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX] = ind_pars_scenes + used_scene_indices[dsrc[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX]]; // replace scene number with the corresponding parameter index
int used_index = used_scene_indices[dsrc[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX]];
if (used_index >= 0) {
used_index += ind_pars_scenes;
}
dsrc[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX] = used_index; // ind_pars_scenes + used_scene_indices[dsrc[DATA_SOURCE_HEAD][DATA_SOURCE_HEAD_SINDEX]]; // replace scene number with the corresponding parameter index
data_source [out_indx++] = dsrc; // data_src[nScene][windx];
}
}
......@@ -2911,15 +3095,16 @@ public class VegetationLMA {
return;
}
@SuppressWarnings("unused")
private void setupParametersVector(
final boolean start_warm_veget,// start with vegetation warmer than terrain
final double terrain_warmest, // warmest terrain (above is initially vegetation)
final double initial_split, // initial alpha: terrain 0.0+, vegetation 1.0-.
final double min_frac, // minimal modality fraction to use split by temperature
final double default_alpha) { // areas where both terrain and vegetation are available
for (int i = 0; i < tvao[TVAO_VEGETATION_ALPHA].length; i++) {
for (int i = 0; i < tvao[TVAO_ALPHA].length; i++) {
if (!Double.isNaN(tvao[TVAO_VEGETATION][i])) {
tvao[TVAO_VEGETATION_ALPHA][i] = default_alpha;
tvao[TVAO_ALPHA][i] = default_alpha;
}
}
boolean split = false;
......@@ -2930,7 +3115,7 @@ public class VegetationLMA {
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int indx = row*full.width + col;
if (par_index[TVAO_VEGETATION_ALPHA][indx] >= 0) {
if (par_index[TVAO_ALPHA][indx] >= 0) {
if (par_index[TVAO_VEGETATION][indx] >= 0) {
if (!Double.isNaN(vegetation_filtered[indx])) {
/// if (tvao[TVAO_VEGETATION][indx] > terrain_warmest) {
......@@ -2968,21 +3153,22 @@ public class VegetationLMA {
System.out.println("parameters_vector[par_index["+TVAO_VEGETATION+"]["+indx+"]]=NaN");
}
}
if (par_index[TVAO_VEGETATION_ALPHA][indx] >= 0) {
if (par_index[TVAO_ALPHA][indx] >= 0) {
if (split && (par_index[TVAO_VEGETATION][indx] >= 0)) {
/// if (tvao[TVAO_VEGETATION][indx] > terrain_warmest) {
if (!Double.isNaN(vegetation_filtered[indx])) {
parameters_vector[par_index[TVAO_VEGETATION_ALPHA][indx]] = 1.0 - initial_split;
parameters_vector[par_index[TVAO_ALPHA][indx]] = 1.0 - initial_split;
} else {
parameters_vector[par_index[TVAO_VEGETATION_ALPHA][indx]] = initial_split;
parameters_vector[par_index[TVAO_ALPHA][indx]] = initial_split;
}
} else {
parameters_vector[par_index[TVAO_VEGETATION_ALPHA][indx]] = default_alpha;
parameters_vector[par_index[TVAO_ALPHA][indx]] = default_alpha;
}
if (Double.isNaN(parameters_vector[par_index[TVAO_VEGETATION_ALPHA][indx]])) {
System.out.println("parameters_vector[par_index["+TVAO_VEGETATION_ALPHA+"]["+indx+"]]=NaN");
if (Double.isNaN(parameters_vector[par_index[TVAO_ALPHA][indx]])) {
System.out.println("parameters_vector[par_index["+TVAO_ALPHA+"]["+indx+"]]=NaN");
}
}
// TODO: process TVAO_ELEVATION
}
}
for (int i = 0; i < par_index[TVAO_SCENE_OFFSET].length; i++) if (par_index[TVAO_SCENE_OFFSET][i] >= 0){
......@@ -2994,99 +3180,239 @@ public class VegetationLMA {
return;
}
private void setupParametersVector() { // areas where both terrain and vegetation are available
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int indx = row*full.width + col;
if (par_index[TVAO_TERRAIN][indx] >= 0) {
parameters_vector[par_index[TVAO_TERRAIN][indx]] = tvao[TVAO_TERRAIN][indx];
if (Double.isNaN(parameters_vector[par_index[TVAO_TERRAIN][indx]])) {
System.out.println("parameters_vector[par_index["+TVAO_TERRAIN+"]["+indx+"]]=NaN");
}
}
if (par_index[TVAO_VEGETATION][indx] >= 0) {
parameters_vector[par_index[TVAO_VEGETATION][indx]] = tvao[TVAO_VEGETATION][indx];
if (Double.isNaN(parameters_vector[par_index[TVAO_VEGETATION][indx]])) {
System.out.println("parameters_vector[par_index["+TVAO_VEGETATION+"]["+indx+"]]=NaN");
}
}
if (par_index[TVAO_ALPHA][indx] >= 0) {
parameters_vector[par_index[TVAO_ALPHA][indx]] = tvao[TVAO_ALPHA][indx];
if (Double.isNaN(parameters_vector[par_index[TVAO_ALPHA][indx]])) {
System.out.println("parameters_vector[par_index["+TVAO_ALPHA+"]["+indx+"]]=NaN");
}
}
// TODO: process TVAO_ELEVATION
}
}
for (int i = 0; i < par_index[TVAO_SCENE_OFFSET].length; i++) if (par_index[TVAO_SCENE_OFFSET][i] >= 0){
parameters_vector[par_index[TVAO_SCENE_OFFSET][i]]= tvao[TVAO_SCENE_OFFSET][i];
if (Double.isNaN(parameters_vector[par_index[TVAO_SCENE_OFFSET][i]])) {
System.out.println("parameters_vector[par_index["+TVAO_SCENE_OFFSET+"]["+i+"]]");
}
}
return;
}
private void setupParametersVector(double alpha_contrast) {
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int indx = row*full.width + col;
if (par_index[TVAO_TERRAIN][indx] >= 0) {
parameters_vector[par_index[TVAO_TERRAIN][indx]] = tvao[TVAO_TERRAIN][indx];
if (Double.isNaN(parameters_vector[par_index[TVAO_TERRAIN][indx]])) {
System.out.println("parameters_vector[par_index["+TVAO_TERRAIN+"]["+indx+"]]=NaN");
}
}
if (par_index[TVAO_VEGETATION][indx] >= 0) {
parameters_vector[par_index[TVAO_VEGETATION][indx]] = tvao[TVAO_VEGETATION][indx];
if (Double.isNaN(parameters_vector[par_index[TVAO_VEGETATION][indx]])) {
System.out.println("parameters_vector[par_index["+TVAO_VEGETATION+"]["+indx+"]]=NaN");
}
}
if (par_index[TVAO_ALPHA][indx] >= 0) {
parameters_vector[par_index[TVAO_ALPHA][indx]] = tvao[TVAO_ALPHA][indx];
if (Double.isNaN(parameters_vector[par_index[TVAO_ALPHA][indx]])) {
System.out.println("parameters_vector[par_index["+TVAO_ALPHA+"]["+indx+"]]=NaN");
}
}
// TODO: process TVAO_ELEVATION
}
}
if (alpha_contrast > 0) {
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int indx = row*full.width + col;
if (par_index[TVAO_ALPHA][indx] >= 0) {
double t = tvao[TVAO_TERRAIN][indx];
double v = tvao[TVAO_VEGETATION][indx];
double vp = vegetation_pull[indx];
double a = par_index[TVAO_ALPHA][indx];
if (!Double.isNaN(t) && !Double.isNaN(v) && !Double.isNaN(vp)) {
if (v < t) {
a = 0.0;
} else if (v > vp) {
a = 1.0;
} else {
a = (v-t)/(vp-t);
if (alpha_contrast != 1.0) {
a = 0.5 + (a - 0.5)*alpha_contrast;
if (a < 0) a = 0;
else if (a > 1.0) a = 1;
}
// if (par_index[TVAO_VEGETATION][indx] >= 0) {
// parameters_vector[par_index[TVAO_VEGETATION][indx]] = vp;
// }
}
if (par_index[TVAO_VEGETATION][indx] >= 0) {
parameters_vector[par_index[TVAO_VEGETATION][indx]] = vp;
}
parameters_vector[par_index[TVAO_ALPHA][indx]] = a;
}
}
}
}
}
for (int i = 0; i < par_index[TVAO_SCENE_OFFSET].length; i++) if (par_index[TVAO_SCENE_OFFSET][i] >= 0){
parameters_vector[par_index[TVAO_SCENE_OFFSET][i]]= tvao[TVAO_SCENE_OFFSET][i];
if (Double.isNaN(parameters_vector[par_index[TVAO_SCENE_OFFSET][i]])) {
System.out.println("parameters_vector[par_index["+TVAO_SCENE_OFFSET+"]["+i+"]]");
}
}
return;
}
private void setupParametersIndices(boolean [][] valid_woi) {
par_index = new int [4][];
private void setupParametersIndices(
boolean [][] valid_woi) { // [0] - valid terrain pixels, [1] valid vegetation and alpha pixels;
par_index = new int [TVAO_TYPES][];
par_index[TVAO_TERRAIN] = new int [image_length];
par_index[TVAO_VEGETATION] = new int [image_length];
par_index[TVAO_VEGETATION_ALPHA] = new int [image_length];
par_index[TVAO_ALPHA] = new int [image_length];
par_index[TVAO_SCENE_OFFSET] = new int [num_scenes];
par_index[TVAO_ELEVATION] = new int [image_length];
if (valid_woi == null) {
valid_woi = new boolean [2][woi.width*woi.height];
Arrays.fill(valid_woi[0],true);
Arrays.fill(valid_woi[1],true);
// mask with available vegetation data? Or later? Same/separate with warp/elevation?
}
int max_pars = 0;
for (int i = 0; i < par_index.length; i++) {
max_pars += par_index[i].length;
}
int [][] par_rindex = new int [max_pars][2];
Arrays.fill(par_index[TVAO_TERRAIN], -1);
Arrays.fill(par_index[TVAO_VEGETATION], -1);
Arrays.fill(par_index[TVAO_VEGETATION_ALPHA], -1);
Arrays.fill(par_index[TVAO_TERRAIN], -1);
Arrays.fill(par_index[TVAO_VEGETATION], -1);
Arrays.fill(par_index[TVAO_ALPHA], -1);
Arrays.fill(par_index[TVAO_ELEVATION], -1);
int par_num = 0;
ind_pars_terrain = 0;
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int windx =dcol + drow * woi.width;
if (valid_woi[0][windx]) {
int indx = row*full.width + col;
par_rindex[par_num][0] = 0;
par_rindex[par_num][1] = indx;
par_index[TVAO_TERRAIN][indx] = par_num++;
if (fits[TVAO_TERRAIN]) {
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int windx =dcol + drow * woi.width;
if (valid_woi[0][windx]) {
int indx = row*full.width + col;
par_rindex[par_num][0] = TVAO_TERRAIN;
par_rindex[par_num][1] = indx;
par_index[TVAO_TERRAIN][indx] = par_num++;
}
}
}
}
num_pars_terrain = par_num;
ind_pars_vegetation = par_num;
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int windx =dcol + drow * woi.width;
if (valid_woi[1][windx]) {
int indx = row*full.width + col;
par_rindex[par_num][0] = 1;
par_rindex[par_num][1] = indx;
par_index[TVAO_VEGETATION][indx] = par_num++;
if (fits[TVAO_VEGETATION]) {
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int windx =dcol + drow * woi.width;
if (valid_woi[1][windx]) {
int indx = row*full.width + col;
if (!Double.isNaN(vegetation_average[indx])) {
par_rindex[par_num][0] = TVAO_VEGETATION;
par_rindex[par_num][1] = indx;
par_index[TVAO_VEGETATION][indx] = par_num++;
}
}
}
}
}
ind_pars_vegetation_alpha = par_num;
ind_pars_alpha = par_num;
num_pars_vegetation = par_num - ind_pars_vegetation;
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int windx =dcol + drow * woi.width;
if (valid_woi[1][windx]) {
int indx = row*full.width + col;
par_rindex[par_num][0] = 2;
par_rindex[par_num][1] = indx;
par_index[TVAO_VEGETATION_ALPHA][indx] = par_num++;
if (fits[TVAO_ALPHA]) {
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int windx =dcol + drow * woi.width;
if (valid_woi[1][windx]) {
int indx = row*full.width + col;
if (!Double.isNaN(vegetation_average[indx])) {
par_rindex[par_num][0] = TVAO_ALPHA;
par_rindex[par_num][1] = indx;
par_index[TVAO_ALPHA][indx] = par_num++;
}
}
}
}
}
num_pars_alpha = par_num - ind_pars_alpha; // normally should be == num_pars_vegetation
ind_pars_elevations = par_num;
if (fits[TVAO_ELEVATION]) {
for (int drow = 0; drow < woi.height; drow++) {
int row = woi.y + drow;
for (int dcol = 0; dcol < woi.width; dcol++) {
int col = woi.x + dcol;
int windx =dcol + drow * woi.width;
if (valid_woi[1][windx]) {
int indx = row*full.width + col;
if (!Double.isNaN(vegetation_average[indx])) {
par_rindex[par_num][0] = TVAO_ELEVATION;
par_rindex[par_num][1] = indx;
par_index[TVAO_ELEVATION][indx] = par_num++;
}
}
}
}
}
num_pars_vegetation_alpha = par_num - ind_pars_vegetation_alpha; // normally should be == num_pars_vegetation
// int num_pars_pixel = par_num;
num_pars_elevations = par_num - ind_pars_elevations; // normally should be == num_pars_vegetation
ind_pars_scenes = par_num;
// temporarily assign for all scenes, later (when know which scenes are used) will be updated, so it must be the last (after TVAO_ELEVATION)
for (int i = 0; i < par_index[TVAO_SCENE_OFFSET].length; i++) {
par_rindex[par_num][0] = 3;
par_rindex[par_num][0] = TVAO_SCENE_OFFSET;
par_rindex[par_num][1] = i;
par_index[TVAO_SCENE_OFFSET][i] = par_num++;
}
// temporarily assign full rindex, will finalize in finishParametersIndices()
this.par_rindex = par_rindex;
/*
this.par_rindex = new int [par_num][2];
System.arraycopy(par_rindex, 0, this.par_rindex, 0, par_num);
num_pars_scenes = par_num - ind_pars_scenes;
parameters_vector = new double [par_num];
*/
return;
}
private void finishParametersIndices() {
int par_num = ind_pars_scenes;
int [][] par_rindex = this.par_rindex;
for (int i = 0; i < par_index[TVAO_SCENE_OFFSET].length; i++) {
if (used_scenes[i]) {
par_rindex[par_num][0] = 3;
if (used_scenes[i] && !isReferenceScene(i)) {
par_rindex[par_num][0] = TVAO_SCENE_OFFSET;
par_rindex[par_num][1] = i;
par_index[TVAO_SCENE_OFFSET][i] = par_num++;
} else {
......@@ -3097,9 +3423,12 @@ public class VegetationLMA {
System.arraycopy(par_rindex, 0, this.par_rindex, 0, par_num);
num_pars_scenes = par_num - ind_pars_scenes;
parameters_vector = new double [par_num];
return;
}
private boolean isReferenceScene(int scene) {
return vegetationModel.reference_index == scene;
}
private boolean [][] filterValidWoi (
int min_scenes_uses,
......
......@@ -57,10 +57,13 @@ public class VegetationModel {
public double [] terrain_filled = null;
public double [] vegetation_full = null; // ~same as vegetation_average_render
public double [] vegetation_filtered = null; // more NaNs than in vegetation_average_render
public double [] vegetation_pull = null; // filled around actual vegetation to use as LMA pull
public double initial_transparent= 0.0;
public double initial_opaque= 1.0;
public double [][][] vegetation_warp = null;
public double [][][] vegetation_warp_md = null; // magnitude+direction
public String [] scene_names = null; // TODO: Implement!
public boolean diff_mode = true;
public Rectangle full;
......@@ -1307,6 +1310,13 @@ public class VegetationModel {
double terr_pull_cold = clt_parameters.imp.terr_pull_cold; // pull vegetations to warm, terrain to cold
double default_alpha = 0.5; // 0.8;
double hifreq_weight = 22.5; // 0 - do not use high-freq. Relative weight of laplacian components
boolean fit_terr = clt_parameters.imp.terr_fit_terr; // true; // adjust terrain pixels
boolean fit_veget = clt_parameters.imp.terr_fit_veget; // true; // adjust vegetation pixels
boolean fit_alpha = clt_parameters.imp.terr_fit_alpha; // true; // adjust vegetation alpha pixels
boolean fit_scenes = clt_parameters.imp.terr_fit_scenes; // true; // adjust scene offsets (start from 0 always?)
boolean fit_elevations = clt_parameters.imp.terr_fit_elevations; // false; // adjust elevation pixels (not yet implemented)
double reg_weights = 0.25; // fraction of the total weight used for regularization
double alpha_loss = 100.0; // 10.0; /// 100.0; // 10.0; // 10000.0; // 1000.0; // 100.; // 10.0; // quadratic loss when alpha reaches -1.0 or 2.0
double alpha_offset = 0.0; // 0.02; // 0.03; // if >0, start losses above 0.0 and below 1.0;
......@@ -1499,6 +1509,11 @@ public class VegetationModel {
terr_pull_cold, // final double terr_pull_cold, // pull vegetations to warm, terrain to cold
default_alpha, // final double default_alpha,
hifreq_weight, //final double hifreq_weight, // 22.5 0 - do not use high-freq. Relative weight of laplacian components
fit_terr, // final boolean adjust_terr,
fit_veget, // final boolean adjust_veget,
fit_alpha, // final boolean adjust_alpha,
fit_scenes, // final boolean adjust_scenes,
fit_elevations, // final boolean adjust_elevations,
reg_weights, // final double reg_weights, // fraction of the total weight used for regularization
alpha_loss, // final double alpha_loss, // quadratic loss when alpha reaches -1.0 or 2.0
alpha_offset, // final double alpha_offset, // quadratic loss when alpha reaches -1.0 or 2.0
......@@ -1621,8 +1636,15 @@ public class VegetationModel {
boolean um_en = clt_parameters.imp.terr_um_en; // true;
double um_sigma = clt_parameters.imp.terr_um_sigma; // 1.0;
double um_weight = clt_parameters.imp.terr_um_weight; // 0.8;
double nan_tolerance = clt_parameters.imp.terr_nan_tolerance; // 0.001;
double nan_tolerance = um_en ? clt_parameters.imp.terr_nan_tolerance : 0; // 0.001;
int nan_grow = clt_parameters.imp.terr_nan_grow; // 20;
int shrink_veget = clt_parameters.imp.terr_shrink_veget; // 20;
int shrink_terrain = clt_parameters.imp.terr_shrink_terrain; // 20;
double vegetation_over = clt_parameters.imp.terr_vegetation_over; // 35;
int filter_veget = clt_parameters.imp.terr_filter_veget; // 10;
boolean tile_woi = clt_parameters.imp.terr_tile_woi; // scan woi_enclosing, false - run a single woi_last
Rectangle woi_enclosing = clt_parameters.imp.terr_woi_enclos; // new Rectangle(0, 0, 200, 160); // will be tiled, using width/height from woi_step;
......@@ -1663,6 +1685,13 @@ public class VegetationModel {
double boost_parallax = clt_parameters.imp.terr_boost_parallax; // 3.0; /// 1.0; /////// 5.0; /// 1.0; // 5;
double max_parallax = clt_parameters.imp.terr_max_parallax; // 10;
double hifreq_weight = clt_parameters.imp.terr_hifreq_weight; // 22.5; // 0 - do not use high-freq. Relative weight of laplacian components double reg_weights = 0.25; // fraction of the total weight used for regularization
boolean fit_terr = clt_parameters.imp.terr_fit_terr; // true; // adjust terrain pixels
boolean fit_veget = clt_parameters.imp.terr_fit_veget; // true; // adjust vegetation pixels
boolean fit_alpha = clt_parameters.imp.terr_fit_alpha; // true; // adjust vegetation alpha pixels
boolean fit_scenes = clt_parameters.imp.terr_fit_scenes; // true; // adjust scene offsets (start from 0 always?)
boolean fit_elevations = clt_parameters.imp.terr_fit_elevations; // false; // adjust elevation pixels (not yet implemented)
double reg_weights = clt_parameters.imp.terr_reg_weights; // 0.25; // fraction of the total weight used for regularization
double lambda = clt_parameters.imp.terr_lambda; // 5.0; // 0.1;
double lambda_scale_good = clt_parameters.imp.terr_lambda_scale_good; // 0.5;
......@@ -1722,12 +1751,8 @@ public class VegetationModel {
boolean last_run = false;
boolean test_laplacian = false;
int nan_grow_render = nan_grow;
int shrink_veget = nan_grow; // same
int filter_vegetation = nan_grow;
double vegetation_over_terrain = 35;
boolean test_laplacian = false;
if (test_laplacian) {
double [][] laplacian_in = new double [2 + terrain_scenes_render.length][];
......@@ -1767,14 +1792,16 @@ public class VegetationModel {
}
double [][] initial_terrain_vegetation = getInitialTerrainVegetation(
terrain_average_render, // double [] terrain_average_zeros,
vegetation_average_render, // double [] vegetation_average_zeros,
vegetation_warp, // final double [][][] vegetation_warp, // to count all pixels used
full.width, // int width,
shrink_veget, // int shrink_veget,
nan_grow + shrink_veget, // int shrink_terrain) {
vegetation_over_terrain, // double vegetation_over_terrain, // trust vegetation that is hotter than filled terrain
filter_vegetation); //int filter_vegetation) // shrink+grow filtered vegetation to remove small clusters
terrain_average_render, // double [] terrain_average_zeros,
vegetation_average_render, // double [] vegetation_average_zeros,
vegetation_warp, // final double [][][] vegetation_warp, // to count all pixels used
full.width, // int width,
shrink_veget, // int shrink_veget,
shrink_terrain + shrink_veget,// int shrink_terrain) {
vegetation_over, // double vegetation_over_terrain, // trust vegetation that is hotter than filled terrain
filter_veget, //int filter_vegetation) // shrink+grow filtered vegetation to remove small clusters
100); // final int extra_pull_vegetation){
if (debugLevel > -2){
ShowDoubleFloatArrays.showArrays(
initial_terrain_vegetation,
......@@ -1782,19 +1809,20 @@ public class VegetationModel {
full.height,
true,
reference_scene+"-terrain_vegetation_conditioned.tiff",
new String[] {"terrain_filled", "vegetation_fiull", "vegetation_filtered"});
new String[] {"terrain_filled", "vegetation_full", "vegetation_filtered", "pull_vegetation"});
}
terrain_filled = initial_terrain_vegetation[0];
vegetation_full = initial_terrain_vegetation[1];
vegetation_filtered = initial_terrain_vegetation[2];
vegetation_pull = initial_terrain_vegetation[3];
zerosToNans (
terrain_scenes_render, // final double [][] data,
full.width, // final int width,
nan_tolerance, // final double tolerance,
false, // final boolean negative_nan,
nan_grow_render); // final int grow)
nan_grow); // final int grow)
......@@ -1813,9 +1841,16 @@ public class VegetationModel {
um_weight); // final double um_weight)
}
double dir_sigma = 16;
for (int nscene = 0; nscene < vegetation_warp.length; nscene++) {
vegetation_warp_md[nscene ]= warpMagnitudeDirection(
vegetation_warp[nscene], // final double [][] warp_dxy,
full.width, // final int width,
dir_sigma ); // final double dir_sigma,// Gaussian blur sigma to smooth direction variations
}
if ((debugLevel > 3) || um_en) {
double [][][] dbg_img = new double[3][vegetation_warp.length][vegetation_warp[0].length];
double [][][] dbg_img = new double[6][vegetation_warp.length][vegetation_warp[0].length];
for (int n = 0; n < dbg_img.length; n++) {
for (int nscene = 0; nscene < vegetation_warp.length; nscene++) {
Arrays.fill(dbg_img[n][nscene], Double.NaN);
......@@ -1826,9 +1861,14 @@ public class VegetationModel {
if (vegetation_warp[nscene][npix] != null) {
double dx =vegetation_warp[nscene][npix][0];
double dy =vegetation_warp[nscene][npix][1];
dbg_img[0][nscene][npix] = Math.sqrt(dx*dx + dy*dy);
dbg_img[1][nscene][npix] = dx;
dbg_img[2][nscene][npix] = dy;
dbg_img[0][nscene][npix] = Math.sqrt(dx*dx + dy*dy);
dbg_img[1][nscene][npix] = Math.atan2(dy, dx);
dbg_img[2][nscene][npix] = dx;
dbg_img[3][nscene][npix] = dy;
}
if (vegetation_warp_md[nscene][npix] != null) {
dbg_img[4][nscene][npix] = vegetation_warp_md[nscene][npix][0]; // magnitude
dbg_img[5][nscene][npix] = vegetation_warp_md[nscene][npix][1]; // direction
}
}
}
......@@ -1838,7 +1878,7 @@ public class VegetationModel {
full.width, // int width,
reference_scene+"-vegetation_offsets.tiff", // String title, "time_derivs-rt"+diff_time_rt+"-rxy"+diff_time_rxy,
scene_names, // String [] titles, // all slices*frames titles or just slice titles or null
new String[] {"dist","dx","dy"}, // String [] frame_titles, // frame titles or null
new String[] {"dist","angle","dx","dy","mag","dir"}, // String [] frame_titles, // frame titles or null
true); // boolean show)
ShowDoubleFloatArrays.showArrays(
......@@ -1857,12 +1897,14 @@ public class VegetationModel {
reference_scene+"-terrain_vegetation_averages.tiff",
new String[] {"terrain","vegetation"});
}
initial_transparent = initial_split;
initial_opaque = 1.0 - initial_split;
VegetationLMA vegetationLMA = new VegetationLMA (this);
if (run_combine) {
VegetationSegment [] segments = vegetationLMA.readAllSegments(
segments_dir); // String dir_path)
......@@ -1913,6 +1955,11 @@ public class VegetationModel {
terr_pull_cold, // final double terr_pull_cold, // pull vegetations to warm, terrain to cold
default_alpha, // final double default_alpha,
hifreq_weight, //final double hifreq_weight, // 22.5 0 - do not use high-freq. Relative weight of laplacian components
fit_terr, // final boolean adjust_terr,
fit_veget, // final boolean adjust_veget,
fit_alpha, // final boolean adjust_alpha,
fit_scenes, // final boolean adjust_scenes,
fit_elevations, // final boolean adjust_elevations,
reg_weights, // final double reg_weights, // fraction of the total weight used for regularization
alpha_loss, // final double alpha_loss, // quadratic loss when alpha reaches -1.0 or 2.0
alpha_offset, // final double alpha_offset, // quadratic loss when alpha reaches -1.0 or 2.0
......@@ -1989,7 +2036,7 @@ public class VegetationModel {
}
/// next_run = true;
vegetationLMA.debug_index = 0;
vegetationLMA.debug_image = new double [num_iter][];
vegetationLMA.debug_image = new double [100][]; // num_iter][];
int lma_rslt= vegetationLMA.runLma( // <0 - failed, >=0 iteration number (1 - immediately)
lambda, // double lambda, // 0.1
......@@ -2023,6 +2070,109 @@ public class VegetationModel {
}
public static double [][] warpMagnitudeDirection(
final double [][] warp_dxy,
final int width,
final double dir_sigma){// Gaussian blur sigma to smooth direction variations
// Assuming objects have positive elevation, so dx, dy sign does not change
final int num_pixels = warp_dxy.length;
final double [][] warp_md = new double [warp_dxy.length][];
final Thread[] threads = ImageDtt.newThreadArray();
final AtomicInteger ai = new AtomicInteger(0);
final AtomicInteger ati = new AtomicInteger(0);
final double [][] sdxy = new double [threads.length][2];
for (int ithread = 0; ithread < threads.length; ithread++) { // first sum for pairs
threads[ithread] = new Thread() {
public void run() {
int nthread = ati.getAndIncrement();
double sdx = 0, sdy = 0;
for (int nPix = ai.getAndIncrement(); nPix < num_pixels; nPix = ai.getAndIncrement()) {
double [] dxy = warp_dxy[nPix];
if (dxy != null) {
sdx += dxy[0];
sdy += dxy[1];
}
}
sdxy[nthread][0] = sdx;
sdxy[nthread][1] = sdy;
}
};
}
ImageDtt.startAndJoin(threads);
double sdx = 0, sdy = 0;
for (int nthread =0; nthread < sdxy.length; nthread++) {
sdx += sdxy[nthread][0];
sdy += sdxy[nthread][1];
}
final double avg_dir = Math.atan2(sdy, sdx); // average angle to avoid crossing full rotation
final double [] mags = new double[num_pixels], dirs = new double[num_pixels], mag_dirs=new double[num_pixels];
// Arrays.fill(dirs, Double.NaN);
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) { // first sum for pairs
threads[ithread] = new Thread() {
public void run() {
for (int nPix = ai.getAndIncrement(); nPix < num_pixels; nPix = ai.getAndIncrement()) {
double [] dxy = warp_dxy[nPix];
if (dxy != null) {
double m = Math.sqrt(dxy[0]*dxy[0]+dxy[1]*dxy[1]);
mags[nPix] = m;
double angle = Math.atan2(dxy[1],dxy[0]) - avg_dir;
if (angle > Math.PI) {
angle -= 2 * Math.PI;
} else if (angle < -Math.PI){
angle += 2 * Math.PI;
}
dirs[nPix] = angle;
mag_dirs[nPix] = m * angle;
}
}
}
};
}
ImageDtt.startAndJoin(threads);
final double [] mag_blur = mags.clone();
(new DoubleGaussianBlur()).blurDouble(
mag_blur, //
width,
num_pixels/width,
dir_sigma, // double sigmaX,
dir_sigma, // double sigmaY,
0.01); // double accuracy)
(new DoubleGaussianBlur()).blurDouble(
mag_dirs, //
width,
num_pixels/width,
dir_sigma, // double sigmaX,
dir_sigma, // double sigmaY,
0.01); // double accuracy)
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) { // first sum for pairs
threads[ithread] = new Thread() {
public void run() {
for (int nPix = ai.getAndIncrement(); nPix < num_pixels; nPix = ai.getAndIncrement()) {
double [] dxy = warp_dxy[nPix];
if (dxy != null) {
double angle = mag_dirs[nPix]/mag_blur[nPix] + avg_dir;
if (angle > Math.PI) {
angle -= 2 * Math.PI;
} else if (angle < -Math.PI){
angle += 2 * Math.PI;
}
if (Double.isNaN(angle)) {
angle = 0;
}
warp_md[nPix] = new double [] {mags[nPix], angle};
}
}
}
};
}
ImageDtt.startAndJoin(threads);
return warp_md;
}
......@@ -3074,7 +3224,8 @@ public class VegetationModel {
final int shrink_veget,
final int shrink_terrain,
final double vegetation_over_terrain, // trust vegetation that is hotter than filled terrain
final int filter_vegetation) { // shrink+grow filtered vegetation to remove small clusters
final int filter_vegetation, // shrink+grow filtered vegetation to remove small clusters
final int extra_pull_vegetation){
boolean debug_img = false; // true;
double [] terrain_filled = terrain_average_zeros.clone();
double [] vegetation_holes = vegetation_average_zeros.clone();
......@@ -3168,7 +3319,22 @@ public class VegetationModel {
for (int i = 0; i < vegetation_mask.length; i++) if (!vegetation_mask[i]){
vegetation_filtered[i] = Double.NaN;
}
double [][] result = {terrain_filled, vegetation_full, vegetation_filtered} ;
// Create pull vegetation
// boolean
//extra_pull_vegetation
double [] pull_vegetation = TileProcessor.fillNaNs(
vegetation_filtered, // final double [] data,
null, // final boolean [] prohibit,
width, // int width,
// CAREFUL ! Remaining NaN is grown by unsharp mask filter ************* !
extra_pull_vegetation, // 100, // 2*width, // 16, // final int grow,
0.7, // double diagonal_weight, // relative to ortho
100, // int num_passes,
0.03, // final double max_rchange, // = 0.01 - does not need to be accurate
ImageDtt.THREADS_MAX); // final int threadsMax) // maximal number of threads to launch
double [][] result = {terrain_filled, vegetation_full, vegetation_filtered, pull_vegetation} ;
return result;
}
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment