Skip to content

I
imagej-elphel
Commit 6a053711 authored by Andrey Filippov's avatar Andrey Filippov
Browse files
Options

before vegetation_lpf modification

parent cc900bac
Show whitespace changes
Inline Side-by-side
Showing with 965 additions and 604 deletions
src/main/java/com/elphel/imagej/tileprocessor/IntersceneMatchParameters.java
View file @ 6a053711
......@@ -803,7 +803,8 @@ min_str_neib_fpn 0.35
public double terr_alpha_loss = 100.0;
public double terr_alpha_loss_lin = 0.0;
public double terr_alpha_offset = 0.1;
public double terr_alpha_offset = 0.1; // alpha offset near alpha ~=1.0
public double terr_alpha_0offset = 0.0; // New - a separate alpha offset nea alpha ~= 0.0
public double terr_alpha_min_veg = -1; // old version 0.5; // Minimal vegetation alpha. If (alpha-alpha_offset)/(1-2*alpha_offset) < alpha_min_veg, pull down to lpha_offset
......@@ -829,6 +830,22 @@ min_str_neib_fpn 0.35
public double terr_elev_alpha_pwr = 2.0; // raise alpha to this power (when alpha > 0)
public double terr_low_veget = 2.0; // (pix) Elevation considered low (lower loss for high alpha)
public double terr_scenes_pull0 = 1.0; // pull average scene offset to zero
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
public double terr_elev_scale_thresh = 1.0; // reduce losses for higher (initial) elevations TODO: consider actual elevations
public boolean terr_elev_scale_pull = false; // scale elevation pull losses for high elevations
public boolean terr_elev_scale_lpf = false; // scale elevation diffusion losses for high elevations
// tree-top removal
public boolean terr_ttop_en = false; // remove tree tops from transparency weights
public double terr_ttop_gb = 1.0; // Elevation Gaussian blur sigma to detect tree tops
public double terr_ttop_min = 3.0; // Minimal tree top elevation
public double terr_ttop_rel_lev = 0.9; // Relative (to the top height) sample level
public double terr_ttop_rel_rad = 0.25; // Relative (to the top height) sample ring radius
public double terr_ttop_frac = 0.5; // Minimal fraction of the ring pixels below sample level
public double terr_ttop_rem_rad = 0.25; // Relative (to the top height) remove transparency radius
// LMA parameters
public double terr_boost_parallax = 3.0; //
......@@ -2134,7 +2151,8 @@ min_str_neib_fpn 0.35
gd.addNumericField("Alpha loss", terr_alpha_loss, 5,7,"", "Alpha quadratic growing loss for when out of [0,1] range");
gd.addNumericField("Alpha loss linear", terr_alpha_loss_lin, 5,7,"", "Alpha linear growing loss for when out of [0,1] range and below minimal vegetation alpha.");
gd.addNumericField("Alpha offset", terr_alpha_offset, 5,7,"", "Start alpha losses above 0.0 and below 1.0 by this value.");
gd.addNumericField("Alpha offset @1.0", terr_alpha_offset, 5,7,"", "Start alpha losses below 1.0 by this value.");
gd.addNumericField("Alpha offset @0.0", terr_alpha_0offset, 5,7,"", "Start alpha losses above 0.0 by this value.");
gd.addNumericField("Minimal vegetation alpha", terr_alpha_min_veg, 5,7,"","Minimal vegetation alpha. If (alpha-alpha_offset)/(1-2*alpha_offset) < alpha_min_veg, pull down to lpha_offset.");
gd.addNumericField("Alpha diffusion", terr_alpha_lpf, 5,7,"", "Alpha diffusion to 4 ortho neighbors.");
......@@ -2159,9 +2177,23 @@ min_str_neib_fpn 0.35
gd.addNumericField("Low vegetation loss", terr_elev_alpha, 5,7,"", "Multiply alpha by under-low elevation for loss.");
gd.addNumericField("Low vegetation power", terr_elev_alpha_pwr, 5,7,"", "Raise alpha to this power for low vegetation loss.");
gd.addNumericField("Low elevation", terr_low_veget, 5,7,"pix", "Elevation considered low (lower loss for high alpha.");
gd.addNumericField("Pull scene offset", terr_scenes_pull0, 5,7,"", "Pull average scene offset to zero.");
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
gd.addNumericField("High elevation", terr_elev_scale_thresh, 5,7,"","Reduce losses for higher (initial) elevations threshold.");
gd.addCheckbox ("High elevation pull", terr_elev_scale_pull, "Scale elevation pull losses for high elevations.");
gd.addCheckbox ("High elevation diffusion",terr_elev_scale_lpf, "Scale elevation diffusion losses for high elevations.");
gd.addMessage ("Tree tops detection/filtering - no (stray) transparency around tree tops");
gd.addCheckbox ("Filter tree tops", terr_ttop_en, "Filter tree tops transparency.");
gd.addNumericField("Elevation blur sigma", terr_ttop_gb, 5,7,"", "Elevation Gaussian blur sigma to detect tree tops.");
gd.addNumericField("Minimal tree top elevation",terr_ttop_min, 5,7,"", "Minimal tree top elevation.");
gd.addNumericField("Relative sample level",terr_ttop_rel_lev, 5,7,"", "Relative (to the top height) sample level.");
gd.addNumericField("Relative ring radius", terr_ttop_rel_rad, 5,7,"", "Relative (to the top height) sample ring radius.");
gd.addNumericField("Minimal fraction lower",terr_ttop_frac, 5,7,"", "Minimal fraction of the ring pixels below sample level.");
gd.addNumericField("Remove transparency radius", terr_ttop_rem_rad, 5,7,"","Relative (to the top height) remove transparency radius.");
gd.addMessage ("LMA parameters");
gd.addNumericField("Boost parallax", terr_boost_parallax, 5,7,"", "Increase weight of scenes that have high parallax to the reference one.");
gd.addNumericField("Limit parallax", terr_max_parallax, 5,7,"", "Parallax limit when evaluating boost parallax.");
......@@ -2899,6 +2931,7 @@ min_str_neib_fpn 0.35
terr_alpha_loss = gd.getNextNumber();// double
terr_alpha_loss_lin = gd.getNextNumber();// double
terr_alpha_offset = gd.getNextNumber();// double
terr_alpha_0offset = gd.getNextNumber();// double
terr_alpha_min_veg = gd.getNextNumber();// double
terr_alpha_lpf = gd.getNextNumber();// double
terr_alpha_piece_linear = gd.getNextBoolean();// boolean
......@@ -2923,6 +2956,20 @@ min_str_neib_fpn 0.35
terr_elev_alpha_pwr = gd.getNextNumber();// double
terr_low_veget = gd.getNextNumber();// double
terr_scenes_pull0 = gd.getNextNumber();// double
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
terr_elev_scale_thresh = gd.getNextNumber();// double
terr_elev_scale_pull = gd.getNextBoolean();// boolean
terr_elev_scale_lpf = gd.getNextBoolean();// boolean
// tree-top removal
terr_ttop_en = gd.getNextBoolean();// boolean
terr_ttop_gb = gd.getNextNumber();// double
terr_ttop_min = gd.getNextNumber();// double
terr_ttop_rel_lev = gd.getNextNumber();// double
terr_ttop_rel_rad = gd.getNextNumber();// double
terr_ttop_frac = gd.getNextNumber();// double
terr_ttop_rem_rad = gd.getNextNumber();// double
terr_boost_parallax = gd.getNextNumber();// double
terr_max_parallax = gd.getNextNumber();// double
terr_hifreq_weight = gd.getNextNumber();// double
......@@ -3625,6 +3672,7 @@ min_str_neib_fpn 0.35
properties.setProperty(prefix+"terr_alpha_loss", terr_alpha_loss+""); // double
properties.setProperty(prefix+"terr_alpha_loss_lin", terr_alpha_loss_lin+""); // double
properties.setProperty(prefix+"terr_alpha_offset", terr_alpha_offset+""); // double
properties.setProperty(prefix+"terr_alpha_0offset", terr_alpha_0offset+""); // double
properties.setProperty(prefix+"terr_alpha_min_veg", terr_alpha_min_veg+""); // double
properties.setProperty(prefix+"terr_alpha_lpf", terr_alpha_lpf+""); // double
properties.setProperty(prefix+"terr_alpha_piece_linear", terr_alpha_piece_linear+""); // boolean
......@@ -3650,6 +3698,20 @@ min_str_neib_fpn 0.35
properties.setProperty(prefix+"terr_low_veget", terr_low_veget+""); // double
properties.setProperty(prefix+"terr_scenes_pull0", terr_scenes_pull0+""); // double
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
properties.setProperty(prefix+"terr_elev_scale_thresh", terr_elev_scale_thresh+""); // double
properties.setProperty(prefix+"terr_elev_scale_pull", terr_elev_scale_pull+""); // boolean
properties.setProperty(prefix+"terr_elev_scale_lpf", terr_elev_scale_lpf+""); // boolean
// tree-top removal
properties.setProperty(prefix+"terr_ttop_en", terr_ttop_en+""); // boolean
properties.setProperty(prefix+"terr_ttop_gb", terr_ttop_gb+""); // double
properties.setProperty(prefix+"terr_ttop_min", terr_ttop_min+""); // double
properties.setProperty(prefix+"terr_ttop_rel_lev", terr_ttop_rel_lev+""); // double
properties.setProperty(prefix+"terr_ttop_rel_rad", terr_ttop_rel_rad+""); // double
properties.setProperty(prefix+"terr_ttop_frac", terr_ttop_frac+""); // double
properties.setProperty(prefix+"terr_ttop_rem_rad", terr_ttop_rem_rad+""); // double
properties.setProperty(prefix+"terr_boost_parallax", terr_boost_parallax+""); // double
properties.setProperty(prefix+"terr_max_parallax", terr_max_parallax+""); // double
properties.setProperty(prefix+"terr_hifreq_weight", terr_hifreq_weight+""); // double
......@@ -4371,6 +4433,7 @@ min_str_neib_fpn 0.35
if (properties.getProperty(prefix+"terr_alpha_loss")!= null) terr_alpha_loss=Double.parseDouble(properties.getProperty(prefix+"terr_alpha_loss"));
if (properties.getProperty(prefix+"terr_alpha_loss_lin")!= null) terr_alpha_loss_lin=Double.parseDouble(properties.getProperty(prefix+"terr_alpha_loss_lin"));
if (properties.getProperty(prefix+"terr_alpha_offset")!= null) terr_alpha_offset=Double.parseDouble(properties.getProperty(prefix+"terr_alpha_offset"));
if (properties.getProperty(prefix+"terr_alpha_0offset")!= null) terr_alpha_0offset=Double.parseDouble(properties.getProperty(prefix+"terr_alpha_0offset"));
if (properties.getProperty(prefix+"terr_alpha_min_veg")!= null) terr_alpha_min_veg=Double.parseDouble(properties.getProperty(prefix+"terr_alpha_min_veg"));
if (properties.getProperty(prefix+"terr_alpha_lpf")!= null) terr_alpha_lpf=Double.parseDouble(properties.getProperty(prefix+"terr_alpha_lpf"));
if (properties.getProperty(prefix+"terr_alpha_piece_linear")!= null) terr_alpha_piece_linear=Boolean.parseBoolean(properties.getProperty(prefix+"terr_alpha_piece_linear"));
......@@ -4395,6 +4458,18 @@ min_str_neib_fpn 0.35
if (properties.getProperty(prefix+"terr_elev_alpha_pwr")!= null) terr_elev_alpha_pwr=Double.parseDouble(properties.getProperty(prefix+"terr_elev_alpha_pwr"));
if (properties.getProperty(prefix+"terr_low_veget")!= null) terr_low_veget=Double.parseDouble(properties.getProperty(prefix+"terr_low_veget"));
if (properties.getProperty(prefix+"terr_scenes_pull0")!= null) terr_scenes_pull0=Double.parseDouble(properties.getProperty(prefix+"terr_scenes_pull0"));
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
if (properties.getProperty(prefix+"terr_elev_scale_thresh")!= null) terr_elev_scale_thresh=Double.parseDouble(properties.getProperty(prefix+"terr_elev_scale_thresh"));
if (properties.getProperty(prefix+"terr_elev_scale_pull")!= null) terr_elev_scale_pull=Boolean.parseBoolean(properties.getProperty(prefix+"terr_elev_scale_pull"));
if (properties.getProperty(prefix+"terr_elev_scale_lpf")!= null) terr_elev_scale_lpf=Boolean.parseBoolean(properties.getProperty(prefix+"terr_elev_scale_lpf"));
// tree-top removal
if (properties.getProperty(prefix+"terr_ttop_en")!= null) terr_ttop_en=Boolean.parseBoolean(properties.getProperty(prefix+"terr_ttop_en"));
if (properties.getProperty(prefix+"terr_ttop_gb")!= null) terr_ttop_gb=Double.parseDouble(properties.getProperty(prefix+"terr_ttop_gb"));
if (properties.getProperty(prefix+"terr_ttop_min")!= null) terr_ttop_min=Double.parseDouble(properties.getProperty(prefix+"terr_ttop_min"));
if (properties.getProperty(prefix+"terr_ttop_rel_lev")!= null) terr_ttop_rel_lev=Double.parseDouble(properties.getProperty(prefix+"terr_ttop_rel_lev"));
if (properties.getProperty(prefix+"terr_ttop_rel_rad")!= null) terr_ttop_rel_rad=Double.parseDouble(properties.getProperty(prefix+"terr_ttop_rel_rad"));
if (properties.getProperty(prefix+"terr_ttop_frac")!= null) terr_ttop_frac=Double.parseDouble(properties.getProperty(prefix+"terr_ttop_frac"));
if (properties.getProperty(prefix+"terr_ttop_rem_rad")!= null) terr_ttop_rem_rad=Double.parseDouble(properties.getProperty(prefix+"terr_ttop_rem_rad"));
if (properties.getProperty(prefix+"terr_boost_parallax")!= null) terr_boost_parallax=Double.parseDouble(properties.getProperty(prefix+"terr_boost_parallax"));
if (properties.getProperty(prefix+"terr_max_parallax")!= null) terr_max_parallax=Double.parseDouble(properties.getProperty(prefix+"terr_max_parallax"));
......@@ -5087,6 +5162,7 @@ min_str_neib_fpn 0.35
imp.terr_alpha_loss = this.terr_alpha_loss;
imp.terr_alpha_loss_lin = this.terr_alpha_loss_lin;
imp.terr_alpha_offset = this.terr_alpha_offset;
imp.terr_alpha_0offset = this.terr_alpha_0offset;
imp.terr_alpha_min_veg = this.terr_alpha_min_veg;
imp.terr_alpha_lpf = this.terr_alpha_lpf;
imp.terr_alpha_piece_linear = this.terr_alpha_piece_linear;
......@@ -5110,6 +5186,18 @@ min_str_neib_fpn 0.35
imp.terr_elev_alpha_pwr = this.terr_elev_alpha_pwr;
imp.terr_low_veget = this.terr_low_veget;
imp.terr_scenes_pull0 = this.terr_scenes_pull0;
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
imp.terr_elev_scale_thresh = this.terr_elev_scale_thresh;
imp.terr_elev_scale_pull = this.terr_elev_scale_pull;
imp.terr_elev_scale_lpf = this.terr_elev_scale_lpf;
// tree-top removal
imp.terr_ttop_en = this.terr_ttop_en;
imp.terr_ttop_gb = this.terr_ttop_gb;
imp.terr_ttop_min = this.terr_ttop_min;
imp.terr_ttop_rel_lev = this.terr_ttop_rel_lev;
imp.terr_ttop_rel_rad = this.terr_ttop_rel_rad;
imp.terr_ttop_frac = this.terr_ttop_frac;
imp.terr_ttop_rem_rad = this.terr_ttop_rem_rad;
imp.terr_boost_parallax = this.terr_boost_parallax;
imp.terr_max_parallax = this.terr_max_parallax;
......
src/main/java/com/elphel/imagej/vegetation/VegetationLMA.java
View file @ 6a053711
......@@ -87,19 +87,6 @@ public class VegetationLMA {
private int [] num_pars = new int [TVAO_TYPES];
private int [] ind_pars = new int [TVAO_TYPES]; // should be in increasing order?
/*
private int num_pars[TVAO_TERRAIN];
private int num_pars[TVAO_VEGETATION];
private int num_pars[TVAO_ALPHA]; //== num_pars[TVAO_VEGETATION]
private int num_pars[TVAO_ELEVATION]; //== num_pars[TVAO_VEGETATION]
private int num_pars[TVAO_SCENE_OFFSET];
private int ind_pars[TVAO_TERRAIN] = 0; // index of the first terrain parameter
private int ind_pars[TVAO_VEGETATION]; // index of the first vegetation parameter
private int ind_pars[TVAO_ALPHA]; // index of the first alpha parameter
private int ind_pars[TVAO_ELEVATION]; // index of the first elevation parameter
private int ind_pars[TVAO_SCENE_OFFSET]; // index of the first scene parameter
*/
private int [][] y_src; // data pointers for y-vector
private int [][] y_src_hf; // subset of pointers that have all 4 neighbors
private int [][] y_src_scene; // [scene]{start_y_src, end_y_src+1} pointers to y_src per scene
......@@ -159,7 +146,8 @@ public class VegetationLMA {
private boolean [] fits = new boolean[TVAO_TYPES];// (+)
private double alpha_loss = 0; // (+) quadratic loss
private double alpha_loss_lin = 0.5; // (+) linear loss
private double alpha_offset = 0; // (+) if >0, start losses above 0.0 and below 1.0;
private double alpha_offset = 0; // (+) if >0, start losses below 1.0;
private double alpha_0offset = 0; // (+) if >0, start losses above 0.0
private double alpha_min_veg = 0.5; // (+) if (alpha-alpha_offset)/(1-2*alpha_offset) < alpha_min_veg, pull down to lpha_offset
private double alpha_lpf = 0; // (+)
......@@ -177,10 +165,6 @@ public class VegetationLMA {
private double veget_lpf = 0; // (+)
private double elevation_lpf = 0; // (*)
private double elevation_lpf_transp = 1000; // lpf for transparent areas
private double elevation_lpf_pow = 2; // lpf for transparent areas power
private boolean [] lpf_fixed = {true, true,true, true, true}; // (*) using Laplacian to lpf parameters, pull to fixed ones too
// when unsharp mask is applied , pulling to 0 (when alpha is 0 (for vegetation) or 1.0 (for terrain) makes sense
......@@ -196,13 +180,37 @@ public class VegetationLMA {
public double low_veget = 2.0; // (+) (pix) Elevation considered low (lower loss for high alpha)
private double scenes_pull0 = 0; // (*) pull average scene offset to 0;
private double scale_scenes_pull = 0; // (*) used in getFxDerivs to scale scene offsets as their weight will be reg_weights / extra_samples
private double scale_scenes_pull = 0; // (?) used in getFxDerivs to scale scene offsets as their weight will be reg_weights / extra_samples
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
public double elev_scale_thresh = 1.0; // (+) reduce losses for higher (initial) elevations TODO: consider actual elevations
public boolean elev_scale_pull = false; // (+) scale elevation pull losses for high elevations
public boolean elev_scale_lpf = false; // (+) scale elevation diffusion losses for high elevations
// tree-top removal
public boolean ttop_en = false; // (+) remove tree tops from transparency weights
public double ttop_gb = 1.0; // (+) Elevation Gaussian blur sigma to detect tree tops
public double ttop_min = 3.0; // (+) Minimal tree top elevation
public double ttop_rel_lev = 0.9; // (+) Relative (to the top height) sample level
public double ttop_rel_rad = 0.25; // (+) Relative (to the top height) sample ring radius
public double ttop_frac = 0.5; // (+) Minimal fraction of the ring pixels below sample level
public double ttop_rem_rad = 0.25; // (+) Relative (to the top height) remove transparency radius
public double boost_parallax = 1; // (+)
public double max_parallax = 10; // (+)
public double max_warp = 3.8; // 1.8 - do not use scenes where distance between vegetation projection exceeds this
public int max_elevation = 22; // maximal "elevation" to consider
public double elevation_radius = 1.5; // Radius of elevation/vegetation influence.
private double elevation_lpf_transp = 1000; // lpf for transparent areas
private double elevation_lpf_pow = 2; // lpf for transparent areas power
// data used to calculate lpf pull of the alpha pixel to average of four neighbors. Below similar (weaker pull) for terrain and vegetation
// to smooth areas where there is no data from available images.
public int [][] alpha_neibs; // corresponds to parameters for alpha (num_pars[TVAO_VEGETATION]_alpha), each has 4 ortho neibs, -1 - border, >= 0
......@@ -250,33 +258,6 @@ public class VegetationLMA {
return vegetationModel;
}
public VegetationLMA (
int width,
boolean diff_mode,
double [] terrain_average, // full image average terrain (no nulls!)
double [] vegetation_average, // full image average vegetation (with nulls)
double [][] terrain_rendered, // terrain rendered for scenes (has nulls)
double [][][] vegetation_offsets // [num_scenes][pixel]{dx,dy} differential offsets of vegetation to terrain grid
) {
full = new Rectangle(0,0,width, terrain_average.length/width);
image_length = terrain_average.length;
// 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[TVAO_TYPES][];
tvao[TVAO_TERRAIN] = this.terrain_average.clone();
tvao[TVAO_VEGETATION] = this.vegetation_average.clone();
tvao[TVAO_ALPHA] = new double [image_length]; // 0 - use terrain
tvao[TVAO_SCENE_OFFSET] = new double [num_scenes];
setupLaplacians(0.0); // double weight_diag)
scene_weights = new double [num_scenes];
}
public VegetationLMA (
VegetationModel vegetationModel,
......@@ -369,43 +350,6 @@ public class VegetationLMA {
weight_diag); // final double weight_diag)
}
}
private void setupInitialTerrainVegetationAlpha(
double alpha_contrast,
double [] terrain_filled,
double [] vegetation_extended,
double [] vegetation, // best guess
int debugLevel) {
double dbg_lim = 0.01;
tvao[TVAO_TERRAIN] = terrain_filled.clone(); // == vegetationModel.terrain_filled
tvao[TVAO_VEGETATION] = vegetation_extended.clone();
tvao[TVAO_ALPHA] = new double [image_length];
for (int npix = 0; npix < image_length; npix++) {
double t = terrain_filled[npix];
double v = vegetation[npix];
double vp = vegetation_extended[npix];
double a = 0; // Double.NaN; // .isNaN(vegetation_filtered[npix])? 0.0:1.0; // not needed
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 (a < dbg_lim) a = dbg_lim;
// if ( a > 1.0-dbg_lim) a = 1.0-dbg_lim;
// tvao[TVAO_ALPHA][npix] = a;
}
if (a < dbg_lim) a = dbg_lim;
if ( a > 1.0-dbg_lim) a = 1.0-dbg_lim;
tvao[TVAO_ALPHA][npix] = a;
}
}
private void setupInitialTerrainVegetationAlpha(
boolean blur_vegetation,
......@@ -475,34 +419,6 @@ public class VegetationLMA {
}
@SuppressWarnings("unused")
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(
final boolean keep_parameters,
final Rectangle woi,
......@@ -522,7 +438,8 @@ public class VegetationLMA {
final double reg_weights, // fraction of the total weight used for regularization
final double alpha_loss, // alpha quadratic growing loss for when out of [0,1] range
final double alpha_loss_lin, // alpha linear growing loss for when out of [0,1] range and below minimal vegetation alpha
final double alpha_offset, // quadratic loss when alpha reaches -1.0 or 2.0
final double alpha_offset, // quadratic loss when alpha > 1 - alpha_offset
final double alpha_0offset, // quadratic loss when alpha < alpha0_offset
final double alpha_min_veg, // 0.5; // if (alpha-alpha_offset)/(1-2*alpha_offset) < alpha_min_veg, pull down to lpha_offset
final double alpha_lpf, // pull vegetation alpha to average of 4 neighbors
final boolean alpha_piece_linear, // true - piece-linear, false - half-cosine
......@@ -546,6 +463,16 @@ public class VegetationLMA {
final double elev_alpha_pwr, // 2.0; // raise alpha to this power (when alpha > 0)
final double low_veget, // 2.0; // (pix) Elevation considered low (lower loss for high alpha)
final double scenes_pull0,
final double elev_scale_thresh, // 1.0;// reduce losses for higher (initial) elevations TODO: consider actual elevations
final boolean elev_scale_pull,// false; // scale elevation pull losses for high elevations
final boolean elev_scale_lpf, // false; // scale elevation diffusion losses for high elevations
final boolean ttop_en, // false; // remove tree tops from transparency weights
final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
final double ttop_min, // 3.0; // Minimal tree top elevation
final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
final double boost_parallax, // increase weight of scene with maximal parallax relative to the reference scene
final double max_parallax, // do not consider maximal parallax above this (consider it a glitch)
// final double um_sigma, // just use in debug image names
......@@ -569,6 +496,7 @@ public class VegetationLMA {
this.alpha_loss = alpha_loss;
this.alpha_loss_lin = alpha_loss_lin;
this.alpha_offset = alpha_offset;
this.alpha_0offset = alpha_0offset;
this.alpha_min_veg = Math.max(alpha_min_veg,0);
this.alpha_lpf = alpha_lpf;
this.alpha_scale_avg = alpha_scale_avg;
......@@ -591,6 +519,18 @@ public class VegetationLMA {
this.elev_alpha_pwr = elev_alpha_pwr;
this.low_veget = low_veget;
this.scenes_pull0 = scenes_pull0;
this.elev_scale_thresh = elev_scale_thresh;
this.elev_scale_pull = elev_scale_pull;
this.elev_scale_lpf = elev_scale_lpf;
this.ttop_en = ttop_en;
this.ttop_gb = ttop_gb;
this.ttop_min = ttop_min;
this.ttop_rel_lev = ttop_rel_lev;
this.ttop_rel_rad = ttop_rel_rad;
this.ttop_frac = ttop_frac;
this.ttop_rem_rad = ttop_rem_rad;
this.boost_parallax = boost_parallax;
this.max_parallax = max_parallax; // parallax limit when evaluating boost parallax
this.max_warp = max_warp;
......@@ -759,22 +699,7 @@ public class VegetationLMA {
"woi_veg_uses_"+woi_veg.x+"-"+woi_veg.y+"-"+woi_veg.width+"-"+woi_veg.height,
titles);
}
// used_veg = valid_scene_pix[2];
/*
if ("RESTORE".equals(parameters_read_path)) {
from_file = false;
String restore_dir = debug_path;
String debug_title = getParametersDebugTitle();
if (!restore_dir.endsWith(Prefs.getFileSeparator())) {
restore_dir += Prefs.getFileSeparator();
}
String restore_path=restore_dir + debug_title+"-live.tiff";
parameters_read_path = restore_path;
}
*/
// final double [] scene_weights = //sets this.scene_weights too
double [] scene_weights0 = //sets this.scene_weights too
setupSceneWeights(
boost_parallax, // double boost_parallax)
......@@ -839,32 +764,16 @@ public class VegetationLMA {
}
from_file = false;
/*
if (parameters_read_path != null) {
readParametersFromImage(
parameters_read_path, // String path,
parameters_vector, // double [] vector,
show_extra, // boolean extra,
1) ;// int gap)
}
*/
boolean use_terr_corr = (terrain_correction>=0); // maybe not needed
setupYVector(
use_terr_corr, // boolean use_terr_corr); //(hifreq_weight > 0)); // boolean use_hf);
scene_weights); //
/*
if (use_y_avg) {
setupYWavg(scene_weights); // double [] scene_weights);
}
*/
setupWeights( // after setupParametersIndices
scene_weights,
reg_weights, // final double reg_weights,
hifreq_weight, // final double hifreq_weight );
debugLevel); // final int debugLevel) {
// false, // final boolean apply_transparency,
// 1.0, // final double alpha_opaque,
// 0.0); // final double alpha_pedestal) {
last_jt = new double [parameters_vector.length][];
return weights.length;
......@@ -1471,7 +1380,7 @@ public class VegetationLMA {
String debug_title = "parameters_vector-x"+woi.x+"-y"+woi.y+"-w"+woi.width+"-h"+woi.height;
if (!short_name) {
debug_title += "-hf"+hifreq_weight+"-tc"+terrain_correction;
debug_title += "-al"+alpha_loss+"-alo"+alpha_offset+"-amv"+alpha_min_veg+"-alp"+alpha_lpf+"-als"+alpha_scale_avg+(alpha_piece_linear?"-alin":"-acos");
debug_title += "-al"+alpha_loss+"-alo"+alpha_offset+"-al0o"+alpha_0offset+"-amv"+alpha_min_veg+"-alp"+alpha_lpf+"-als"+alpha_scale_avg+(alpha_piece_linear?"-alin":"-acos");
debug_title += "-ap"+alpha_push+"-apn"+alpha_push_neutral+"-apc"+alpha_push_center+(alpha_en_holes?("-amm"+alpha_mm_hole):"");
debug_title += "-tl"+terr_lpf+"-vl"+veget_lpf+"-el"+elevation_lpf+"-tp"+terr_pull0+"-tu"+terr_pull_up+"-ta"+terr_pull_avg+"-vp"+veget_pull0+"-ep"+elevation_pull0;
if (elev_alpha_en) {
......@@ -1850,78 +1759,6 @@ public class VegetationLMA {
}
/*
@Deprecated
private double [][] getTerrainWeights( // FIXME !
final double alpha_threshold, // discard images with too low transparency
final double [] vector)
{
final int woi_length = woi_veg.width*woi_veg.height;
final double [][] terrain_scene_weights = new double [num_scenes][woi_length];
final double [][] terrain_scene_max = new double [num_scenes][woi_length];
final double [][] full_scene_weights = new double [num_scenes][woi_length];
final double [][] terrain_weights = new double [2][woi_length];
final Thread[] threads = ImageDtt.newThreadArray();
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
double [] alpha = new double [4];
for (int n = ai.getAndIncrement(); n < y_vector.length; n = ai.getAndIncrement()) {
int nscene = y_src[n][YSRC_SCENE];
int indx = y_src[n][YSRC_FINDEX];
int x = indx % full.width - woi_veg.x;
int y = indx / full.width - woi_veg.y;
int windx = x + y * woi_veg.width;
double [] cw = corners_weights[n];
int [] indx_alpha = data_source[n][2];
double sum_a =0;
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_ALPHA][-1-ia]);
sum_a += alpha[i];
}
double tw = Math.max(0, 1.0 - sum_a);
terrain_scene_max[nscene][windx] = tw;
if (tw >= alpha_threshold) {
terrain_scene_weights[nscene][windx] = weights[n] * tw;
}
full_scene_weights[nscene][windx] = weights[n];
}
}
}
};
}
ImageDtt.startAndJoin(threads);
// sum for all scenes, divide (undefined - also 0)
Arrays.fill(terrain_weights[0], Double.NaN);
Arrays.fill(terrain_weights[1], Double.NaN);
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int windx = ai.getAndIncrement(); windx < woi_length; windx = ai.getAndIncrement()) {
double sum_w = 0, sum_wt = 0, max_t = 0;
for (int nscene = 0; nscene < num_scenes;nscene++) {
sum_wt += terrain_scene_weights[nscene][windx];
sum_w += full_scene_weights[nscene][windx];
max_t = Math.max(max_t, terrain_scene_max[nscene][windx]);
}
if (sum_w > 0) {
terrain_weights[0][windx] = sum_wt/sum_w;
terrain_weights[1][windx] = max_t;
}
}
}
};
}
ImageDtt.startAndJoin(threads);
// regularization weights and derivatives
return terrain_weights;
}
*/
private void setupElevationLMA( // return pivot?
final double [] vector,
......@@ -1985,35 +1822,22 @@ public class VegetationLMA {
int wvy = wvindex / woi_veg.width; // relative to woi_veg
int x = wvx + woi_veg.x; // relative to full
int y = wvy + woi_veg.y; // relative to full
int wx = x - woi.x; // relative to woi
int wy = x - woi.y; // relative to woi
int npix = x + y * full.width;
// double elevation = tvao[TVAO_ELEVATION][npix];
int fnpix = x + y * full.width;
double elevation;
if (elevation_parameters) {
int npar = par_index[TVAO_ELEVATION][npix];
int npar = par_index[TVAO_ELEVATION][fnpix];
if (npar >= 0) {
elevation = vector[npar];
} else {
System.out.println("setupElevationLMA(): BUG: par_index["+TVAO_ELEVATION+"]["+npix+"] ="+npar);
System.out.println("setupElevationLMA(): BUG: par_index["+TVAO_ELEVATION+"]["+fnpix+"] ="+npar);
continue;
}
} else {
elevation = tvao[TVAO_ELEVATION][npix];
elevation = tvao[TVAO_ELEVATION][fnpix];
}
double [] scales = scales_xy[nScene][npix];
double radius = elev_radius[npix];
/*
int diameter = (int) Math.floor(2*radius);
int diameter2 = diameter * diameter;
if (wnd_x.length < (diameter+1)) {
wnd_x = new double[diameter+1];
if (need_derivs) {
dwnd_x = new double[diameter+1];
}
}
*/
double [] scales = scales_xy[nScene][fnpix];
double radius = elev_radius[fnpix];
Arrays.fill(wnd_x, Double.NaN);
double px = x + scales[0] * elevation; // ELEV-SIGN
double py = y + scales[1] * elevation; // ELEV-SIGN
......@@ -2069,7 +1893,7 @@ public class VegetationLMA {
if (need_derivs) {
double dww = wnd_y * dwnd_x[dx] + dwnd_y*wnd_x[dx];
elev_dweights4[nScene][wvindex][dindx] = dww;
int par_indx_elev = par_index[TVAO_ELEVATION][npix]; // ipx+ipy*full.width];
int par_indx_elev = par_index[TVAO_ELEVATION][fnpix]; // ipx+ipy*full.width];
elev_dsum_weights[nScene][par_indx_elev-ind_pars[TVAO_ELEVATION]][wpindex] += dww;
contrib_thread.get(wpindex).add(par_indx_elev); // full parameter index
// pivot_thread[nthread][wvindex][par_indx_elev] = true;
......@@ -2919,13 +2743,13 @@ public class VegetationLMA {
int nthread = ati.getAndIncrement();
for (int yIndx = ai.getAndIncrement(); yIndx < y_src.length; yIndx = ai.getAndIncrement()){
int nscene = y_src[yIndx][YSRC_SCENE];
int indx = y_src[yIndx][YSRC_FINDEX];
int windx = getWindexFromFindex(indx);
int findx = y_src[yIndx][YSRC_FINDEX];
int windx = getWoiIndex(null, woi, full, findx); // getWindexFromFindex(findx);
if (windx == dbg_windx) {
System.out.println("getFxDerivs() SAMPLES_Y_AVG-1 windx="+windx);
}
if (debug_print && (windx == dbg_windx)) {
System.out.println("yIndx="+yIndx+", nscene="+nscene+", indx="+indx+", windx="+windx+", fX[yIndx]="+fX[yIndx]);
System.out.println("yIndx="+yIndx+", nscene="+nscene+", indx="+findx+", windx="+windx+", fX[yIndx]="+fX[yIndx]);
}
double w = recalc_average ? weights[yIndx]/weight_to_scene_weight: scene_weights[nscene]; // norm_scene_weights[nscene];
w *= terrain_correction;
......@@ -2971,7 +2795,8 @@ public class VegetationLMA {
if (samples_pointers[SAMPLES_ALPHA_PULL][1] > 0) { // fits[TVAO_ALPHA] && ((alpha_loss > 0) || (alpha_push > 0))) {
int dbg_nx = -120155; // 30578; // -120114; // -76340;
final boolean dbg_old = false;
final double alpha_threshold = alpha_offset + alpha_min_veg *(1 - 2*alpha_offset);
// final double alpha_threshold = alpha_offset + alpha_min_veg *(1 - 2*alpha_offset);
final double alpha_threshold = alpha_0offset + alpha_min_veg *(1 - alpha_offset - alpha_0offset);
final double alpha_plateau = alpha_threshold; // * alpha_loss;
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
......@@ -2991,14 +2816,14 @@ public class VegetationLMA {
fX[nx] = 0.0;
if (alpha_loss > 0) {
double alpha = vector[np];
if (alpha < alpha_offset) {
double d = alpha - alpha_offset; // negative
if (alpha < alpha_0offset) {
double d = alpha - alpha_0offset; // negative
dquad = - d * d; // negative
ddquad = - 2 * d; // positive
dlin = d; // negative
ddlin = 1; // positive
} else if (alpha < alpha_threshold) {
dlin = (alpha - alpha_offset); // positive
dlin = (alpha - alpha_0offset); // positive
ddlin = 1; // positive
} else if (alpha <= (1 - alpha_offset)) {
dlin = alpha_plateau;
......@@ -3017,8 +2842,8 @@ public class VegetationLMA {
if (dbg_old) {
// old version
double d = 0;
if (alpha < alpha_offset) {
d = alpha- alpha_offset;
if (alpha < alpha_0offset) {
d = alpha- alpha_0offset;
} else if (alpha > (1 - alpha_offset)) {
d = alpha - (1.0 - alpha_offset);
}
......@@ -3030,50 +2855,15 @@ public class VegetationLMA {
}
}
/*
if (alpha_min_veg >= 0) {
if (alpha < alpha_threshold) {
fX[nx] = (alpha - alpha_offset) * alpha_loss;
if (jt != null) {
jt[np][nx] = alpha_loss; // d/dalpha[i]
}
} else if (alpha > (1 - alpha_offset)) {
fX[nx] = (alpha - (1.0 - alpha_offset)) * alpha_loss + alpha_plateau;
if (jt != null) {
jt[np][nx] = alpha_loss; // d/dalpha[i]
}
} else {
fX[nx] = alpha_plateau;
// jt[np][nx] = 0.0; // d/dalpha[i]
}
} else { // old version, requires alpha_min_veg < 0
if (alpha < alpha_offset) {
d = alpha- alpha_offset;
} else if (alpha > (1 - alpha_offset)) {
d = alpha - (1.0 - alpha_offset);
}
if (d != 0) {
fX[nx] = d * d * alpha_loss;
if (jt != null) {
jt[np][nx] = 2 * alpha_loss * d; // d/dalpha[i]
}
}
}
*/
}
if (alpha_push > 0) {
// average including center:
/// double sum_w = (nn + alpha_push_center);
/// double avg_c = (avg * nn + vector[np] * alpha_push_center) / sum_w;
double sum_w = 1; // (nn + alpha_push_center);
double avg_c = vector[np];
double one_minus_n = 1.0 - alpha_push_neutral;
/// double one_minus_p = 1.0 - alpha_push;
double x_minus_p = avg_c - alpha_push;
// double x_minus_p = avg_c - alpha_push;
if ((avg_c > 0) || (avg_c < 1.0)) {
if (avg_c <= alpha_push_neutral) {
fX[nx] += alpha_push / (2 * alpha_push_neutral) * avg_c * (1 - avg_c / (2 * alpha_push_neutral));
......@@ -3100,8 +2890,6 @@ public class VegetationLMA {
if (samples_pointers[SAMPLES_ALPHA_LPF][1] > 0) { // fits[TVAO_ALPHA] && (alpha_lpf >= 0)) {
int dbg_nx = -76340;
// final int ind_y_alpha_lpf = samples_pointers[SAMPLES_ALPHA_LPF][0]; // ind_next;
// ind_next += num_pars[TVAO_ALPHA];
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
......@@ -3264,8 +3052,6 @@ public class VegetationLMA {
}
if (samples_pointers[SAMPLES_TERRAIN_LPF][1]>0) { // fits[TVAO_TERRAIN] && (terr_lpf >= 0)) {
// final int ind_y_terr = samples_pointers[SAMPLES_TERRAIN_LPF][0]; // ind_next;
// ind_next += num_pars[TVAO_TERRAIN];
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
......@@ -3310,8 +3096,6 @@ public class VegetationLMA {
}
if (samples_pointers[SAMPLES_VEGETATION_PULL][1] > 0) { // 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 = samples_pointers[SAMPLES_VEGETATION_PULL][0]; // ind_next;
// ind_next += num_pars[TVAO_VEGETATION];
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
......@@ -3339,54 +3123,6 @@ public class VegetationLMA {
ImageDtt.startAndJoin(threads);
}
/*
if (samples_pointers[SAMPLES_VEGETATION_LPF][1] > 0) { // fits[TVAO_VEGETATION] && (veget_lpf >= 0)) { // should be positive for pull0 and terr_pull_cold (difference between vegetation and terrain)
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[TVAO_VEGETATION]; n = ai.getAndIncrement()) {
int np = n + samples_pointers[SAMPLES_VEGETATION_LPF][2];
int nx = n + samples_pointers[SAMPLES_VEGETATION_LPF][0];
double d = 0;
if (veget_lpf > 0) {
double avg = 0;
int nn = 0;
for (int i = 0; i < veget_neibs[n].length; i++) { // now 4, may be increased
int di = veget_neibs[n][i];
d=0;
if (di >= 0) {
d = vector[di]; // d - full parameter index
avg+=d;
nn++;
} else if ((di < -1) && lpf_fixed[TVAO_VEGETATION]) {
d = tvao[TVAO_VEGETATION][-di - 2];
avg+=d;
nn++;
}
}
avg /= nn; // average
fX[nx] += veget_lpf * (vector[np] - avg); // + veget_pull0 * vector[np];
if (jt != null) {
jt[np][nx] += veget_lpf;
for (int i = 0; i < veget_neibs[n].length; i++) { // now 4, may be increased
int di = veget_neibs[n][i];
if (di >= 0) {
jt[di][nx] -= veget_lpf/nn;
}
}
}
}
}
}
};
}
ImageDtt.startAndJoin(threads);
}
*/
if (samples_pointers[SAMPLES_VEGETATION_LPF][1] > 0) { // should be positive for pull0 and terr_pull_cold (difference between vegetation and terrain)
if (num_pars[TVAO_ALPHA] <= 0) {
System.out.println("getFxDerivs(): BUG: Asjusting SAMPLES_VEGETATION_LPF requires alpha parameters aslo to be fitted");
......@@ -3433,9 +3169,6 @@ public class VegetationLMA {
if (veget_lap >=0) { // only for cold spots
continue;
}
// fX[nx] += veget_lpf * (vector[np] - veget_avg4); // + veget_pull0 * vector[np];
// double w_alpha_neib = 0.7;
// double holes_pwr = 2.0;
d = vector[npa];
double alpha_avg5;
if (alpha_piece_linear) {
......@@ -3519,16 +3252,20 @@ public class VegetationLMA {
for (int n = ai.getAndIncrement(); n < num_pars[TVAO_ELEVATION]; n = ai.getAndIncrement()) {
int np = n + samples_pointers[SAMPLES_ELEVATION_PULL][2];
int nx = n + samples_pointers[SAMPLES_ELEVATION_PULL][0];
double d = 0;
// double d = 0;
if (elevation_pull0 > 0) {
double eff_elevation_pull0 = elevation_pull0;
int findx = par_rindex[np][1];
double elevation_pull = tvao[TVAO_ELEVATION][findx]; // vegetation_pull[findx]; // maybe use tvao[TVAO_TERRAIN][findx] ?
if (Double.isNaN(elevation_pull)) {
elevation_pull= 0; // should not happen unless too far from the vegetation
}
fX[nx] += elevation_pull0 * (vector[np] - elevation_pull);
if (elev_scale_pull && (elevation_pull > elev_scale_thresh)) {
eff_elevation_pull0 *= elev_scale_thresh /elevation_pull;
}
fX[nx] += eff_elevation_pull0 * (vector[np] - elevation_pull);
if (jt != null) {
jt[np][nx] += elevation_pull0;
jt[np][nx] += eff_elevation_pull0;
}
}
}
......@@ -3547,7 +3284,6 @@ public class VegetationLMA {
for (int n = ai.getAndIncrement(); n < num_pars[TVAO_ELEVATION]; n = ai.getAndIncrement()) {
int np = n + samples_pointers[SAMPLES_ELEVATION_LPF][2];
int nx = n + samples_pointers[SAMPLES_ELEVATION_LPF][0];
int findx = par_rindex[np][1];
double d = 0;
if (elevation_lpf > 0) {
double avg = 0;
......@@ -3566,51 +3302,24 @@ public class VegetationLMA {
}
}
avg /= nn; // average vegetation
// if (elevation_lpf_pow <=0) { // without transparency correction
fX[nx] += elevation_lpf * (vector[np] - avg); // + veget_pull0 * vector[np];
if (jt != null) {
jt[np][nx] += elevation_lpf;
for (int i = 0; i < elevation_neibs[n].length; i++) { // now 4, may be increased
int di = elevation_neibs[n][i];
if (di >= 0) {
jt[di][nx] -= elevation_lpf/nn;
}
}
}
/*
} else {
int npa = par_index[TVAO_ALPHA][findx]; // parameter for same-pixel alpha
int na = npa - ind_pars[TVAO_ALPHA];
if (npa < 0) {
System.out.println("getFxDerivs(): BUG: no alpha for vegetation, n="+n+", np="+np+", nx ="+nx+
" findx="+findx+", npa="+npa);
continue;
double eff_elevation_lpf = elevation_lpf;
if (elev_scale_lpf) {
int findx = par_rindex[np][1];
double elevation_pull = tvao[TVAO_ELEVATION][findx]; // vegetation_pull[findx]; // maybe use tvao[TVAO_TERRAIN][findx] ?
if (elevation_pull > elev_scale_thresh) {
eff_elevation_lpf *= elev_scale_thresh /elevation_pull;
}
double alpha_center;
if (alpha_piece_linear) {
alpha_center = (d < 0)? 0: ((d > 1) ? 1.0: d);
} else {
alpha_center = (d < 0)? 0: ((d > 1) ? 1.0: 0.5 * (1.0 - Math.cos(d*Math.PI)));
}
double transparency_pwr = Math.pow(1.0-alpha_center, elevation_lpf_pow); // 0.. 1
double lpf = (1.0 - transparency_pwr) * elevation_lpf + transparency_pwr * elevation_lpf_transp;
fX[nx] += lpf * (vector[np] - avg); // + veget_pull0 * vector[np];
fX[nx] += eff_elevation_lpf * (vector[np] - avg); // + veget_pull0 * vector[np];
if (jt != null) {
jt[np][nx] += lpf;
jt[np][nx] += eff_elevation_lpf;
for (int i = 0; i < elevation_neibs[n].length; i++) { // now 4, may be increased
int di = elevation_neibs[n][i];
if (di >= 0) {
jt[di][nx] -= lpf/nn;
}
jt[di][nx] -= eff_elevation_lpf/nn;
}
// new derivative with respect to alpha
}
}
*/
}
}
}
......@@ -4283,6 +3992,7 @@ public class VegetationLMA {
imp.setProperty("ALPHA_LOSS", ""+alpha_loss);
imp.setProperty("ALPHA_LOSS_LIN", ""+alpha_loss_lin);
imp.setProperty("ALPHA_OFFSET", ""+alpha_offset);
imp.setProperty("ALPHA_0OFFSET", ""+alpha_0offset);
imp.setProperty("ALPHA_MIN_VEG", ""+alpha_min_veg);
imp.setProperty("ALPHA_LPF", ""+alpha_lpf);
imp.setProperty("ALPHA_PIECE_LINEAR", ""+alpha_piece_linear);
......@@ -4316,6 +4026,19 @@ public class VegetationLMA {
imp.setProperty("ELEV_LOW_VEGET", ""+low_veget);
imp.setProperty("SCENES_PULL0", ""+scenes_pull0);
imp.setProperty("SCALE_SCENES_PULL", ""+scale_scenes_pull);
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
imp.setProperty("ELEV_SCALE_THRESH", ""+elev_scale_thresh);
imp.setProperty("ELEV_SCALE_PULL", ""+elev_scale_pull);
imp.setProperty("ELEV_SCALE_LPF", ""+elev_scale_lpf);
// tree-top removal
imp.setProperty("TTOP_EN", ""+ttop_en);
imp.setProperty("TTOP_GB", ""+ttop_gb);
imp.setProperty("TTOP_MIN", ""+ttop_min);
imp.setProperty("TTOP_REL_LEV", ""+ttop_rel_lev);
imp.setProperty("TTOP_REL_RAD", ""+ttop_rel_rad);
imp.setProperty("TTOP_FRAC", ""+ttop_frac);
imp.setProperty("TTOP_REM_RAD", ""+ttop_rem_rad);
imp.setProperty("BOOST_PARALLAX", ""+boost_parallax);
imp.setProperty("MAX_PARALLAX", ""+max_parallax);
......@@ -4471,7 +4194,8 @@ public class VegetationLMA {
alpha_loss = Double.parseDouble((String) imp_pars.getProperty("ALPHA_LOSS"));
alpha_loss_lin = getProperty(imp_pars,"ALPHA_LOSS_LIN", alpha_loss_lin);
alpha_offset = Double.parseDouble((String) imp_pars.getProperty("ALPHA_OFFSET"));
alpha_offset = getProperty(imp_pars,"ALPHA_OFFSET", alpha_offset);
alpha_0offset = getProperty(imp_pars,"ALPHA_0OFFSET", alpha_0offset);
alpha_min_veg = getProperty(imp_pars,"ALPHA_MIN_VEG",alpha_min_veg);
alpha_lpf = Double.parseDouble((String) imp_pars.getProperty("ALPHA_LPF"));
alpha_piece_linear = Boolean.parseBoolean((String) imp_pars.getProperty("ALPHA_PIECE_LINEAR"));
......@@ -4498,18 +4222,31 @@ public class VegetationLMA {
terr_pull_avg = Double.parseDouble((String) imp_pars.getProperty("TERR_PULL_AVG"));
veget_pull0 = Double.parseDouble((String) imp_pars.getProperty("VEGET_PULL0"));
elevation_pull0 = Double.parseDouble((String) imp_pars.getProperty("ELEVATION_PULL0"));
if (imp_pars.getProperty("ELEV_ALPHA_EN") != null) elev_alpha_en = Boolean.parseBoolean((String) imp_pars.getProperty("ELEV_ALPHA_EN"));
if (imp_pars.getProperty("ELEV_ALPHA") != null) elev_alpha = Double.parseDouble((String) imp_pars.getProperty("ELEV_ALPHA"));
if (imp_pars.getProperty("ELEV_ALPHA_PWR") != null) elev_alpha_pwr =Double.parseDouble((String) imp_pars.getProperty("ELEV_ALPHA_PWR"));
if (imp_pars.getProperty("ELEV_LOW_VEGET") != null) low_veget = Double.parseDouble((String) imp_pars.getProperty("ELEV_LOW_VEGET"));
scenes_pull0 = Double.parseDouble((String) imp_pars.getProperty("SCENES_PULL0"));
scale_scenes_pull = Double.parseDouble((String) imp_pars.getProperty("SCALE_SCENES_PULL"));
boost_parallax = Double.parseDouble((String) imp_pars.getProperty("BOOST_PARALLAX"));
// renew above using getProperty(,,);
elev_alpha_en= getProperty(imp_pars,"ELEV_ALPHA_EN", elev_alpha_en);
elev_alpha = getProperty(imp_pars,"ELEV_ALPHA", elev_alpha);
elev_alpha_pwr = getProperty(imp_pars,"ELEV_ALPHA_PWR", elev_alpha_pwr);
low_veget = getProperty(imp_pars,"ELEV_LOW_VEGET", low_veget);
scenes_pull0 = getProperty(imp_pars,"SCENES_PULL0", scenes_pull0);
scale_scenes_pull = getProperty(imp_pars,"SCALE_SCENES_PULL", scale_scenes_pull);
elev_scale_thresh = getProperty(imp_pars,"ELEV_SCALE_THRESH", elev_scale_thresh);
elev_scale_pull = getProperty(imp_pars,"ELEV_SCALE_PULL", elev_scale_pull);
elev_scale_lpf = getProperty(imp_pars,"ELEV_SCALE_LPF", elev_scale_lpf);
ttop_en = getProperty(imp_pars,"TTOP_EN", ttop_en);
ttop_gb = getProperty(imp_pars,"TTOP_GB", ttop_gb);
ttop_min = getProperty(imp_pars,"TTOP_MIN", ttop_min);
ttop_rel_lev = getProperty(imp_pars,"TTOP_REL_LEV", ttop_rel_lev);
ttop_rel_rad = getProperty(imp_pars,"TTOP_REL_RAD", ttop_rel_rad);
ttop_frac = getProperty(imp_pars,"TTOP_FRAC", ttop_frac);
ttop_rem_rad = getProperty(imp_pars,"TTOP_REM_RAD", ttop_rem_rad);
boost_parallax = getProperty(imp_pars,"BOOST_PARALLAX", boost_parallax);
max_parallax = getProperty(imp_pars,"MAX_PARALLAX", max_parallax);
max_warp = getProperty(imp_pars,"MAX_WARP", max_warp);
max_elevation = getProperty(imp_pars,"MAX_ELEVATION", max_elevation);
elevation_radius = getProperty(imp_pars,"ELEVATION_RADIUS", elevation_radius);
terrain_offset = getProperty(imp_pars,"TERRAIN_OFFSET", terrain_offset);
// if (imp_pars.getProperty("TERRAIN_OFFSET") != null) terrain_offset = Double.parseDouble((String) imp_pars.getProperty("TERRAIN_OFFSET"));
......@@ -4687,7 +4424,8 @@ public class VegetationLMA {
getProperty(imp_pars, "REG_WEIGHTS", reg_weights), // final double reg_weights, // fraction of the total weight used for regularization
getProperty(imp_pars, "ALPHA_LOSS", alpha_loss), // final double alpha_loss, // alpha quadratic growing loss for when out of [0,1] range
getProperty(imp_pars, "ALPHA_LOSS_LIN", alpha_loss_lin), // final double alpha_loss_lin, // alpha linear growing loss for when out of [0,1] range and below minimal vegetation alpha
getProperty(imp_pars, "ALPHA_OFFSET", alpha_offset), // final double alpha_offset, // quadratic loss when alpha reaches -1.0 or 2.0
getProperty(imp_pars, "ALPHA_OFFSET", alpha_offset), // final double alpha_offset, // quadratic loss when alpha > 1 - alpha_offset
getProperty(imp_pars, "ALPHA_0OFFSET", alpha_0offset), // final double alpha_0offset, // quadratic loss when alpha < alpha0_offset
getProperty(imp_pars, "ALPHA_MIN_VEG", alpha_min_veg), // final double alpha_min_veg, // 0.5; // if (alpha-alpha_offset)/(1-2*alpha_offset) < alpha_min_veg, pull down to lpha_offset
getProperty(imp_pars, "ALPHA_LPF", alpha_lpf), // final double alpha_lpf, // pull to average of 4 neighbors
......@@ -4712,6 +4450,18 @@ public class VegetationLMA {
getProperty(imp_pars, "ELEV_ALPHA_PWR", elev_alpha_pwr), // final double elev_alpha_pwr, // 2.0; // raise alpha to this power (when alpha > 0)
getProperty(imp_pars, "ELEV_LOW_VEGET", low_veget), // final double low_veget, // 2.0; // (pix) Elevation considered low (lower loss for high alpha)
getProperty(imp_pars, "SCENES_PULL0", scenes_pull0), // final double scenes_pull0,
getProperty(imp_pars, "ELEV_SCALE_THRESH", elev_scale_thresh), // final double elev_scale_thresh,
getProperty(imp_pars, "ELEV_SCALE_PULL", elev_scale_pull), // final boolean elev_scale_pull,
getProperty(imp_pars, "ELEV_SCALE_LPF", elev_scale_lpf), // final boolean elev_scale_lpf,
getProperty(imp_pars, "TTOP_EN", ttop_en), // final boolean ttop_en,
getProperty(imp_pars, "TTOP_GB", ttop_gb), // final double ttop_gb,
getProperty(imp_pars, "TTOP_MIN", ttop_min), // final double ttop_min,
getProperty(imp_pars, "TTOP_REL_LEV", ttop_rel_lev), // final double ttop_rel_lev,
getProperty(imp_pars, "TTOP_REL_RAD", ttop_rel_rad), // final double ttop_rel_rad,
getProperty(imp_pars, "TTOP_FRAC", ttop_frac), // final double ttop_frac,
getProperty(imp_pars, "TTOP_REM_RAD", ttop_rem_rad), // final double ttop_rem_rad,
getProperty(imp_pars, "BOOST_PARALLAX", boost_parallax), // final double boost_parallax, // increase weight of scene with maximal parallax relative to the reference scene
getProperty(imp_pars, "MAX_PARALLAX", max_parallax), // final double max_parallax, // do not consider maximal parallax above this (consider it a glitch)
debugLevel, // final int debugLevel);
......@@ -4985,14 +4735,23 @@ public class VegetationLMA {
public VegetationSegment [] readAllSegments(
String dir_path,
String suffix) {
String suffix,
double transparency_opaque,
double transparency_pedestal,
double transparency_frac,
double transparency_dist,
double transparency_pow,
double transparency_gb,
// TODO: See if pass ttop_* parameters, common to all instead of the ones from the parameters files.
int debugLevel,
String debug_path,
boolean debug_save_improved, // Save debug image after successful LMA step.");
boolean debug_save_worsened) { // Save debug image after unsuccessful LMA step.");
FileFilter parFileFilter = new FileFilter()
{
//Override accept method
public boolean accept(File file) {
//if the file extension is .log return true, else false
public boolean accept(File file) { //if the file extension is .log return true, else false
if (file.isFile() && file.getName().endsWith(suffix+PAR_EXT)) {
return true;
}
......@@ -5002,73 +4761,12 @@ public class VegetationLMA {
File dir = new File(dir_path);
File [] par_files = dir.listFiles(parFileFilter);
VegetationSegment [] segments = new VegetationSegment [par_files.length];
double [] other_pars = new double [1];
for (int n = 0; n < par_files.length; n++) {
Rectangle [] wois = new Rectangle[2];
double [][] tvao = restoreParametersFile( // already trimmed to actual lengths FIXME: Not finished for real import !
par_files[n].getPath(), // String path,
false, // boolean keep_settings,
wois, // Rectangle file_woi)
other_pars); //double [] other_pars)
Rectangle woi_veg = wois[0];
Rectangle woi = wois[1];
segments[n] = new VegetationSegment(
par_files[n].getPath(), // String path,
woi_veg, // Rectangle woi_veg,
woi, // Rectangle woi,
other_pars[0], // terrain_offset,
tvao[TVAO_SCENE_OFFSET], // double [] scene_offsets,
new double [][] {tvao[0],tvao[1],tvao[2],tvao[2]}, // double [][] tva) {} {}
null, // double [] confidence
overlaid); //boolean [] overlaid)
int dbg_n = 13;
boolean bug_check = false; // true;
}
// sort in line-scan order
Arrays.sort(segments, new Comparator<VegetationSegment>() {
@Override
public int compare(VegetationSegment lhs, VegetationSegment rhs) {
return (rhs.woi_veg.y > lhs.woi_veg.y) ? -1 : (rhs.woi_veg.y < lhs.woi_veg.y) ? 1 :
((rhs.woi_veg.x > lhs.woi_veg.x) ? -1 : (rhs.woi_veg.x < lhs.woi_veg.x) ? 1 : 0); // increasing
}
});
return segments;
}
public VegetationSegment [] readAllSegments(
String dir_path,
String suffix,
double transparency_opaque,
double transparency_pedestal,
double transparency_frac,
double transparency_dist,
double transparency_pow,
double transparency_gb,
int debugLevel,
String debug_path,
boolean debug_save_improved, // Save debug image after successful LMA step.");
boolean debug_save_worsened) { // Save debug image after unsuccessful LMA step.");
FileFilter parFileFilter = new FileFilter()
{
//Override accept method
public boolean accept(File file) { //if the file extension is .log return true, else false
if (file.isFile() && file.getName().endsWith(suffix+PAR_EXT)) {
return true;
}
return false;
}
};
File dir = new File(dir_path);
File [] par_files = dir.listFiles(parFileFilter);
VegetationSegment [] segments = new VegetationSegment [par_files.length];
int dbg_n = 13;
boolean bug_check = false; // true;
for (int n = 0; n < par_files.length; n++) {
if (debugLevel > -3) {
System.out.println("====== readAllSegments(): n="+n+ ", "+par_files[n].getPath());
for (int n = 0; n < par_files.length; n++) {
if (debugLevel > -3) {
System.out.println("====== readAllSegments(): n="+n+ ", "+par_files[n].getPath());
}
if (n==dbg_n) {
System.out.println("readAllSegments(): n="+n);
......@@ -5094,6 +4792,13 @@ public class VegetationLMA {
transparency_dist, // final double transparency_dist,
transparency_pow, // final double transparency_pow,
transparency_gb, // final double transparency_gb,
ttop_en, // final boolean ttop_en, // false; // Remove tree tops from transparency weights
ttop_gb, // final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
ttop_min, // final double ttop_min, // 3.0; // Minimal tree top elevation
ttop_rel_lev, // final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
ttop_rel_rad, // final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
ttop_frac, // final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
ttop_rem_rad, // final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
null); // transparency_data); // final double [][][] transparency_data) { // double [3][][]
segments[n] = new VegetationSegment(
par_files[n].getPath(), // String path,
......@@ -5370,8 +5075,130 @@ public class VegetationLMA {
return;
}
public double [] getElevation ( // 0.0 in undefined pixels
double [] vector_in,
boolean only_parameters) {
final double [] vector = (vector_in != null) ? vector_in : parameters_vector;
final int woi_veg_length = woi_veg.width*woi_veg.height;
final double [] elevations = new double [woi_veg_length];
final boolean elevation_parameters = (num_pars[TVAO_ELEVATION] > 0);
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int wvindex = ai.getAndIncrement(); wvindex < woi_veg_length; wvindex = ai.getAndIncrement()) {
int wvx = wvindex % woi_veg.width; // relative to woi_veg
int wvy = wvindex / woi_veg.width; // relative to woi_veg
int x = wvx + woi_veg.x; // relative to full
int y = wvy + woi_veg.y; // relative to full
int fpix = x + y * full.width;
if (elevation_parameters) {
int npar = par_index[TVAO_ELEVATION][fpix];
if (npar >= 0) {
elevations[wvindex] = vector[npar];
} else if (!only_parameters) {
elevations[wvindex] = tvao[TVAO_ELEVATION][fpix]; // use default
}
} else {
elevations[wvindex] = tvao[TVAO_ELEVATION][fpix];
}
if (Double.isNaN(elevations[wvindex])) {
elevations[wvindex] = 0;
}
}
}
};
}
ImageDtt.startAndJoin(threads);
return elevations;
}
public double [] getAlpha ( // NaN in undefined pixels
double [] vector_in,
boolean only_parameters) {
final double [] vector = (vector_in != null) ? vector_in : parameters_vector;
final int woi_veg_length = woi_veg.width*woi_veg.height;
final double [] alphas = new double [woi_veg_length];
Arrays.fill(alphas, Double.NaN);
final boolean alpha_parameters = (num_pars[TVAO_ALPHA] > 0);
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int wvindex = ai.getAndIncrement(); wvindex < woi_veg_length; wvindex = ai.getAndIncrement()) {
int wvx = wvindex % woi_veg.width; // relative to woi_veg
int wvy = wvindex / woi_veg.width; // relative to woi_veg
int x = wvx + woi_veg.x; // relative to full
int y = wvy + woi_veg.y; // relative to full
int fpix = x + y * full.width;
if (alpha_parameters) {
int npar = par_index[TVAO_ALPHA][fpix];
if (npar >= 0) {
alphas[wvindex] = Math.min(Math.max(vector[npar],0),1.0);
} else if (!only_parameters) {
alphas[wvindex] = tvao[TVAO_ALPHA][fpix]; // use default
}
} else {
alphas[wvindex] = tvao[TVAO_ALPHA][fpix];
}
if (!only_parameters && Double.isNaN(alphas[wvindex])) {
alphas[wvindex] = 0;
}
}
}
};
}
ImageDtt.startAndJoin(threads);
return alphas;
}
public double [] getVegetation ( // NaN in undefined pixels
double [] vector_in,
boolean only_parameters) {
final double [] vector = (vector_in != null) ? vector_in : parameters_vector;
final int woi_veg_length = woi_veg.width*woi_veg.height;
final double [] vegetations = new double [woi_veg_length];
Arrays.fill(vegetations, Double.NaN);
final boolean vegetation_parameters = (num_pars[TVAO_VEGETATION] > 0);
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int wvindex = ai.getAndIncrement(); wvindex < woi_veg_length; wvindex = ai.getAndIncrement()) {
int wvx = wvindex % woi_veg.width; // relative to woi_veg
int wvy = wvindex / woi_veg.width; // relative to woi_veg
int x = wvx + woi_veg.x; // relative to full
int y = wvy + woi_veg.y; // relative to full
int fpix = x + y * full.width;
if (vegetation_parameters) {
int npar = par_index[TVAO_VEGETATION][fpix];
if (npar >= 0) {
vegetations[wvindex] = vector[npar];
} else if (!only_parameters) {
vegetations[wvindex] = tvao[TVAO_VEGETATION][fpix]; // use default
}
} else {
vegetations[wvindex] = tvao[TVAO_VEGETATION][fpix];
}
if (!only_parameters && Double.isNaN(vegetations[wvindex])) {
vegetations[wvindex] = 0;
}
}
}
};
}
ImageDtt.startAndJoin(threads);
return vegetations;
}
public double [][] getTransparency(// [scene][windx]
final double [] vector){
double [] vector_in,
boolean only_parameters) {
final double [] vector = (vector_in != null) ? vector_in : parameters_vector;
final int woi_length = woi.width*woi.height;
final double [][] transparency = new double [num_scenes][woi_length];
for (int nscene = 0; nscene < num_scenes; nscene++) {
......@@ -5392,7 +5219,7 @@ public class VegetationLMA {
for (int ny = ai.getAndIncrement(); ny < y_src.length; ny = ai.getAndIncrement()) {
int nscene = y_src[ny][YSRC_SCENE];
int findx = y_src[ny][YSRC_FINDEX];
int windx = getWindexFromFindex(findx);
int windx = getWoiIndex(null, woi, full, findx); //getWindexFromFindex(findx);
double t = 1.0 - alpha[ny];
transparency[nscene][windx] = t;
}
......@@ -5440,7 +5267,7 @@ public class VegetationLMA {
if (!Double.isNaN(t_max)) {
double t_max_offs = t_max - transparency_opaque;
if (t_max_offs > 0) {
int fpix = getFindexFromWindex(wPix);
int fpix = getWoiIndex(woi, null, full, wPix); // getFindexFromWindex(wPix);
double [] scale_xy_max = scales_xy[n_max][fpix];
for (int nscene = 0; nscene < transp.length; nscene++) if (!Double.isNaN(transp[nscene])){
double tf = (transp[nscene] - transparency_opaque) / t_max_offs;
......@@ -5511,8 +5338,45 @@ public class VegetationLMA {
final double transparency_dist,
final double transparency_pow,
final double transparency_gb,
final boolean ttop_en, // false; // Remove tree tops from transparency weights
final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
final double ttop_min, // 3.0; // Minimal tree top elevation
final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
final double [][][] transparency_data) { // double [3][][]
double [][] transparency = getTransparency((vector==null)? parameters_vector:vector);
final double [] border_weights = {0.01, 0.03, 0.1, 0.3};
final boolean apply_transparency_fade = true;
// final int woi_veg_length = woi_veg.width * woi_veg.height;
final int woi_length = woi.width * woi.height;
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
final double [][] transparency = getTransparency((vector==null)? parameters_vector:vector, false);
final double [][] transparency_raw = transparency.clone(); // shallow
double [][] transparency_fade = null;
if (apply_transparency_fade) {
transparency_fade = getSceneTransparencyFade(
vector, // final double [] vector, // parameters vector or null
border_weights); // final double [] border_weights);
for (int i = 0; i < transparency.length; i++) {
transparency[i] = new double [woi_length];
}
final double [][] transparency_fade_final =transparency_fade;
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int nScene = ai.getAndIncrement(); nScene < num_scenes; nScene = ai.getAndIncrement()) {
for (int windx = 0; windx < woi_length; windx++) {
transparency[nScene][windx] = transparency_raw[nScene][windx] * transparency_fade_final[nScene][windx];
}
}
}
};
}
ImageDtt.startAndJoin(threads);
}
double [][] transparency_weights_preblur = getTransparenctWeights( // should be no NaNs!
transparency, // final double [][] transparency,
transparency_opaque, // final double transparency_opaque,
......@@ -5521,8 +5385,6 @@ public class VegetationLMA {
transparency_dist, // final double transparency_dist); //0 - do not apply
transparency_pow); // final double transparency_pow){
final double [][] transparency_weights = new double [num_scenes][];
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
if (transparency_gb <= 0) {
for (int n = 0; n < transparency_weights.length; n++) {
transparency_weights[n] = transparency_weights_preblur[n];
......@@ -5552,12 +5414,335 @@ public class VegetationLMA {
transparency_data[0] = transparency;
transparency_data[1] = transparency_weights_preblur;
transparency_data[2] = transparency_weights;
if (transparency_data.length > 4) {
transparency_data[3] = transparency_fade;
transparency_data[4] = transparency_raw;
}
}
double [] confidence = getTerrainConfidence(transparency_weights); // blurred
return confidence;
}
/**
* Calculate weights of vegetation fade near the edges to use for weight reduction
* @param vector parameters vector or null to use current
* @param border_weights array such as {0.01,0.03,0.1,0.3} to specify border fading. All outside
* the defined vegetation pixels will have weight 0.0, their ortho neighbors
* border_weights[0], their ortho neighbors - border_weights[1], etc. The
* inside area weight is 1.0
* @return weight in linescan order inside woi_veg.
*/
public double [] getVegetationFade(
final double [] vector, // parameters vector or null
final double [] border_weights) {
final int woi_veg_length = woi_veg.width * woi_veg.height;
double [] alphas = getAlpha (vector, true);
double [] wnd = new double [woi_veg_length];
boolean [] mask = new boolean [woi_veg_length];
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int n = ai.getAndIncrement(); n < woi_veg_length; n = ai.getAndIncrement()) {
mask[n] = !Double.isNaN(alphas[n]);
if (mask[n]) {
wnd[n] = Double.NaN;
}
}
}
};
}
ImageDtt.startAndJoin(threads);
TileNeibs tn = new TileNeibs(woi_veg.width,woi_veg.height);
for (int step = 0; step < border_weights.length; step++) {
final double w = border_weights[step];
final boolean last = step == border_weights.length -1;
tn.shrinkSelection(
1, // final int shrink, // grow tile selection by 1 over non-background tiles 1: 4 directions, 2 - 8 directions, 3 - 8 by 1, 4 by 1 more
mask, // final boolean [] tiles,
null); // final boolean [] prohibit)
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int n = ai.getAndIncrement(); n < woi_veg_length; n = ai.getAndIncrement()) {
if (Double.isNaN(wnd[n])) {
if (!mask[n]) {
wnd[n] = w;
} else if (last) {
wnd[n] = 1.0;
}
}
}
}
};
}
ImageDtt.startAndJoin(threads);
}
return wnd;
}
double [][] getSceneTransparencyFade(
final double [] vector, // parameters vector or null
final double [] border_weights) {
final int woi_veg_length = woi_veg.width * woi_veg.height;
final double [] vector_fade = (vector == null) ? parameters_vector.clone() : vector.clone();
double [] vegetation_fade = getVegetationFade(
vector, // final double [] vector, // parameters vector or null
border_weights); // final double [] border_weights)
double [] alpha_fade = new double [woi_veg_length];
for (int i = 0; i < woi_veg_length; i++) {
alpha_fade[i] = 1.0 - vegetation_fade[i];
}
// replace alpha with vegetation_fade in both vectro and tvao, saving tvao before and restoring it after
double [] tvao_alpha = copyWoi ( // save tvao[TVAO_VEGETATION] inside woi_veg
full, // final Rectangle woi_src_in,
woi_veg, // final Rectangle woi_dst_in,
full, // final Rectangle woi_full,
tvao[TVAO_VEGETATION], // final double [] data_src,
null); // final double [] data_dst_in)
copyWoi ( // save tvao[TVAO_VEGETATION] inside woi_veg
woi_veg, // final Rectangle woi_src_in,
full, // final Rectangle woi_dst_in,
full, // final Rectangle woi_full,
alpha_fade, // final double [] data_src,
tvao[TVAO_VEGETATION]); // final double [] data_dst_in)
// replace alpha values in vector_fade
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int n = ai.getAndIncrement(); n < num_pars[TVAO_ALPHA]; n = ai.getAndIncrement()) {
int np = n + samples_pointers[SAMPLES_ALPHA_PULL][2]; // index of the alpha parameter
int findx = par_rindex[np][1]; // full pixel index
if (findx >= 0) {
vector_fade[np] = tvao[TVAO_VEGETATION][findx];
}
}
}
};
}
ImageDtt.startAndJoin(threads);
final double [][] transparency_fade = getTransparency(vector_fade, false);
// restore tvao[TVAO_VEGETATION]
copyWoi ( // save tvao[TVAO_VEGETATION] inside woi_veg
woi_veg, // final Rectangle woi_src_in,
full, // final Rectangle woi_dst_in,
full, // final Rectangle woi_full,
tvao_alpha, // final double [] data_src,
tvao[TVAO_VEGETATION]);// final double [] data_dst_in)
return transparency_fade;
}
/**
* Get tree tops mask for woi_veg rectangle or null if there are no tree tops detected in this window
* @param vector parameters vector or null (will use current parameters_vector
* @param ttop_en if false, return null immediately
* @param ttop_gb Elevation Gaussian blur sigma to detect tree tops
* @param ttop_min Minimal tree top elevation
* @param ttop_rel_lev Relative (to the top height) sample level
* @param ttop_rel_rad Relative (to the top height) sample ring radius
* @param ttop_frac Minimal fraction of the ring pixels below sample level
* @param ttop_rem_rad Relative (to the top height) remove transparency radius
* @param dbg_prefix - Generate image if not null
* @param debugLevel - text debug verbosity
* @return array [woi_veg.width * woi_veg.height] for tree tops pixels (to be removed from transparency)
*/
final boolean [] getTreeTops (
final double [] vector, // parameters vector or null
final boolean ttop_en, // false; // Remove tree tops from transparency weights
final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
final double ttop_min, // 3.0; // Minimal tree top elevation
final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
final boolean ttop_no_border,// true; // No maximums on the border allowed
final String dbg_prefix,
final int debugLevel) {
if (!ttop_en) {
return null;
}
final int woi_veg_length = woi_veg.width*woi_veg.height;
final boolean [] tree_tops = new boolean [woi_veg_length];
final double [] elevations = getElevation (vector,false);
TileNeibs tn = new TileNeibs(woi_veg.width, woi_veg.height);
(new DoubleGaussianBlur()).blurDouble(
elevations, //
woi_veg.width,
woi_veg.height,
ttop_gb, // double sigmaX,
ttop_gb, // double sigmaY,
0.01); // double accuracy)
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
final AtomicInteger anmax = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int wvindex = ai.getAndIncrement(); wvindex < woi_veg_length; wvindex = ai.getAndIncrement()) {
find_max: {
double elev_center = elevations[wvindex];
if (elev_center > ttop_min) {
for (int dir = 0; dir < TileNeibs.DIRS; dir++) {
int wvindex1 = tn.getNeibIndex(wvindex, dir);
if (wvindex1 >= 0) {
if (elevations[wvindex1] > elev_center) {
break find_max;
}
}
}
int wx = wvindex % woi_veg.width; // relative to woi_veg
int wy = wvindex / woi_veg.width; // relative to woi_veg
if (ttop_no_border && ((wx == 0) || (wy == 0) || (wx == (woi_veg.width - 1)) || (wy == (woi_veg.height - 1)))) {
break find_max;
}
double level = elev_center * ttop_rel_lev;
double radius = elev_center * ttop_rel_rad;
double remove_radius = elev_center * ttop_rem_rad;
int irad = (int) Math.ceil(radius);
int y0 = Math.max(wy - irad, 0);
int y1 = Math.min(wy + irad, woi_veg.height - 1);
int x0 = Math.max(wx - irad, 0);
int x1 = Math.min(wx + irad, woi_veg.width - 1);
int n = 0, n_lower = 0;
for (int y = y0; y <= y1; y++) {
double dy = y - wy;
double dy2 = dy*dy;
for (int x = x0; x <= x1; x++) {
double dx = x - wx;
double dr = Math.sqrt(dy2 + dx*dx) - radius;
if (Math.round(dr) == 0) {
n ++;
if (elevations[y * woi_veg.width + x] < level) {
n_lower++;
}
}
}
}
double frac = 1.0 * n_lower / n;
if (frac > ttop_frac) { // mark
int iremove_radius = (int) Math.ceil(remove_radius);
int yr0 = Math.max(wy - iremove_radius, 0);
int yr1 = Math.min(wy + iremove_radius, woi_veg.height - 1);
int xr0 = Math.max(wx - iremove_radius, 0);
int xr1 = Math.min(wx + iremove_radius, woi_veg.width - 1);
for (int y = yr0; y <= yr1; y++) {
double dy = y - wy;
double dy2 = dy*dy;
for (int x = xr0; x <= xr1; x++) {
double dx = x - wx;
double r = Math.sqrt(dy2 + dx*dx);
if (r <= remove_radius) {
tree_tops[y * woi_veg.width + x] = true;
}
}
}
anmax.getAndIncrement();
if (debugLevel > -2) {
int findx = getWoiIndex(woi_veg, null, full, wvindex);
int fx = findx % full.width;
int fy = findx / full.width;
System.out.println ("getTreeTops(): found a tree top at xy=("+fx+","+fy+"), elevation = "+
elev_center+ ">= "+ttop_min+ ", ring fraction="+frac+" >= "+ttop_frac+"." );
}
// debug print. Show results?
} else {
}
}
}
}
}
};
}
ImageDtt.startAndJoin(threads);
if (dbg_prefix != null) {
double [] elevations_nonblurred = getElevation (vector, false);
double [] alphas = getAlpha (vector, true); //
double [] vegetations = getVegetation (vector, false);
double [] transparency = new double [woi_veg_length];
double [] ttops = new double [woi_veg_length];
double max_elev = 0;
for (int i = 0; i < elevations.length; i++) {
if (elevations[i] > max_elev) {
max_elev = elevations[i];
}
}
for (int i = 0; i < elevations.length; i++) {
elevations_nonblurred[i] /= max_elev;
elevations[i] /= max_elev;
transparency[i] = 1.0 - alphas[i];
if (tree_tops[i]) {
ttops[i] = 1.0;
transparency[i] = 0;
}
}
// elevations - blurred
String [] dbg_titles = {"elevation","blurred","tree_tops","transparency","alpha","vegetation"};
double [][] dbg_img = {elevations_nonblurred, elevations, ttops, transparency, alphas, vegetations};
ShowDoubleFloatArrays.showArrays(
dbg_img,
woi_veg.width,
woi_veg.height,
true,
dbg_prefix+"_TREE-TOPS-"+anmax.get(),
dbg_titles);
}
if (anmax.get() == 0) {
return null;
}
return tree_tops;
}
/**
* Clone vector and modify it's alpha making tree tops completely opaque (alpha = 1.0)
* @param vector parameters vector or null (will use current parameters vector)
* @param tree_tops boolean array of tree tops pixels, matching Rectangle woi_veg
* @return modified vector with tree tops opaque or the unmodified and not cloned vector if tree_tops in null
*/
public double [] applyTreeTops(
double [] vector,
boolean [] tree_tops) {
if (tree_tops == null) {
return (vector == null) ? parameters_vector : vector;
} else {
final double [] vector_applied = (vector == null) ? parameters_vector.clone() : vector.clone();
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int n = ai.getAndIncrement(); n < num_pars[TVAO_ALPHA]; n = ai.getAndIncrement()) {
int np = n + samples_pointers[SAMPLES_ALPHA_PULL][2]; // index of the alpha parameter
int findx = par_rindex[np][1]; // full pixel index
if (findx >= 0) {
int wvindx = getWoiIndex(null, woi_veg, full, findx); // convert to woi_veg index
if (tree_tops[wvindx]) {
vector_applied[np] = 1.0; // make opaque
}
}
}
}
};
}
ImageDtt.startAndJoin(threads);
return vector_applied;
}
}
public double [] applyTransparency(
final double [] vector,
final double transparency_opaque,
......@@ -5568,25 +5753,59 @@ public class VegetationLMA {
final double transparency_gb,
final double transparency_boost_in,
final boolean recalc_average,
String prefix) {
final boolean ttop_en, // false; // Remove tree tops from transparency weights
final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
final double ttop_min, // 3.0; // Minimal tree top elevation
final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
final boolean ttop_no_border,// true; // No maximums on the border allowed
final String dbg_prefix,
final int debugLevel) {
this.recalc_average = recalc_average;
final int woi_length = woi.width*woi.height;
final double transparency_boost = Math.max(transparency_boost_in, 1.0);
boolean use_hf = (y_src_hf != null);
int y_avg_len = y_vector.length - (y_src.length + ((y_src_hf != null) ? y_src_hf.length : 0));
double [][][] transparency_data = new double [3][][];
double [][][] transparency_data = new double [5][][];
boolean [] tree_tops = getTreeTops (
vector, // final double [] vector, // parameters vector or null
ttop_en, // final boolean ttop_en, // false; // Remove tree tops from transparency weights
ttop_gb, // final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
ttop_min, // final double ttop_min, // 3.0; // Minimal tree top elevation
ttop_rel_lev, // final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
ttop_rel_rad, // final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
ttop_frac, // final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
ttop_rem_rad, // final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
ttop_no_border,//final boolean ttop_no_border,// true; // No maximums on the border allowed
dbg_prefix, // final String dbg_prefix,
debugLevel); // debugLevel); // final int debugLevel) {
double [] vector_mod = applyTreeTops(
vector, // double [] vector,
tree_tops); // boolean [] tree_tops)
double [] confidence = getTransparencyConfidence(
vector, // final double [] vector,
vector_mod, // final double [] vector,
transparency_opaque, // final double transparency_opaque,
transparency_pedestal, // final double transparency_pedestal,
transparency_frac, // final double transparency_frac,
transparency_dist, // final double transparency_dist,
transparency_pow, // final double transparency_pow,
transparency_gb, // final double transparency_gb,
ttop_en, // final boolean ttop_en, // false; // Remove tree tops from transparency weights
ttop_gb, // final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
ttop_min, // final double ttop_min, // 3.0; // Minimal tree top elevation
ttop_rel_lev, // final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
ttop_rel_rad, // final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
ttop_frac, // final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
ttop_rem_rad, // final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
transparency_data); // final double [][][] transparency_data) { // double [3][][]
double [][] transparency = transparency_data[0];
double [][] transparency_weights_preblur = transparency_data[1];
double [][] transparency_weights = transparency_data[2];
double [][] transparency_fade = transparency_data[3];
double [][] transparency_raw = transparency_data[4];
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
......@@ -5619,29 +5838,38 @@ public class VegetationLMA {
}
ImageDtt.startAndJoin(threads);
if (prefix != null) {
String title =prefix+String.format("-opq%3f-pdst%5f-frac%3f-dst%2f-pwr%4f-gb%4f-boost%4f-avg%1d.tiff",
if (dbg_prefix != null) {
String title =dbg_prefix+String.format("-opq%3f-pdst%5f-frac%3f-dst%2f-pwr%4f-gb%4f-boost%4f-avg%1d.tiff",
transparency_opaque,transparency_pedestal,transparency_frac,
transparency_dist,transparency_pow,transparency_gb, transparency_boost, recalc_average?1:0);
int nscenes = transparency.length;
String [] titles_top = {"transparency","pre-blur","weights","normalized"};
String [] titles_top = {"transparency_fade", "transparity_raw", "transparency","pre-blur","weights","normalized"};
String [] titles = new String[nscenes + 1];
double [][][] dbg_img = new double[titles_top.length][titles.length][]; // woi_length];
dbg_img[0][nscenes] = new double[woi_length];
dbg_img[1][nscenes] = getTerrainConfidence(transparency_weights_preblur);
dbg_img[2][nscenes] = getTerrainConfidence(transparency_weights);
dbg_img[3][nscenes] = scales;
dbg_img[1][nscenes] = new double[woi_length];
dbg_img[2][nscenes] = new double[woi_length];
dbg_img[3][nscenes] = getTerrainConfidence(transparency_weights_preblur);
dbg_img[4][nscenes] = getTerrainConfidence(transparency_weights);
dbg_img[5][nscenes] = scales;
for (int nscene = 0; nscene < transparency.length; nscene++) {
titles[nscene] = "scene-"+nscene;
dbg_img[0][nscene] = transparency[nscene];
dbg_img[1][nscene] = transparency_weights_preblur[nscene];
dbg_img[2][nscene] = transparency_weights[nscene];
dbg_img[3][nscene] = scaled_weights[nscene];
dbg_img[0][nscene] = (transparency_fade != null) ? transparency_fade[nscene] : new double[woi_length];
dbg_img[1][nscene] = transparency_raw[nscene];
dbg_img[2][nscene] = transparency[nscene];
dbg_img[3][nscene] = transparency_weights_preblur[nscene];
dbg_img[4][nscene] = transparency_weights[nscene];
dbg_img[5][nscene] = scaled_weights[nscene];
for (int wpix = 0;wpix < woi_length; wpix++) {
if (dbg_img[0][nscenes][wpix] > transparency[nscene][wpix]) { // takes care of NaN
transparency[nscene][wpix] = dbg_img[0][nscenes][wpix];
if (dbg_img[2][nscenes][wpix] < transparency[nscene][wpix]) { // takes care of NaN
dbg_img[2][nscenes][wpix] = transparency[nscene][wpix];
}
if (dbg_img[1][nscenes][wpix] < transparency_raw[nscene][wpix]) { // takes care of NaN
dbg_img[1][nscenes][wpix] = transparency_raw[nscene][wpix];
}
if ((transparency_fade != null) && (dbg_img[0][nscenes][wpix] < transparency_fade[nscene][wpix])) { // takes care of NaN
dbg_img[0][nscenes][wpix] = transparency_fade[nscene][wpix] ;
}
// dbg_img[0][nscenes][wpix] = Math.max(transparency[nscene][wpix],dbg_img[0][nscenes][wpix]);
}
}
titles[nscenes] = "max/conf";
......@@ -5658,7 +5886,7 @@ public class VegetationLMA {
for (int ny = 0; ny < y_src.length; ny++) {
int nscene = y_src[ny][YSRC_SCENE];
int findx = y_src[ny][YSRC_FINDEX];
int windx = getWindexFromFindex(findx);
int windx = getWoiIndex(null, woi, full, findx); // getWindexFromFindex(findx);
double w = scaled_weights[nscene][windx];
if (Double.isNaN(w)) {
System.out.println("applyTransparency() 1: ny="+ny+", nscene="+nscene+", findx="+findx+", windx="+windx+", w="+w);
......@@ -5673,7 +5901,7 @@ public class VegetationLMA {
for (int ny = 0; ny < y_src_hf.length; ny++) {
int nscene = y_src_hf[ny][YSRC_SCENE];
int findx = y_src_hf[ny][YSRC_FINDEX];
int windx = getWindexFromFindex(findx);
int windx = getWoiIndex(null, woi, full, findx); // getWindexFromFindex(findx);
double w = scaled_weights[nscene][windx] * hifreq_weight ;
if (Double.isNaN(w)) {
System.out.println("applyTransparency() 2: ny="+ny+", nscene="+nscene+", findx="+findx+", windx="+windx+", w="+w);
......@@ -5710,7 +5938,7 @@ public class VegetationLMA {
int ny = nw;
int nscene = y_src[ny][YSRC_SCENE];
int findx = y_src[ny][YSRC_FINDEX];
int windx = getWindexFromFindex(findx);
int windx = getWoiIndex(null, woi, full, findx); // getWindexFromFindex(findx);
double w = scaled_weights[nscene][windx] * k;
if (!Double.isNaN(w)) {
weights[nw] = w;
......@@ -5722,7 +5950,7 @@ public class VegetationLMA {
int ny = nw-y_src.length;
int nscene = y_src_hf[ny][YSRC_SCENE];
int findx = y_src_hf[ny][YSRC_FINDEX];
int windx = getWindexFromFindex(findx);
int windx = getWoiIndex(null, woi, full, findx); //; // getWindexFromFindex(findx);
double w = scaled_weights[nscene][windx] * k * hifreq_weight; // scene_weights[nscene] * k * hifreq_weight;
if (!Double.isNaN(w)) {
weights[nw] = w;
......@@ -5752,7 +5980,7 @@ public class VegetationLMA {
if (y_src_hf != null) { // set average offset y - either old way with scene_weights[nscene] or new way with weights[ny]
setupYVectorTransparency();
}
if (prefix != null) {
if (dbg_prefix != null) {
double sw = 0;
for (int n = 0; n < (samples_pointers.length -1); n++) {
double swn = 0;
......@@ -5822,8 +6050,8 @@ public class VegetationLMA {
public void run() {
for (int ny = ai.getAndIncrement(); ny < y_src.length; ny = ai.getAndIncrement()) {
int nscene = y_src[ny][YSRC_SCENE];
int indx = y_src[ny][YSRC_FINDEX];
int windx = getWindexFromFindex(indx);
int findx = y_src[ny][YSRC_FINDEX];
int windx = getWoiIndex(null, woi, full, findx); //getWindexFromFindex(findx);
double d =y_vector[ny] * scene_weights[nscene] * terrain_correction;
if (Double.isNaN(d)) {
System.out.println("y_vector["+ny+"]*scene_weights["+nscene+"]=NaN");
......@@ -5882,8 +6110,8 @@ public class VegetationLMA {
public void run() {
for (int ny = ai.getAndIncrement(); ny < y_src.length; ny = ai.getAndIncrement()) {
int nscene = y_src[ny][YSRC_SCENE];
int indx = y_src[ny][YSRC_FINDEX];
int windx = getWindexFromFindex(indx);
int findx = y_src[ny][YSRC_FINDEX];
int windx = getWoiIndex(null, woi, full, findx); // getWindexFromFindex(findx);
// double d =y_vector[ny] * scene_weights[nscene];
double w = recalc_average ? (weights[ny] / weight_to_scene_weight): scene_weights[nscene];
w *= terrain_correction;
......@@ -5929,16 +6157,6 @@ public class VegetationLMA {
private int getWindexFromFindex(int findex) {
int wx = (findex % full.width) - woi.x;
int wy = (findex / full.width) - woi.y;
return wx + wy*woi.width;
}
private int getFindexFromWindex(int windex) {
int x = (windex % woi.width) + woi.x;
int y = (windex / woi.width) + woi.y;
return x + y*full.width;
}
/*
private void setupYWavg( double [] scene_weights) {
......@@ -6002,16 +6220,16 @@ public class VegetationLMA {
public void run() {
for (int i = ai.getAndIncrement(); i < num_samples; i = ai.getAndIncrement()) {
int pindx = ind_samples + i; // full parameter index
int indx = par_rindex[pindx][1]; // full pixel index
int findx = 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);
int nindx = tn.getNeibIndex(findx, dir);
if (nindx < 0) {
alpha_neibs[i][dir4] = -1;
} else if (par_index[tvao][nindx] >= 0) {
alpha_neibs[i][dir4] = par_index[tvao][nindx];
} else { // pull to fixed alpha value
alpha_neibs[i][dir4] = -2 -indx; // OFFSET by 2 !
alpha_neibs[i][dir4] = -2 -findx; // OFFSET by 2 !
}
}
}
......@@ -7127,6 +7345,104 @@ public class VegetationLMA {
return par_files;
}
/*
private int getWindexFromFindex(int findex) {
int wx = (findex % full.width) - woi.x;
int wy = (findex / full.width) - woi.y;
return wx + wy*woi.width;
}
private int getFindexFromWindex(int windex) {
int x = (windex % woi.width) + woi.x;
int y = (windex / woi.width) + woi.y;
return x + y*full.width;
}
*/
/**
* Convert between indices inside rectangular WOIs
* @param woi_src source WOI (or null if it is the full image WOI)
* @param woi_dst destination WOI (or null if it is the full image WOI)
* @param indx source index inside woi_src
* @param woi_full full WOI (not null)
* @return destination index inside woi_dst
*/
public static int getWoiIndex(
Rectangle woi_src,
Rectangle woi_dst,
Rectangle woi_full,
int indx) {
int x,y;
if (woi_src != null) {
x = indx % woi_src.width + woi_src.x;
y = indx / woi_src.width + woi_src.y;
} else {
x = indx % woi_full.width;
y = indx / woi_full.width;
}
if (woi_dst != null) {
return (x - woi_dst.x) + (y - woi_dst.y) * woi_dst.width;
} else {
return x + y * woi_full.width;
}
}
/**
* Copy data beween WOIs
* @param woi_src_in source WOI (or null if it is the full image WOI)
* @param woi_dst_in destination WOI (or null if it is the full image WOI)
* @param woi_full full WOI (not null)
* @param data_src - contents of woi_src in linescan order
* @param data_dst_in - current contents of woi_dst or null (will use NaN) to copy data to
* @return dst_data
*/
public static double [] copyWoi (
final Rectangle woi_src_in,
final Rectangle woi_dst_in,
final Rectangle woi_full,
final double [] data_src,
final double [] data_dst_in) {
final Rectangle woi_src = (woi_src_in != null) ? woi_src_in : woi_full;
final Rectangle woi_dst = (woi_dst_in != null) ? woi_dst_in : woi_full;
final int woi_dst_length = woi_dst.width * woi_dst.height;
final double [] data_dst = (data_dst_in != null) ? data_dst_in : new double[woi_dst_length];
if (data_dst_in == null) {
Arrays.fill (data_dst, Double.NaN);
}
final Rectangle woi_intersect = woi_src.intersection(woi_dst);
if (!woi_intersect.isEmpty()) {
final Thread[] threads = ImageDtt.newThreadArray();
final AtomicInteger ai = new AtomicInteger(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int irow = ai.getAndIncrement(); irow < woi_intersect.height; irow = ai.getAndIncrement()) {
int inter_indx = irow * woi_intersect.width;
int src_indx = getWoiIndex(
woi_intersect, // Rectangle woi_src,
woi_src, // Rectangle woi_dst,
woi_full, // Rectangle woi_full,
inter_indx); // int indx)
int dst_indx = getWoiIndex(
woi_intersect, // Rectangle woi_src,
woi_dst, // Rectangle woi_dst,
woi_full, // Rectangle woi_full,
inter_indx); // int indx)
System.arraycopy(
data_src,
src_indx,
data_dst,
dst_indx,
woi_intersect.width);
}
}
};
}
ImageDtt.startAndJoin(threads);
}
return data_dst;
}
}
src/main/java/com/elphel/imagej/vegetation/VegetationModel.java
View file @ 6a053711
......@@ -1512,7 +1512,8 @@ public class VegetationModel {
// double default_alpha = clt_parameters.imp.terr_alpha_dflt; // 0.5; // 0.8;
double alpha_loss = clt_parameters.imp.terr_alpha_loss; //100.0; // alpha quadratic growing loss for when out of [0,1] range
double alpha_loss_lin = clt_parameters.imp.terr_alpha_loss_lin; // alpha linear growing loss for when out of [0,1] range and below minimal vegetation alpha
double alpha_offset = clt_parameters.imp.terr_alpha_offset; // 0.0; // 0.02; // 0.03; // if >0, start losses above 0.0 and below 1.0;
double alpha_offset = clt_parameters.imp.terr_alpha_offset; // 0.0; // 0.02; // 0.03; // if >0, start losses below 1.0;
double alpha_0offset = clt_parameters.imp.terr_alpha_0offset; // 0.0; // if >0, start losses above 0.0
double alpha_min_veg = clt_parameters.imp.terr_alpha_min_veg; // 0.5 // Minimal vegetation alpha. If (alpha-alpha_offset)/(1-2*alpha_offset) < alpha_min_veg, pull down to lpha_offset
double alpha_lpf = clt_parameters.imp.terr_alpha_lpf; // 2.5; // 5; /// 2; // 5; /// 10; /// 15; // 10.0; // 5.0; // 10.0; // 3; // 10; // 20; // 6.0; // 3.0; // 2.0; // 1.5; // 5.0; // 0.5; // pull to average of 4 neighbors
......@@ -1538,7 +1539,19 @@ public class VegetationModel {
double elev_alpha_pwr = clt_parameters.imp.terr_elev_alpha_pwr; // 1.0; // multiply alpha by under-low elevation for loss
double low_veget = clt_parameters.imp.terr_low_veget; // 2.0; // (pix) Elevation considered low (lower loss for high alpha)
double scenes_pull0 = clt_parameters.imp.terr_scenes_pull0; // 1.0
// scaling elevation losses for high elevations (decrease pull and/or lpf for high elevations)
double elev_scale_thresh = clt_parameters.imp.terr_elev_scale_thresh; // 1.0; // reduce losses for higher (initial) elevations TODO: consider actual elevations
boolean elev_scale_pull = clt_parameters.imp.terr_elev_scale_pull; // false; // scale elevation pull losses for high elevations
boolean elev_scale_lpf = clt_parameters.imp.terr_elev_scale_lpf; // false; // scale elevation diffusion losses for high elevations
// tree-top removal
boolean ttop_en = clt_parameters.imp.terr_ttop_en; // false; // remove tree tops from transparency weights
double ttop_gb = clt_parameters.imp.terr_ttop_gb; // 1.0; // Elevation Gaussian blur sigma to detect tree tops
double ttop_min = clt_parameters.imp.terr_ttop_min; // 3.0; // Minimal tree top elevation
double ttop_rel_lev = clt_parameters.imp.terr_ttop_rel_lev; // 0.9; // Relative (to the top height) sample level
double ttop_rel_rad = clt_parameters.imp.terr_ttop_rel_rad; // 0.25; // Relative (to the top height) sample ring radius
double ttop_frac = clt_parameters.imp.terr_ttop_frac; // 0.5; // Minimal fraction of the ring pixels below sample level
double ttop_rem_rad = clt_parameters.imp.terr_ttop_rem_rad; // 0.25; // Relative (to the top height) remove transparency radius
boolean ttop_no_border = true; // No maximums on the border allowed
// LMA parameters
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;
......@@ -1706,37 +1719,6 @@ public class VegetationModel {
}
*/
double dir_sigma = 16;
/*
vegetation_warp_md = new double [vegetation_warp.length][][];
for (int nscene = 0; nscene < vegetation_warp.length; nscene++) {
// if (nscene == reference_index) {
// vegetation_warp_md[nscene] = new double [full.width * full.height][];
// } else {
boolean after_ref = (nscene > reference_index);
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
after_ref); // final boolean after_ref);
// }
}
vegetation_inv_warp_md = new double [vegetation_warp.length][][];
for (int nscene = 0; nscene < vegetation_warp.length; nscene++) {
// if (nscene == reference_index) {
// vegetation_inv_warp_md[nscene] = new double [full.width * full.height][];
// } else {
boolean after_ref = (nscene > reference_index);
vegetation_inv_warp_md[nscene ]= warpMagnitudeDirection(
vegetation_inv_warp[nscene], // final double [][] warp_dxy,
full.width, // final int width,
dir_sigma, // final double dir_sigma,// Gaussian blur sigma to smooth direction variations
after_ref); // final boolean after_ref);
// }
}
*/
// TODO: remove some of the above
if ((elevations == null) || (scale_dirs == null)){
//Moving it here to generate needed vegetation_inv_warp_md
......@@ -1744,29 +1726,21 @@ public class VegetationModel {
// probably will not use these, but as optional
vegetation_warp_md = new double [vegetation_warp.length][][];
for (int nscene = 0; nscene < vegetation_warp.length; nscene++) {
// if (nscene == reference_index) {
// vegetation_warp_md[nscene] = new double [full.width * full.height][];
// } else {
boolean after_ref = (nscene > reference_index);
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
after_ref); // final boolean after_ref);
// }
}
vegetation_inv_warp_md = new double [vegetation_warp.length][][];
for (int nscene = 0; nscene < vegetation_warp.length; nscene++) {
// if (nscene == reference_index) {
// vegetation_inv_warp_md[nscene] = new double [full.width * full.height][];
// } else {
boolean after_ref = (nscene > reference_index);
vegetation_inv_warp_md[nscene ]= warpMagnitudeDirection(
vegetation_inv_warp[nscene], // final double [][] warp_dxy,
full.width, // final int width,
dir_sigma, // final double dir_sigma,// Gaussian blur sigma to smooth direction variations
after_ref); // final boolean after_ref);
// }
}
double [][][] dbg_img = new double[12][vegetation_warp.length][vegetation_warp[0].length];
......@@ -1910,10 +1884,8 @@ public class VegetationModel {
}
for (int npix = 0; npix < num_pixels; npix++) {
if (vegetation_inv_warp_md[nscene][npix] != null) {
// compare_offsets[0][nscene][npix] = vegetation_inv_warp_md[nscene][npix][1]; // direction
compare_offsets[1][nscene][npix] = vegetation_inv_warp_md[nscene][npix][0];
}
// if (!Double.isNaN(elevation_scales[nscene][npix])) {
if (elevation_scale_dirs[nscene][npix] != null) {
compare_offsets[0][nscene][npix] = elevation_scale_dirs[nscene][npix][1]; // direction
compare_offsets[2][nscene][npix] = elevation_scale_dirs[nscene][npix][0] * elevations[npix];
......@@ -1931,34 +1903,25 @@ public class VegetationModel {
}
// if ((debugLevel > 3) || um_en) {
if (debugLevel > 3) { //-2) {
// probably will not use these, but as optional
vegetation_warp_md = new double [vegetation_warp.length][][];
for (int nscene = 0; nscene < vegetation_warp.length; nscene++) {
// if (nscene == reference_index) {
// vegetation_warp_md[nscene] = new double [full.width * full.height][];
// } else {
boolean after_ref = (nscene > reference_index);
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
after_ref); // final boolean after_ref);
// }
}
vegetation_inv_warp_md = new double [vegetation_warp.length][][];
for (int nscene = 0; nscene < vegetation_warp.length; nscene++) {
// if (nscene == reference_index) {
// vegetation_inv_warp_md[nscene] = new double [full.width * full.height][];
// } else {
boolean after_ref = (nscene > reference_index);
vegetation_inv_warp_md[nscene ]= warpMagnitudeDirection(
vegetation_inv_warp[nscene], // final double [][] warp_dxy,
full.width, // final int width,
dir_sigma, // final double dir_sigma,// Gaussian blur sigma to smooth direction variations
after_ref); // final boolean after_ref);
// }
}
double [][][] dbg_img = new double[12][vegetation_warp.length][vegetation_warp[0].length];
......@@ -1994,8 +1957,6 @@ public class VegetationModel {
dbg_img[10][nscene][npix] = vegetation_inv_warp_md[nscene][npix][0]; // magnitude
dbg_img[11][nscene][npix] = vegetation_inv_warp_md[nscene][npix][1]; // direction
}
}
}
......@@ -2032,12 +1993,6 @@ public class VegetationModel {
} // if ((elevations == null) || (scale_dirs == null)){
/*
VegetationLMA vegetationLMA = new VegetationLMA (
this,
alpha_initial_contrast, // double alpha_initial_contrast)
debugLevel); // int debugLevel)
*/
VegetationLMA vegetationLMA = new VegetationLMA (
this,
alpha_init_offs, // 0.01; // double alpha_offs,
......@@ -2063,11 +2018,6 @@ public class VegetationModel {
if (run_combine) {
/*
VegetationSegment [] segments = vegetationLMA.readAllSegments(
segments_dir, // String dir_path)
segments_suffix); // String suffix);
*/
VegetationSegment [] segments = vegetationLMA.readAllSegments(
segments_dir, // String dir_path)
segments_suffix, // String suffix);
......@@ -2100,6 +2050,11 @@ public class VegetationModel {
}
Rectangle woi = woi_last_done;
while (true) {
if (debugLevel>-2) {
Runtime runtime = Runtime.getRuntime();
runtime.gc();
System.out.println("--- Free memory="+runtime.freeMemory()+" (of "+runtime.totalMemory()+")");
}
step_restore= par_restore? clt_parameters.imp.terr_step_restore : 0;
woi = nextTileWoi(
woi_enclosing, // Rectangle enclosing,
......@@ -2111,8 +2066,6 @@ public class VegetationModel {
System.out.println ("===== No WOIs to process left, exiting");
break;
}
// String par_path = restore_mode? "RESTORE": (read_pars ? parameters_path : null);
// String par_path = read_pars ? parameters_path : null;
if (tile_woi) {
System.out.println("===== Will process WOI ("+woi.x+", "+woi.y+", "+woi.width+", "+woi.height+") of the enclosing WOI ("+
+woi_enclosing.x+", "+woi_enclosing.y+", "+woi_enclosing.width+", "+woi_enclosing.height+").");
......@@ -2139,7 +2092,8 @@ public class VegetationModel {
reg_weights, // final double reg_weights, // fraction of the total weight used for regularization
alpha_loss, // final double alpha_loss, // alpha quadratic growing loss for when out of [0,1] range
alpha_loss_lin, // final double alpha_loss_lin, // alpha linear growing loss for when out of [0,1] range and below minimal vegetation alpha
alpha_offset, // final double alpha_offset, // quadratic loss when alpha reaches -1.0 or 2.0
alpha_offset, // final double alpha_offset, // quadratic loss when alpha > 1.0 - alpha_offset
alpha_0offset, // final double alpha_0ffset, // quadratic loss when alpha < alpha0_offset
alpha_min_veg, // final double alpha_min_veg, // 0.5; // if (alpha-alpha_offset)/(1-2*alpha_offset) < alpha_min_veg, pull down to lpha_offset
alpha_lpf, // final double alpha_lpf, // pull to average of 4 neighbors
alpha_piece_linear, // final boolean alpha_piece_linear, // true - piece-linear, false - half-cosine
......@@ -2163,11 +2117,20 @@ public class VegetationModel {
elev_alpha_pwr, // final double elev_alpha_pwr, // 2.0; // raise alpha to this power (when alpha > 0)
low_veget, // final double low_veget, // 2.0; // (pix) Elevation considered low (lower loss for high alpha)
scenes_pull0, // final double scenes_pull0,
elev_scale_thresh, // final double elev_scale_thresh, // 1.0;// reduce losses for higher (initial) elevations TODO: consider actual elevations
elev_scale_pull,// final boolean elev_scale_pull,// false; // scale elevation pull losses for high elevations
elev_scale_lpf, // final boolean elev_scale_lpf, // false; // scale elevation diffusion losses for high elevations
ttop_en, // final boolean ttop_en, // false; // remove tree tops from transparency weights
ttop_gb, // final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
ttop_min, // final double ttop_min, // 3.0; // Minimal tree top elevation
ttop_rel_lev, // final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
ttop_rel_rad, // final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
ttop_frac, // final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
ttop_rem_rad, // final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
boost_parallax,// final double boost_parallax, // increase weight of scene with maximal parallax relative to the reference scene
max_parallax, //final double max_parallax, // do not consider maximal parallax above this (consider it a glitch)
// um_sigma, // final double um_sigma, // just use in debug image names
// (um_en? um_weight: 0.0), // final double um_weight,
// par_path, // final String parameters_read_path,
debugLevel, // final int debugLevel);
debug_path, // final String debug_path,
debug_save_improved, // final boolean debug_save_improved, // Save debug image after successful LMA step.");
......@@ -2178,15 +2141,7 @@ public class VegetationModel {
continue;
}
// if (save_par_files && skip_existing_woi) { // check that segment already exists
if (save_par_files && skip_existing_woi) { // check that segment already exists
// public boolean segmentFileExist(String dir_path) {
/*
String save_path = VegetationLMA.getSavePath(
segments_dir, // String dir,
vegetationLMA.getParametersDebugTitle()); // String title)
if (new File (save_path).exists()) {
*/
File[] segment_files = vegetationLMA.getSegmentFiles(segments_dir);
if (tile_woi) {
if (segment_files.length > 0) {
......@@ -2238,26 +2193,19 @@ public class VegetationModel {
transparency_gb, // final double transparency_gb,
transparency_boost, // final double transparency_boost,
recalc_average, // final boolean recalc_average);
dbg_title); // String title)
ttop_en, // final boolean ttop_en, // false; // Remove tree tops from transparency weights
ttop_gb, // final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
ttop_min, // final double ttop_min, // 3.0; // Minimal tree top elevation
ttop_rel_lev, // final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
ttop_rel_rad, // final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
ttop_frac, // final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
ttop_rem_rad, // final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
ttop_no_border, // final boolean ttop_no_border,// true; // No maximums on the border allowed
dbg_title, // String title)
debugLevel); // final int debugLevel);
}
// num_iter = num_iters[num_iters.length-1]; // use last, may be 3-rd
//num_iters.length
}
// old
/*
if ("RESTORE".equals(par_path)) {
System.out.println("Reading fitted parameters from file");
if (save_par_files) {
String restore_dir = vegetationLMA.debug_path;
String debug_title = vegetationLMA.getParametersDebugTitle();
vegetationLMA.saveParametersFile(
restore_dir, // String dir,
debug_title, // String title, // no .par-tiff
null); // double [] vector)
continue;
}
}
*/
if ((show_final_result) && (debugLevel > -2)) { // 0)) {
String reconstructed_title = reference_scene+"-reconstructed-pre-step"+step_restore; //
vegetationLMA.showYfX(
......@@ -2324,7 +2272,16 @@ public class VegetationModel {
transparency_gb, // final double transparency_gb,
transparency_boost, // final double transparency_boost,
recalc_average, // final boolean recalc_average);
dbg_title); // String title)
ttop_en, // final boolean ttop_en, // false; // Remove tree tops from transparency weights
ttop_gb, // final double ttop_gb, // 1.0; // Elevation Gaussian blur sigma to detect tree tops
ttop_min, // final double ttop_min, // 3.0; // Minimal tree top elevation
ttop_rel_lev, // final double ttop_rel_lev, // 0.9; // Relative (to the top height) sample level
ttop_rel_rad, // final double ttop_rel_rad, // 0.25; // Relative (to the top height) sample ring radius
ttop_frac, // final double ttop_frac, // 0.5; // Minimal fraction of the ring pixels below sample level
ttop_rem_rad, // final double ttop_rem_rad, // 0.25; // Relative (to the top height) remove transparency radius
ttop_no_border, // final boolean ttop_no_border,// true; // No maximums on the border allowed
dbg_title, // final String title
debugLevel); // final int debugLevel);
}
if (save_par_files) {
......
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