Skip to content

I
imagej-elphel
Commit 81c431ad authored by Andrey Filippov's avatar Andrey Filippov
Browse files
Options

DATI and pattern matching

parent 1c62dfa5
Show whitespace changes
Inline Side-by-side
Showing with 1633 additions and 135 deletions
src/main/java/com/elphel/imagej/common/DoubleFHT.java
View file @ 81c431ad
...@@ -418,6 +418,53 @@ public class DoubleFHT { ...@@ -418,6 +418,53 @@ public class DoubleFHT {
return first; return first;
} }
/**
* Transform pattern once, and then use FD of the pattern for each overlapping image tile
* @param second pattern to be transform
* @return transformed pattern ( the original is also modified)
*/
public double [] transformPattern ( // new
double [] second ){ //null-OK
updateMaxN(second);
swapQuadrants(second);
if (!transform(second,false)) return null; // direct FHT
return second;
}
/**
* Phase correlate with pattern, that is transformed separately with transformPattern()
* @param first square array to be phase-correlated, modified but not with result!
* @param secondFD pattern already transformed to FD
* @param phaseCoeff 0 - regular correlation, 1.0 - pure phase correlation
* @param filter frequency filter or null
* @param first_save null or array to accommodate FD of the first array before updating it with correlation
* @return phase correlation, first still contains original transformed!
*/
public double [] phaseCorrelatePattern ( // new
double [] first,
double [] secondFD,
double phaseCoeff,
double [] filter, // high/low pass filtering
double [] first_save ){ //null-OK
if (first.length!=secondFD.length) {
IJ.showMessage("Error","Correlation arrays should be the same size");
return null;
}
updateMaxN(first);
swapQuadrants(first);
if (!transform(first,false)) return null; // direct FHT
if (first_save != null) System.arraycopy(first, 0, first_save, 0, first.length);
first= phaseMultiplyNorm(first, secondFD, phaseCoeff); // correlation, not convolution
if (filter!=null) multiplyByReal(first, filter);
transform(first,true) ; // inverse transform
swapQuadrants(first);
return first;
}
// //
......
src/main/java/com/elphel/imagej/correction/Eyesis_Correction.java
View file @ 81c431ad
...@@ -847,9 +847,12 @@ public class Eyesis_Correction implements PlugIn, ActionListener { ...@@ -847,9 +847,12 @@ public class Eyesis_Correction implements PlugIn, ActionListener {
addButton("Set pair", panelLWIRWorld, color_process); addButton("Set pair", panelLWIRWorld, color_process);
addButton("Warp pair", panelLWIRWorld, color_process); addButton("Warp pair", panelLWIRWorld, color_process);
addButton("Read Tiff", panelLWIRWorld, color_process); addButton("Read Tiff", panelLWIRWorld, color_process);
addButton("Set pair GPS", panelLWIRWorld, color_process);
addButton("Test video", panelLWIRWorld, color_process); addButton("Test video", panelLWIRWorld, color_process);
addButton("Deconvolve Slices", panelLWIRWorld, color_process); addButton("Ortho Pairs", panelLWIRWorld, color_process);
addButton("Mismatched resolutions", panelLWIRWorld, color_process);
addButton("Generate DATI", panelLWIRWorld, color_process);
addButton("Create mine", panelLWIRWorld, color_process);
addButton("Pattern correlate", panelLWIRWorld, color_process);
plugInFrame.add(panelLWIRWorld); plugInFrame.add(panelLWIRWorld);
} }
...@@ -5705,7 +5708,7 @@ public class Eyesis_Correction implements PlugIn, ActionListener { ...@@ -5705,7 +5708,7 @@ public class Eyesis_Correction implements PlugIn, ActionListener {
e.printStackTrace(); e.printStackTrace();
} }
// ComboMatch.testReadTiff(); // ComboMatch.testReadTiff();
} else if (label.equals("Set pair GPS")) { } else if (label.equals("Ortho Pairs")) {
DEBUG_LEVEL = MASTER_DEBUG_LEVEL; DEBUG_LEVEL = MASTER_DEBUG_LEVEL;
EYESIS_CORRECTIONS.setDebug(DEBUG_LEVEL); EYESIS_CORRECTIONS.setDebug(DEBUG_LEVEL);
if (GPU_TILE_PROCESSOR == null) { if (GPU_TILE_PROCESSOR == null) {
...@@ -5734,24 +5737,44 @@ public class Eyesis_Correction implements PlugIn, ActionListener { ...@@ -5734,24 +5737,44 @@ public class Eyesis_Correction implements PlugIn, ActionListener {
return; return;
} }
OrthoMap.testVideo(imp_sel); OrthoMap.testVideo(imp_sel);
} else if (label.equals("Deconvolve Slices")) { } else if (label.equals("Mismatched resolutions")) {
ImagePlus imp_sel = WindowManager.getCurrentImage(); ImagePlus imp_sel = WindowManager.getCurrentImage();
if (imp_sel == null) { if (imp_sel == null) {
IJ.showMessage("Error", "No images selected"); IJ.showMessage("Error", "No images selected");
return; return;
} }
OrthoMap.testDeconvolveSlices( OrthoMap.createMismatchedResolutionKernel(
imp_sel, // ImagePlus imp, imp_sel, // ImagePlus imp,
512, // int [] slices, 512, // int [] slices,
new int [] {1,2}, // int [] slices, deconvolve slice 2 with slice1 new int [] {1,2}, // int [] slices, deconvolve slice 2 with slice1
4, // int kernel_radius, 6, // 4, // int kernel_radius,
true, // boolean hor_sym, true, // boolean hor_sym,
true, // boolean vert_sym, true, // boolean vert_sym,
true, // false, // boolean all_sym, true, // false, // boolean all_sym,
DEBUG_LEVEL); // int debugLevel) { // >0 DEBUG_LEVEL); // int debugLevel) { // >0
} else if (label.equals("Generate DATI")) {
ImagePlus imp_sel = WindowManager.getCurrentImage();
if (imp_sel == null) {
IJ.showMessage("Error", "No images selected");
return;
}
OrthoMap.generateDATI(
imp_sel, // ImagePlus imp,
512, // int [] slices,
new int [] {1,2}, // int [] slices, deconvolve slice 2 with slice1
DEBUG_LEVEL); // int debugLevel) { // >0
} else if (label.equals("Create mine")) {
OrthoMap.testPatternGenerate();
} else if (label.equals("Pattern correlate")) {
ImagePlus imp_sel = WindowManager.getCurrentImage();
if (imp_sel == null) {
IJ.showMessage("Error", "No images selected");
return;
} }
OrthoMap.testPatternCorrelate(imp_sel);
} }
}
//
public boolean debugInitOneScene() { public boolean debugInitOneScene() {
DEBUG_LEVEL = MASTER_DEBUG_LEVEL; DEBUG_LEVEL = MASTER_DEBUG_LEVEL;
if (EYESIS_CORRECTIONS_AUX == null) { if (EYESIS_CORRECTIONS_AUX == null) {
......
src/main/java/com/elphel/imagej/orthomosaic/ComboMatch.java
View file @ 81c431ad
...@@ -82,8 +82,15 @@ public class ComboMatch { ...@@ -82,8 +82,15 @@ public class ComboMatch {
// String orthoMapsCollection_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_08_november.data"; // String orthoMapsCollection_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_08_november.data";
// String files_list_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_10_short.list"; // String files_list_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_10_short.list";
// String orthoMapsCollection_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_10_short.data"; // String orthoMapsCollection_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_10_short.data";
String files_list_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_19_sep13_25-50-75-100m.list"; // String files_list_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_19_sep13_25-50-75-100m.list";
String orthoMapsCollection_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_19_sep13_25-50-75-100m.data"; // String orthoMapsCollection_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_19_sep13_25-50-75-100m.data";
String [] files_lists_paths = {
"/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_19_sep13_25-50-75-100m",
"/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_nov3_50-75"};
String files_list_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_nov3_50-75.list";
String orthoMapsCollection_path = "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/maps_nov3_50-75.data";
//maps_nov3_50-75
//maps_19_sep13.list //maps_19_sep13.list
//maps_09_short.list //maps_09_short.list
//maps_08_november.list //maps_08_november.list
...@@ -129,26 +136,20 @@ public class ComboMatch { ...@@ -129,26 +136,20 @@ public class ComboMatch {
boolean restore_temp = true; boolean restore_temp = true;
double frac_remove = 0.15; double frac_remove = 0.15;
double metric_error = 0.05; // 0.02;// 2 cm double metric_error = 0.05; // 0.02;// 2 cm
boolean update_lla = false; // re-read file metadata
if (!use_marked_image) { if (!use_marked_image) {
process_correlation=false; // use already adjusted by defualt process_correlation=false; // use already adjusted by default
} }
GenericJTabbedDialog gd = new GenericJTabbedDialog("Set image pair",1200,900); GenericJTabbedDialog gd = new GenericJTabbedDialog("Set image pair",1200,900);
gd.addStringField ("Image list full path", files_list_path, 180, "Image list full path."); gd.addChoice ("Files list/data path (w/o extension):", files_lists_paths, files_lists_paths[files_lists_paths.length-1]);
gd.addStringField ("Maps collection save path", orthoMapsCollection_path, 180, "Save path for serialized map collection data.");
gd.addCheckbox ("Use saved maps collection", use_saved_collection, "If false - use files list."); gd.addCheckbox ("Use saved maps collection", use_saved_collection, "If false - use files list.");
gd.addCheckbox ("Save maps collection", save_collection, "Save maps collection to be able to restore."); gd.addCheckbox ("Save maps collection", save_collection, "Save maps collection to be able to restore.");
gd.addCheckbox ("Process correlations", process_correlation, "false to skip to just regenerate new save file."); gd.addCheckbox ("Process correlations", process_correlation, "false to skip to just regenerate new save file.");
gd.addCheckbox ("Update match if calculated", update_match, "Will update correlation match for a pair if found."); gd.addCheckbox ("Update match if calculated", update_match, "Will update correlation match for a pair if found.");
gd.addCheckbox ("Render match", render_match, "Render a pair of matched images."); gd.addCheckbox ("Render match", render_match, "Render a pair of matched images.");
/*
//
// for (int n = 0; n < image_paths_pre.length; n++) {
// gd.addStringField ("Image path "+n, image_paths_pre[n], 180, "Image "+n+" full path w/o ext");
// }
// gd.addStringField ("First image path", image_paths_pre[0], 180, "First image full path w/o ext");
// gd.addStringField ("Second image path", image_paths_pre[1], 180, "Second image full path w/o ext");
for (int n = 0; n < image_enuatr.length; n++) { for (int n = 0; n < image_enuatr.length; n++) {
gd.addMessage("image["+n+"] pose"); gd.addMessage("image["+n+"] pose");
gd.addNumericField("East", image_enuatr[n][0][0], 3,7,"m", "Move image "+n+" East."); gd.addNumericField("East", image_enuatr[n][0][0], 3,7,"m", "Move image "+n+" East.");
...@@ -158,6 +159,7 @@ public class ComboMatch { ...@@ -158,6 +159,7 @@ public class ComboMatch {
gd.addNumericField("Tilt", image_enuatr[n][1][1], 3,7,"deg", "Rotate image "+n+" around East."); gd.addNumericField("Tilt", image_enuatr[n][1][1], 3,7,"deg", "Rotate image "+n+" around East.");
gd.addNumericField("Roll", image_enuatr[n][1][2], 3,7,"deg", "Rotate image "+n+" around North."); gd.addNumericField("Roll", image_enuatr[n][1][2], 3,7,"deg", "Rotate image "+n+" around North.");
} }
*/
gd.addNumericField("Zoom level", zoom_lev, 0,4,"", gd.addNumericField("Zoom level", zoom_lev, 0,4,"",
"Zoom level: +1 - zoom in twice, -1 - zoom out twice"); "Zoom level: +1 - zoom in twice, -1 - zoom out twice");
gd.addNumericField("GPU image width", gpu_width, 0,4,"", gd.addNumericField("GPU image width", gpu_width, 0,4,"",
...@@ -167,23 +169,24 @@ public class ComboMatch { ...@@ -167,23 +169,24 @@ public class ComboMatch {
gd.addCheckbox ("Show transformation centers", show_centers, "Mark verticals from the UAS on the ground."); gd.addCheckbox ("Show transformation centers", show_centers, "Mark verticals from the UAS on the ground.");
gd.addCheckbox ("Show combo image map", show_combo_map, "Load and process altitude maps."); gd.addCheckbox ("Show combo image map", show_combo_map, "Load and process altitude maps.");
gd.addCheckbox ("Show altitude combo image", use_alt, "Load and process altitude maps."); gd.addCheckbox ("Show altitude combo image", use_alt, "Load and process altitude maps.");
// gd.addStringField ("First matching timestamp", gpu_spair[0], 20, "First GPU-matching timestamp with '_' for decimal point.");
// gd.addStringField ("Second matching timestamp", gpu_spair[1], 20, "First GPU-matching timestamp with '_' for decimal point.");
gd.addNumericField("Remove fraction of worst matches", frac_remove, 3,7,"", "When fitting scenes remove this fraction of worst match."); gd.addNumericField("Remove fraction of worst matches", frac_remove, 3,7,"", "When fitting scenes remove this fraction of worst match.");
gd.addNumericField("Maximal metric error", metric_error, 3,7,"m", "Maximal tolerable fitting error caused by elevation variations."); gd.addNumericField("Maximal metric error", metric_error, 3,7,"m", "Maximal tolerable fitting error caused by elevation variations.");
if (use_marked_image ) { if (use_marked_image ) {
gd.addCheckbox ("Use marked image data", true, "Use markes from the selected image"); gd.addCheckbox ("Use marked image data", true, "Use markes from the selected image");
} }
gd.addCheckbox ("Update files metadata", update_lla, "Re-read files metadata (ifd it was modified)");
gd.showDialog(); gd.showDialog();
if (gd.wasCanceled()) return false; if (gd.wasCanceled()) return false;
files_list_path = gd.getNextString(); int choice_index = gd.getNextChoiceIndex();
files_list_path = files_lists_paths[choice_index]+".list";
orthoMapsCollection_path = gd.getNextString(); orthoMapsCollection_path =files_lists_paths[choice_index]+".data";
use_saved_collection = gd.getNextBoolean(); use_saved_collection = gd.getNextBoolean();
save_collection = gd.getNextBoolean(); save_collection = gd.getNextBoolean();
process_correlation= gd.getNextBoolean(); process_correlation= gd.getNextBoolean();
update_match= gd.getNextBoolean(); update_match= gd.getNextBoolean();
render_match= gd.getNextBoolean(); render_match= gd.getNextBoolean();
/*
for (int n = 0; n < image_enuatr.length; n++) { for (int n = 0; n < image_enuatr.length; n++) {
image_enuatr[n][0][0] = gd.getNextNumber(); image_enuatr[n][0][0] = gd.getNextNumber();
image_enuatr[n][0][1] = gd.getNextNumber(); image_enuatr[n][0][1] = gd.getNextNumber();
...@@ -192,6 +195,7 @@ public class ComboMatch { ...@@ -192,6 +195,7 @@ public class ComboMatch {
image_enuatr[n][1][1] = gd.getNextNumber(); image_enuatr[n][1][1] = gd.getNextNumber();
image_enuatr[n][1][2] = gd.getNextNumber(); image_enuatr[n][1][2] = gd.getNextNumber();
} }
*/
zoom_lev = (int) gd.getNextNumber(); zoom_lev = (int) gd.getNextNumber();
gpu_width = (int) gd.getNextNumber(); gpu_width = (int) gd.getNextNumber();
gpu_height = (int) gd.getNextNumber(); gpu_height = (int) gd.getNextNumber();
...@@ -209,6 +213,7 @@ public class ComboMatch { ...@@ -209,6 +213,7 @@ public class ComboMatch {
if (use_marked_image ) { // will only be used if found and asked if (use_marked_image ) { // will only be used if found and asked
use_marked_image= gd.getNextBoolean(); use_marked_image= gd.getNextBoolean();
} }
update_lla= gd.getNextBoolean();
OrthoMapsCollection maps_collection=null; OrthoMapsCollection maps_collection=null;
if (use_saved_collection) { if (use_saved_collection) {
try { try {
...@@ -245,8 +250,11 @@ public class ComboMatch { ...@@ -245,8 +250,11 @@ public class ComboMatch {
//getTemperature() //getTemperature()
// get all temperatures // get all temperatures
maps_collection. getAllTemperatures(); maps_collection.getAllTemperatures();
if (update_lla) {
System.out.println("Updating map files metadata");
maps_collection.updateLLa();
}
// which pair to compare // which pair to compare
// String [] gpu_spair = {names[gpu_ipair[0]],names[gpu_ipair[1]]}; // String [] gpu_spair = {names[gpu_ipair[0]],names[gpu_ipair[1]]};
...@@ -310,7 +318,7 @@ public class ComboMatch { ...@@ -310,7 +318,7 @@ public class ComboMatch {
int max_zoom_lev = Math.max( int max_zoom_lev = Math.max(
maps_collection.ortho_maps[gpu_pair[0]].getOriginalZoomLevel(), maps_collection.ortho_maps[gpu_pair[0]].getOriginalZoomLevel(),
maps_collection.ortho_maps[gpu_pair[1]].getOriginalZoomLevel()); maps_collection.ortho_maps[gpu_pair[1]].getOriginalZoomLevel());
int initial_zoom = max_zoom_lev - 3; // another algorithm? int initial_zoom = max_zoom_lev - 4; // another algorithm?
System.out.println("Setting up GPU"); System.out.println("Setting up GPU");
......
src/main/java/com/elphel/imagej/orthomosaic/OrthoMap.java
View file @ 81c431ad
...@@ -2,6 +2,7 @@ package com.elphel.imagej.orthomosaic; ...@@ -2,6 +2,7 @@ package com.elphel.imagej.orthomosaic;
import java.awt.Color; import java.awt.Color;
import java.awt.Font; import java.awt.Font;
import java.awt.Point;
import java.awt.Rectangle; import java.awt.Rectangle;
import java.io.File; import java.io.File;
import java.io.IOException; import java.io.IOException;
...@@ -45,6 +46,7 @@ import ij.gui.Roi; ...@@ -45,6 +46,7 @@ import ij.gui.Roi;
import ij.plugin.filter.AVI_Writer; import ij.plugin.filter.AVI_Writer;
import ij.plugin.filter.GaussianBlur; import ij.plugin.filter.GaussianBlur;
import ij.process.ColorProcessor; import ij.process.ColorProcessor;
import ij.process.FloatPolygon;
import ij.process.FloatProcessor; import ij.process.FloatProcessor;
import ij.process.ImageConverter; import ij.process.ImageConverter;
import ij.process.ImageProcessor; import ij.process.ImageProcessor;
...@@ -57,13 +59,22 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -57,13 +59,22 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
public static int gpu_height = 4096; public static int gpu_height = 4096;
public static final String [] KEY_DIRS= {"rootDirectory", // from EyesisCorrectionParameters public static final String [] KEY_DIRS= {"rootDirectory", // from EyesisCorrectionParameters
"sourceDirectory","linkedModels","videoDirectory","x3dDirectory","resultsDirectory"}; "sourceDirectory","linkedModels","videoDirectory","x3dDirectory","resultsDirectory"};
public static final String [] kernel_paths = {
"--- no convolution (same altitude) ---",
"/media/elphel/NVME/lwir16-proc/ortho_videos/kernel_25_50.tiff",
"/media/elphel/NVME/lwir16-proc/ortho_videos/kernel_50_75.tiff",
"/media/elphel/NVME/lwir16-proc/ortho_videos/kernel_50_100.tiff"
};
public static final String ALT_SUFFIX = "-ALT"; public static final String ALT_SUFFIX = "-ALT";
// public transient final String name; // timestamp // public transient final String name; // timestamp
// public transient final double ts; // public transient final double ts;
public transient String name; // timestamp public transient String name; // timestamp
public transient double ts; public transient double ts;
public transient String path; // full path to the model directory (including /vXX?) public transient String path; // full path to the model directory (including /vXX?)
public transient String scenes_path; // full path to the model directory (including /vXX?) // public transient String scenes_path; // full path to the model directory (including /vXX?)
public double [] lla; // lat/long/alt public double [] lla; // lat/long/alt
public LocalDateTime dt; public LocalDateTime dt;
// affine convert (input) rectified coordinates (meters) relative to vert_meters to source image // affine convert (input) rectified coordinates (meters) relative to vert_meters to source image
...@@ -82,11 +93,12 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -82,11 +93,12 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
public transient double averageImagePixel = Double.NaN; // average image pixel value (to combine with raw) public transient double averageImagePixel = Double.NaN; // average image pixel value (to combine with raw)
transient HashMap <Integer, FloatImageData> images; transient HashMap <Integer, FloatImageData> images;
HashMap <String, PairwiseOrthoMatch> pairwise_matches; HashMap <String, PairwiseOrthoMatch> pairwise_matches;
public transient SensorTemperatureData[] temp_data;
private void writeObject(ObjectOutputStream oos) throws IOException { private void writeObject(ObjectOutputStream oos) throws IOException {
oos.defaultWriteObject(); oos.defaultWriteObject();
oos.writeObject(path); oos.writeObject(path);
oos.writeObject(scenes_path); // oos.writeObject(scenes_path);
// lla is not transient // lla is not transient
// dt is not transient // dt is not transient
// affine is not transient // affine is not transient
...@@ -101,6 +113,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -101,6 +113,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
// need_extra_zoom is not transient // need_extra_zoom is not transient
// images is not saved // images is not saved
// pairwise_matches is not transient // pairwise_matches is not transient
oos.writeObject(temp_data); // temporary, while transient
} }
private void readObject(ObjectInputStream ois) throws ClassNotFoundException, IOException { private void readObject(ObjectInputStream ois) throws ClassNotFoundException, IOException {
...@@ -108,7 +121,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -108,7 +121,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
path = (String) ois.readObject(); path = (String) ois.readObject();
name = getNameFromPath(path); name = getNameFromPath(path);
ts= Double.parseDouble(name.replace("_", ".")); ts= Double.parseDouble(name.replace("_", "."));
scenes_path= (String) ois.readObject(); // scenes_path= (String) ois.readObject();
// lla is not transient // lla is not transient
// dt is not transient // dt is not transient
// affine is not transient // affine is not transient
...@@ -118,6 +131,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -118,6 +131,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
// orig_height is not transient // orig_height is not transient
// orig_image was not saved // orig_image was not saved
// alt_image was not saved // alt_image was not saved
temp_data = (SensorTemperatureData[]) ois.readObject();
images = new HashMap <Integer, FloatImageData>(); // field images was not saved images = new HashMap <Integer, FloatImageData>(); // field images was not saved
averageImagePixel = Double.NaN; // average image pixel value (to combine with raw) averageImagePixel = Double.NaN; // average image pixel value (to combine with raw)
...@@ -169,6 +183,26 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -169,6 +183,26 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
return path.substring(p1+1, p2); return path.substring(p1+1, p2);
} }
public static String getModelPathFromPath(String path) {
int p1 = path.lastIndexOf(Prefs.getFileSeparator());
if (p1 < 0) return null;
return path.substring(0,p1);
}
public String getScenesPath() {
return getModelPathFromPath(this.path)+"/jp4";
}
public void updateLLA() {
Properties imp_prop = null;
try {
imp_prop = ElphelTiffReader.getTiffMeta(path);
} catch (IOException e) {
// TODO Auto-generated catch block
e.printStackTrace();
}
lla = ElphelTiffReader.getLLA(imp_prop);
}
/** /**
* Gets average pixel value of all sensors of the reference scene * Gets average pixel value of all sensors of the reference scene
* @return * @return
...@@ -176,7 +210,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -176,7 +210,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
public double getTemperature() { public double getTemperature() {
if (Double.isNaN(averageRawPixel)) { if (Double.isNaN(averageRawPixel)) {
final String TIFF_EXT=".tiff"; final String TIFF_EXT=".tiff";
final File [] raw_files = (new File(scenes_path)).listFiles(); final File [] raw_files = (new File(getScenesPath())).listFiles();
final Thread[] threads = ImageDtt.newThreadArray(); final Thread[] threads = ImageDtt.newThreadArray();
final AtomicInteger ai = new AtomicInteger(0); final AtomicInteger ai = new AtomicInteger(0);
final AtomicInteger ati = new AtomicInteger(0); final AtomicInteger ati = new AtomicInteger(0);
...@@ -222,12 +256,74 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -222,12 +256,74 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
return averageRawPixel; return averageRawPixel;
} }
public double getTemperatureAndTelemetry() {
if (this.temp_data == null) {
final String TIFF_EXT=".tiff";
final File [] raw_files = (new File(getScenesPath())).listFiles();
final Thread[] threads = ImageDtt.newThreadArray();
final AtomicInteger ai = new AtomicInteger(0);
final AtomicInteger ati = new AtomicInteger(0);
final double [] avg_arr = new double [threads.length];
final double [] npix_arr = new double [threads.length];
final SensorTemperatureData [] temp_data = new SensorTemperatureData[raw_files.length];
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
int thread_num = ati.getAndIncrement();
ImagejJp4Tiff imagejJp4Tiff = new ImagejJp4Tiff();
for (int iFile = ai.getAndIncrement(); iFile < raw_files.length; iFile = ai.getAndIncrement()) {
if (raw_files[iFile].getName().endsWith(TIFF_EXT)) {
ImagePlus imp = null;
String spath = raw_files[iFile].toString();
try {
imp= imagejJp4Tiff.readTiffJp4(spath);
} catch (IOException e) {
System.out.println("getImagesMultithreaded IOException " + spath);
} catch (FormatException e) {
System.out.println("getImagesMultithreaded FormatException " + spath);
}
if (imp != null) {
int channel_sep = spath.lastIndexOf("_");
int channel_dot = spath.lastIndexOf(".");
int chn = Integer.parseInt(spath.substring(channel_sep+1, channel_dot));
temp_data[chn] = new SensorTemperatureData();
Properties properties = imp.getProperties();
temp_data[chn].TLM_FPA_KELV = Double.parseDouble(properties.getProperty("TLM_FPA_KELV"));
temp_data[chn].TLM_FFC_KELV = Double.parseDouble(properties.getProperty("TLM_FFC_KELV")); temp_data[chn].TLM_FPA_KELV = Double.parseDouble(properties.getProperty("TLM_FPA_KELV"));
temp_data[chn].TLM_CORE_TEMP = Double.parseDouble(properties.getProperty("TLM_CORE_TEMP"));
temp_data[chn].TLM_FRAME = Integer.parseInt(properties.getProperty("TLM_FRAME"));
temp_data[chn].TLM_FRAME_FFC = Integer.parseInt(properties.getProperty("TLM_FRAME_FFC"));
npix_arr[thread_num] = imp.getWidth() * imp.getHeight();
float [] fpixels = (float[]) imp.getProcessor().getPixels();
for (int i = 0; i < fpixels.length; i++) {
avg_arr[thread_num] += fpixels[i];
}
avg_arr[thread_num] /= npix_arr[thread_num];
temp_data[chn].averagePixelValue = avg_arr[thread_num];
}
}
}
}
};
}
ImageDtt.startAndJoin(threads);
double avg=0, num=0;
for (int i = 0; i < avg_arr.length; i++) {
avg+=avg_arr[i]*npix_arr[i];
num+=npix_arr[i];
}
averageRawPixel = avg/num;
this.temp_data = temp_data;
}
return averageRawPixel;
}
public OrthoMap ( public OrthoMap (
String path, String path) {
String scenes_path) { // String scenes_path) {
this.path = path; this.path = path;
this.scenes_path = scenes_path; // this.scenes_path = scenes_path; // will not be used
name = getNameFromPath(path); name = getNameFromPath(path);
ts= Double.parseDouble(name.replace("_", ".")); ts= Double.parseDouble(name.replace("_", "."));
Properties imp_prop = null; Properties imp_prop = null;
...@@ -762,12 +858,502 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -762,12 +858,502 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
ImageDtt.startAndJoin(threads); ImageDtt.startAndJoin(threads);
return padded_gpu; return padded_gpu;
} }
public static void generateDATI(
ImagePlus imp,
final int fft_size,
int [] slices,
int debugLevel) { // >0
double [] kernel = getConvolutionKernel();
debugLevel = 1;
boolean flat_high_res = false;
boolean temp_high_res = false; // false - use convolved with kernel
double neib_radius = 2.0; // 1.5;
double neib_radius_other= 1.5; // 1.1;
double max_heib_diff = 20; // 30; // 50; // not much filtering
double max_neib_diff_other= 30;// 50; // 60;
boolean four_sides = true;//false;
int num_temp = 256; // 1024; // height
int num_dati = 256; // 1024; // width
double frac_hi_temp = 0.005;
double frac_lo_temp = 0.0001;
double frac_hi_dati = 0.005;
double frac_lo_dati = 0.005;
boolean invert_y = true; // high temperature - up
int mine_number = 2; // 100; // all which of the first min_markers are mines (can be > length)
double falsepos_rel_radius= 0.05; // 0.04; // relative (to selected ranges) radius of false positives
final boolean generate_dati_histogram = true;
GenericJTabbedDialog gd = new GenericJTabbedDialog("Set image pair",1200,600);
gd.addCheckbox ("Use high-res for flat points", flat_high_res, "If false - use low-res convolved version.");
gd.addCheckbox ("Use high-res for histogram", temp_high_res, "If false - - use low-res convolved version.");
gd.addMessage("--- Finding locally flat points to use for the 2D histogram ---");
gd.addCheckbox ("4-sides flat", four_sides, "Require local min-max in both horizontal and vertical direction (false - one pair is enough).");
gd.addNumericField("Flatness check in first image", neib_radius, 3,7,"pix", "Low/high res as selected above.");
gd.addNumericField("Flatness check in second image", neib_radius_other, 3,7,"pix", "In second image.");
gd.addNumericField("Maximal difference in first image", max_heib_diff, 3,7,"counts", "Maximal difference from the center pixel within the specified radius.");
gd.addNumericField("Maximal difference in first image", max_neib_diff_other, 3,7,"counts", "Maximal difference from the center pixel within the specified radius.");
gd.addMessage("--- Generating 2D histogram ---");
gd.addNumericField("Histogram height (temperature)", num_temp, 0,4,"pix", "Number of temperature bins,");
gd.addNumericField("Histogram width (DATI)", num_dati, 0,4,"pix", "Number of DATI bins,");
gd.addNumericField("Discard high temperature samples", frac_hi_temp, 4,8,"", "Fraction of all samples.");
gd.addNumericField("Discard low temperature samples", frac_lo_temp, 4,8,"", "Fraction of all samples.");
gd.addNumericField("Discard high DATI samples", frac_hi_dati, 4,8,"", "Fraction of all samples.");
gd.addNumericField("Discard low DATI samples", frac_lo_dati, 4,8,"", "Fraction of all samples.");
gd.addCheckbox ("High temperature - up", invert_y, "Plot high-temperature up (y-axis). High DATI is always to the right (x-axis).");
gd.addNumericField("Number of marked mines", mine_number, 0,4,"", "First marks are considered mines, other - just samples to see in the histogram.");
gd.addNumericField("Cluster radius (relative)", falsepos_rel_radius, 3,7,"","Relative to temperature and DATI plotted ranges.");
gd.addMessage("Histogram may be manually blurred with sigma=2.0" );
gd.showDialog();
if (gd.wasCanceled()) return;
flat_high_res= gd.getNextBoolean();
temp_high_res= gd.getNextBoolean();
four_sides= gd.getNextBoolean();
neib_radius= gd.getNextNumber();
neib_radius_other= gd.getNextNumber();
max_heib_diff= gd.getNextNumber();
max_neib_diff_other= gd.getNextNumber();
num_temp= (int) gd.getNextNumber();
num_dati= (int) gd.getNextNumber();
frac_hi_temp= gd.getNextNumber();
frac_lo_temp= gd.getNextNumber();
frac_hi_dati= gd.getNextNumber();
frac_lo_dati= gd.getNextNumber();
invert_y= gd.getNextBoolean();
mine_number= (int) gd.getNextNumber();
falsepos_rel_radius= gd.getNextNumber();
final int slice = flat_high_res? 0: 1; // 1; // convolved
final int temp_index = temp_high_res? 0: 1;
ImageStack stack = imp.getStack();
final int width = stack.getWidth();
final int height = stack.getHeight();
// Roi roi= imp.getRoi();
// boolean good_ROI = false;
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
final double [][] full_img = new double [4][width*height];
final int [] out_slices = {0,2,1,3}; // {high-res, low-res, convoled_high_res,dati}
for (int n = 0; n < 2; n++) {
final int fn = out_slices[n];
final float [] fpixels = (float[]) stack.getPixels(slices[n]);
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int ipix = ai.getAndIncrement(); ipix < full_img[fn].length; ipix = ai.getAndIncrement()) {
full_img[fn][ipix] = fpixels[ipix];
}
}
};
}
ImageDtt.startAndJoin(threads);
}
if (kernel != null) {
full_img[out_slices[2]] = convolveWithKernel(
full_img[out_slices[0]], // final double [] data,
kernel, // final double [] kernel,
width); // final int width)
} else {
full_img[out_slices[2]] = full_img[out_slices[0]].clone();
}
// calculate DATI
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int ipix = ai.getAndIncrement(); ipix < full_img[0].length; ipix = ai.getAndIncrement()) {
full_img[out_slices[3]][ipix] = full_img[out_slices[1]][ipix]-full_img[out_slices[2]][ipix];
}
}
};
}
ImageDtt.startAndJoin(threads);
if (!generate_dati_histogram) { // show just that, otherwise add more layers
if (debugLevel>0) {
String [] test_titles= {"high_res","convolved", "low_res", "DATI"};
ShowDoubleFloatArrays.showArrays(
full_img,
width,
height,
true,
removeKnownExtension(imp.getTitle())+"-DATI",
test_titles);
}
} else {
int [] nonNaNs = {2}; // low_res
int slice_other = 2; // -1 to skip testing for flatness
boolean [] flat = extractFlatPixels(
full_img, // final double [][] data,
width, // final int width,
slice, // final int slice,
nonNaNs, // final int [] nonNaNs,
slice_other, // final int slice_other, // <0 if not controlled
four_sides, // final boolean four_sides, // false - only opposite
neib_radius, // final double neib_radius,
max_heib_diff, // final double max_heib_diff,
neib_radius_other, // final double neib_radius_other,
max_neib_diff_other, // final double max_neib_diff_other,
debugLevel); // final int debugLevel)
double [] data_temp = full_img[temp_index]; // [1]; // [0]; // may be one of the others : [1], [2]
if (debugLevel >-4) {
String sres = (new String[] {"high res","low res (convolved with kernel)"})[temp_index];
System.out.println("Using reference layer "+temp_index+" ("+sres+").");
}
double [] data_dati = full_img[3]; // may be one of the others : [1], [2]
// get mine markers
double [][] mine_markers = null;
PointRoi pRoi = (PointRoi) imp.getRoi();
if (pRoi != null) {
FloatPolygon fp = pRoi.getContainedFloatPoints();
mine_markers = new double [fp.npoints][2];
for (int np = 0; np< mine_markers.length; np++) {
mine_markers[np][0] = fp.xpoints[np];
mine_markers[np][1] = fp.ypoints[np];
}
}
final ArrayList<Integer> falsepos_list=new ArrayList<Integer>();
ImagePlus imp_dati = mapFlatPixels(
data_temp, // final double [] temp, // may try different
data_dati, // final double [] dati,
flat, // final boolean[] flat,
width, // final int width,
num_temp, // final int num_temp, // height
num_dati, // final int num_dati, // width
invert_y, // final boolean invert_y,
frac_hi_temp, // final double frac_hi_temp,
frac_lo_temp, // final double frac_lo_temp,
frac_hi_dati, // final double frac_hi_dati,
frac_lo_dati, // final double frac_lo_dati,
mine_markers, // final double [][] mine_markers, // may be null of the original image
mine_number, // final int mine_number, // which of the first min_markers are mines (can be > length)
falsepos_rel_radius, // final double falsepos_rel_radius, // relative (to selected ranges) radius of false positives
falsepos_list, // final ArrayList<Integer> falsepos_list, // initialized list of false positive indices in temp[],dati[], flat
debugLevel); // final int debugLevel)
String title = removeKnownExtension(imp.getTitle())+"-"+imp_dati.getTitle();
imp_dati.setTitle(title);
imp_dati.show();
if ((falsepos_list != null) && (falsepos_list.size() > 0)) {
System.out.println("generateDATI(): Found "+falsepos_list.size()+" potential false+real positives, marked in the image");
}
if (debugLevel>-4) {
String [] test_titles= {"high_res","convolved", "low_res", "DATI"};
String impfp_title = removeKnownExtension(imp.getTitle())+"-FALSE";
ImageStack stack_impfp = ShowDoubleFloatArrays.makeStack(
full_img,
width,
height,
test_titles,
false);
ImagePlus impfp = new ImagePlus(impfp_title, stack_impfp);
PointRoi roifp = new PointRoi();
for (int ipix:falsepos_list) {
roifp.addPoint(ipix % width, ipix / width, 4);
}
roifp.setOptions("label, circle");
impfp.setRoi(roifp);
impfp.show();
}
}
/*
*/
System.out.println("generateDATI() Done");
}
/**
* Read one of the convolution kernels (or nun if no convolution is needed)
* @return square kernel as a 1-d double array or null
*/
public static double [] getConvolutionKernel() {
String kernel_path = null;
GenericJTabbedDialog gds = new GenericJTabbedDialog("Select kernel path",1200,400);
gds.addChoice("Kernel path:", kernel_paths, kernel_paths[kernel_paths.length-1]);
gds.showDialog();
if (gds.wasCanceled()) return null;
int kernel_index = gds.getNextChoiceIndex();
double [] kernel = null;
if (kernel_index > 0) {
kernel_path= kernel_paths[kernel_index];
ImagePlus imp_kernel = new ImagePlus(kernel_path);
if (imp_kernel.getWidth() == 0) {
System.out.println("generateDATI(): precomputed kernel "+kernel_path+
" is not found");
return null;
}
float [] kernel_pixels = (float[]) imp_kernel.getProcessor().getPixels();
kernel = new double[kernel_pixels.length];
for (int i = 0; i < kernel.length; i++) {
kernel[i] = kernel_pixels[i];
}
}
return kernel;
}
public static ImagePlus mapFlatPixels(
final double [] temp, // may try different
final double [] dati,
final boolean[] flat,
final int width,
final int num_temp, // height
final int num_dati, // width
final boolean invert_y,
final double frac_hi_temp,
final double frac_lo_temp,
final double frac_hi_dati,
final double frac_lo_dati,
final double [][] mine_markers, // may be null of the original image
final int mine_number_in, // which of the first min_markers are mines (can be > length)
final double falsepos_rel_radius, // relative (to selected ranges) radius of false positives
final ArrayList<Integer> falsepos_list, // initialized list of false positive indices in temp[],dati[], flat
final int debugLevel) {
final int mine_number = (mine_markers != null)? Math.min(mine_number_in, mine_markers.length):0;
final int num_bins = 1000;
final double [][] min_max = new double [2][];
min_max[0] = getMinMax(
temp, // final double [] data,
flat, // final boolean[] mask, // may be null
null, // final double [] weights, // may be null
frac_hi_temp, // final double frac_hi,
frac_lo_temp, // final double frac_lo,
num_bins, // int num_bins,
debugLevel); // final int debug_level)
min_max[1] = getMinMax(
dati, // final double [] data,
flat, // final boolean[] mask, // may be null
null, // final double [] weights, // may be null
frac_hi_temp, // final double frac_hi,
frac_lo_temp, // final double frac_lo,
num_bins, // int num_bins,
debugLevel); // final int debug_level)
if (debugLevel > -4) {
System.out.println("mapFlatPixels(): Minimal temperature: "+min_max[0][0]+", maximal: "+min_max[0][1]);
System.out.println(" Minimal difference: "+min_max[1][0]+", maximal: "+min_max[1][1]);
}
final double [] scales = {
num_temp/( min_max[0][1]- min_max[0][0]),
num_dati/( min_max[1][1]- min_max[1][0])
};
final double [] map = new double [num_dati*num_temp];
for (int i = 0; i < flat.length; i++) if (flat[i]) {
// TODO: average dati and temp over circles?
// TODO: use high-res (raw 50m) for temp?
// TODO: add markers
// TODO: read temperatures
int y = (int) ((temp[i] - min_max[0][0]) * scales[0]);
int x = (int) ((dati[i] - min_max[1][0]) * scales[1]);
if (invert_y) y = num_temp - y -1;
if ((y >= 0) && (y < num_temp) && (x >=0) && (x <= num_dati)) {
map[y*num_dati + x] += 1.0;
}
}
String [] map_names = {"map"}; // maybe add more later
String map_title = String.format("DATI_map_t%f06.1:%f06.1_d%f06.1:%f06.1",
min_max[0][0],min_max[0][1],min_max[1][0],min_max[1][1]);
ImageStack stack = ShowDoubleFloatArrays.makeStack(
new double [][] {map},
num_dati,
num_temp,
map_names,
false);
ImagePlus imp = new ImagePlus(map_title, stack);
final double [][] pix_mine_centers = new double [mine_number][2];
if (mine_markers != null) {
int slice_map = 1;
PointRoi roi = new PointRoi();
for (int i = 0; i < mine_markers.length; i++) {
int indx = ((int) Math.round(mine_markers[i][1])) * width + ((int) Math.round(mine_markers[i][0]));
int x = (int) ((dati[indx] - min_max[1][0]) * scales[1]);
int y = (int) ((temp[indx] - min_max[0][0]) * scales[0]);
if (i < mine_number) {
pix_mine_centers[i][0] = x; // before limit
pix_mine_centers[i][1] = y; // before invert and limit
}
if (invert_y) y = num_temp - y -1;
if (y < 0) y = 0;
else if (y >= num_temp) y= num_temp - 1;
if (x < 0) x = 0;
else if (x >= num_dati) x= num_dati - 1;
roi.addPoint(x,y, slice_map);
}
roi.setOptions("label");
imp.setRoi(roi);
}
if (falsepos_list != null) {
final double pix_radius_temp = falsepos_rel_radius * (min_max[0][1]-min_max[0][0]) * scales[0];
final double pix_radius_dati = falsepos_rel_radius * (min_max[1][1]-min_max[1][0]) * scales[1];
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 iPix = ai.getAndIncrement(); iPix < flat.length; iPix = ai.getAndIncrement()) if (flat[iPix]){
int x = (int) ((dati[iPix] - min_max[1][0]) * scales[1]);
int y = (int) ((temp[iPix] - min_max[0][0]) * scales[0]);
for (int imine = 0; imine < mine_number; imine++) {
double dx = (x - pix_mine_centers[imine][0])/pix_radius_dati;
double dy = (y - pix_mine_centers[imine][1])/pix_radius_temp;
if ((dx*dx+dy*dy) < 1.0) {
falsepos_list.add(iPix);
break;
}
}
}
}
};
}
ImageDtt.startAndJoin(threads);
}
return imp;
}
public static double [] getMinMax(
final double [] data,
final boolean[] mask, // may be null
final double [] weights, // may be null
final double frac_hi,
final double frac_lo,
final int num_bins,
final int debug_level) {
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
final AtomicInteger ati = new AtomicInteger(0);
final double [][] min_max_thread = new double [threads.length][2];
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
int thread_num = ati.getAndIncrement();
min_max_thread[thread_num][0] = Double.POSITIVE_INFINITY;
min_max_thread[thread_num][1] = Double.NEGATIVE_INFINITY;
for (int iPix = ai.getAndIncrement(); iPix < data.length; iPix = ai.getAndIncrement()) {
if (!Double.isNaN(data[iPix]) && ((mask==null) || (mask[iPix]))){
if (data[iPix] < min_max_thread[thread_num][0]) {
min_max_thread[thread_num][0] = data[iPix];
}
if (data[iPix] > min_max_thread[thread_num][1]) {
min_max_thread[thread_num][1] = data[iPix];
}
}
}
}
};
}
ImageDtt.startAndJoin(threads);
final double [] min_max = {Double.POSITIVE_INFINITY, Double.NEGATIVE_INFINITY};
for (int i = 0; i < min_max_thread.length; i++) {
if (min_max_thread[i][0] < min_max[0]) {
min_max[0] = min_max_thread[i][0];
}
if (min_max_thread[i][1] > min_max[1]) {
min_max[1] = min_max_thread[i][1];
}
}
if ((frac_hi==0) && (frac_lo == 0)) {
return min_max; // no need for the histograms
}
final double [][] histograms = new double [threads.length][num_bins];
final double scale = num_bins / (min_max[1]-min_max[0]);
ai.set(0);
ati.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
int thread_num = ati.getAndIncrement();
for (int iPix = ai.getAndIncrement(); iPix < data.length; iPix = ai.getAndIncrement()) {
if (!Double.isNaN(data[iPix]) && ((mask==null) || (mask[iPix]))){
int bin = (int) Math.round((data[iPix] - min_max[0])*scale);
if (bin < 0) bin=0;
else if (bin >= num_bins) bin=num_bins-1;
double w = (weights == null)? 1.0: weights[iPix];
histograms[thread_num][bin] += w;
}
}
}
};
}
ImageDtt.startAndJoin(threads);
ai.set(0);
final double [] histogram = new double [num_bins];
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
for (int bin = ai.getAndIncrement(); bin < num_bins; bin = ai.getAndIncrement()) {
for (int i = 0; i < histograms.length; i++) {
histogram[bin]+= histograms[i][bin];
}
}
}
};
}
ImageDtt.startAndJoin(threads);
double sw = 0;
for (int bin = 0; bin < num_bins; bin++) {
sw += histogram[bin];
}
double trlow = sw * frac_lo;
double trhigh = sw * frac_hi;
double sh = 0, shp = 0;
double threshold_high = min_max[1];
if (frac_hi > 0) {
for (int bin = num_bins-1; bin>=0; bin--) {
shp = sh;
sh += histogram[bin];
if (sh > trhigh) {
double r = (sh-trhigh)/(sh-shp);
threshold_high = min_max[0] + (bin + r)/scale;
break;
}
}
}
sh = 0;
shp = 0;
double threshold_low = min_max[0];
if (frac_lo > 0) {
for (int bin = 0; bin < num_bins; bin++) {
shp = sh;
sh += histogram[bin];
if (sh > trlow) {
double r = (trlow-shp)/(sh-shp);
threshold_low = min_max[0] + (bin + r)/scale;
break;
}
}
}
return new double[] {threshold_low,threshold_high};
}
/** /**
* Extracts ROI, rounds to fft_size and deconvolves 2-nd slice (lower resolution) with * Extracts ROI, rounds to fft_size and deconvolves 2-nd slice (lower resolution) with
* the first one (higher resolution) * the first one (higher resolution)
* @param imp * @param imp
*/ */
public static void testDeconvolveSlices( public static void createMismatchedResolutionKernel(
ImagePlus imp, ImagePlus imp,
final int fft_size, final int fft_size,
int [] slices, int [] slices,
...@@ -776,10 +1362,8 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -776,10 +1362,8 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
boolean vert_sym, boolean vert_sym,
boolean all_sym, boolean all_sym,
int debugLevel) { // >0 int debugLevel) { // >0
ImageStack stack = imp.getStack(); ImageStack stack = imp.getStack();
final int width = stack.getWidth(); final int width = stack.getWidth();
final int height = stack.getHeight();
Roi roi= imp.getRoi(); Roi roi= imp.getRoi();
boolean good_ROI = false; boolean good_ROI = false;
Rectangle rroi = null; Rectangle rroi = null;
...@@ -795,28 +1379,10 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -795,28 +1379,10 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
if (slices == null) { if (slices == null) {
slices = new int [] {1,2}; slices = new int [] {1,2};
} }
// String kernel_path = "/media/elphel/NVME/lwir16-proc/ortho_videos/kernel_50_75.tiff";
String kernel_path = "/media/elphel/NVME/lwir16-proc/ortho_videos/kernel_50_100.tiff";
double [] kernel = null;
if (!good_ROI) { if (!good_ROI) {
ImagePlus imp_kernel = new ImagePlus(kernel_path);
if (imp_kernel.getWidth() == 0) {
System.out.println("testDeconvolveSlices(): precomputed kernel "+kernel_path+
" is not found, and for calculation");
System.out.println("it needs rectangular selection of "+fft_size+"x"+fft_size); System.out.println("it needs rectangular selection of "+fft_size+"x"+fft_size);
return; return;
} else {
float [] kernel_pixels = (float[]) imp_kernel.getProcessor().getPixels();
kernel = new double[kernel_pixels.length];
for (int i = 0; i < kernel.length; i++) {
kernel[i] = kernel_pixels[i];
}
}
} }
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
if (kernel==null) {
roi= imp.getRoi(); // retry roi roi= imp.getRoi(); // retry roi
rroi=roi.getBounds(); rroi=roi.getBounds();
if (rroi.width != fft_size) { if (rroi.width != fft_size) {
...@@ -826,10 +1392,9 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -826,10 +1392,9 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
rroi.height = fft_size; rroi.height = fft_size;
} }
final double [][] dpixels = new double [2][fft_size*fft_size]; final double [][] dpixels = new double [2][fft_size*fft_size];
// ImageStack stack = imp.getStack();
// final int width = stack.getWidth();
// final int height = stack.getHeight();
final int tl = rroi.y*width+rroi.x; // top left corner index final int tl = rroi.y*width+rroi.x; // top left corner index
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
for (int n = 0; n < 2; n++) { for (int n = 0; n < 2; n++) {
final int fn = n; final int fn = n;
final float [] fpixels = (float[]) stack.getPixels(slices[fn]); final float [] fpixels = (float[]) stack.getPixels(slices[fn]);
...@@ -864,63 +1429,235 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -864,63 +1429,235 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
test_titles); test_titles);
} }
kernel = extractKernel( double [] kernel = extractKernel(
deconvolved, // double [] data, deconvolved, // double [] data,
kernel_radius, // int kernel_radius, kernel_radius, // int kernel_radius,
hor_sym, // boolean hor_sym, hor_sym, // boolean hor_sym,
vert_sym, // boolean vert_sym, vert_sym, // boolean vert_sym,
all_sym, // boolean all_sym, all_sym, // boolean all_sym,
debugLevel); // int debugLevel) debugLevel); // int debugLevel)
System.out.println("Manually duplicate and save the last slice of teh displayed kernel stack");
System.out.println("createMismatchedResolutionKernel() Done");
// TODO - add automatic saving of the kernel
} }
final double [][] full_img = new double [4][width*height];
final int [] out_slices = {0,2,1,3}; // {high-res, low-res, convoled_high_res,data}
for (int n = 0; n < 2; n++) {
final int fn = out_slices[n];
final float [] fpixels = (float[]) stack.getPixels(slices[n]); public static boolean [] extractFlatPixels(
ai.set(0); final double [][] data,
final int width,
final int slice,
final int [] nonNaNs,
final int slice_other, // <0 if not controlled
final boolean four_sides, // false - only opposite
final double neib_radius,
final double max_heib_diff,
final double neib_radius_other,
final double max_neib_diff_other,
final int debugLevel) {
final boolean dati_subtract_regression=true;
final int regression_temp = slice;
final int height = data[slice].length/width;
final boolean [] flat = new boolean [data[slice].length];
final double [] data_slice = data[slice];
final double [] data_slice_other = (slice_other >=0) ? data[slice_other]: null;
final double [] dati = data[3];
final double [] temp = data[regression_temp];
final int [][] neibs= getCirclePoints(neib_radius);
final int [][] neibs_other= (slice_other >=0) ? getCirclePoints(neib_radius_other):null;
final TileNeibs tn = new TileNeibs(width, height);
final int [][] opposite_pairs= {{1,3},{0,2}}; // may be modified to 8 dirs
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
final int dbg_pix = 1123*2288+1156; // 2037*3060+1631;
for (int ithread = 0; ithread < threads.length; ithread++) { for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() { threads[ithread] = new Thread() {
public void run() { public void run() {
for (int ipix = ai.getAndIncrement(); ipix < full_img[fn].length; ipix = ai.getAndIncrement()) { int [] neibs_ind = new int [4];
full_img[fn][ipix] = fpixels[ipix]; for (int ipix = ai.getAndIncrement(); ipix < data_slice.length; ipix = ai.getAndIncrement()) if (!Double.isNaN(data[slice][ipix])){
if (ipix == dbg_pix) {
System.out.println("extractFlatPixels(): ipix="+ipix);
}
pix_process: {
double d = data_slice[ipix];
for (int s:nonNaNs) {
if (Double.isNaN(data[s][ipix])) {
break pix_process;
}
}
for (int i = 0; i < TileNeibs.DIRS/2; i++) {
neibs_ind[i]= tn.getNeibIndex(ipix, i*2);
}
int num_pairs = 0;
for (int i = 0; i < opposite_pairs.length; i++) {
int p0= neibs_ind[opposite_pairs[i][0]];
int p1= neibs_ind[opposite_pairs[i][1]];
if ((p0>=0) && (p1>=0) && ((data_slice[p0] - d) * (data_slice[p1] - d) > 0)) { // correct NaN
num_pairs++;
}
}
if ((num_pairs > 0) && (!four_sides || (num_pairs > 1))) { // first condition met, try neib_radius
for (int[] ip : neibs) {
int p = tn.getNeibIndex(ipix, ip[0], ip[1]);
if (p >= 0) {
double d1 = data_slice[p];
if (!Double.isNaN(d1)) { // Skip NaNs, they are OK
if (Math.abs(d1-d) > max_heib_diff) {
break pix_process;
}
}
}
}
if (data_slice_other != null) {
double d_other = data_slice_other[ipix];
// Here it should be already tested for NaN as slice_other should be in nonNaNs array to test
for (int[] ip : neibs_other) {
int p = tn.getNeibIndex(ipix, ip[0], ip[1]);
if (p >= 0) {
double d1 = data_slice_other[p];
if (!Double.isNaN(d1)) { // Skip NaNs, they are OK
if (Math.abs(d1-d_other) > max_neib_diff_other) {
break pix_process;
}
}
}
}
} }
flat[ipix] = true;
}
} // pix_process
} // for (int ipix
} }
}; };
} }
ImageDtt.startAndJoin(threads); ImageDtt.startAndJoin(threads);
double [] ab = getDatiRegression(
temp, // final double [] temp,
dati, // final double [] dati,
flat); // final boolean [] mask);
double a = ab[0];
double b = ab[1];
if (debugLevel > -4) {
System.out.println("extractFlatPixels(): regression a="+a+", b="+b);
} }
double [] orig_dati = dati.clone();
if (dati_subtract_regression) {
full_img[out_slices[2]] = convolveWithKernel( // final double [] corr_dati = dati.clone();
full_img[out_slices[0]], // final double [] data,
kernel, // final double [] kernel,
width); // final int width)
// calculate DATI
ai.set(0); ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) { for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() { threads[ithread] = new Thread() {
public void run() { public void run() {
for (int ipix = ai.getAndIncrement(); ipix < full_img[0].length; ipix = ai.getAndIncrement()) { for (int ipix = ai.getAndIncrement(); ipix < data_slice.length; ipix = ai.getAndIncrement())
full_img[out_slices[3]][ipix] = full_img[out_slices[1]][ipix]-full_img[out_slices[2]][ipix]; if (!Double.isNaN(dati[ipix])){
dati[ipix] -= a * temp[ipix] + b;
} }
} }
}; };
} }
ImageDtt.startAndJoin(threads); ImageDtt.startAndJoin(threads);
}
if (debugLevel>0) { if (debugLevel>0) {
String [] test_titles= {"high_res","convolved", "low_res", "DATI"}; String [] test_titles = {"high_res","convolved", "low_res", "DATI", "orig_dati","flat"};
double [][] dbg_img = new double [test_titles.length][];
for (int i = 0; i < data.length; i++) {
dbg_img[i] = data[i];
}
dbg_img[4] = orig_dati;
dbg_img[5] = new double [flat.length];
int num_flat=0;
for (int i = 0; i < flat.length; i++) {
dbg_img[5][i] = flat[i]? 100:Double.NaN;
if (flat[i]) num_flat++;
}
ShowDoubleFloatArrays.showArrays( ShowDoubleFloatArrays.showArrays(
full_img, dbg_img,
width, width,
height, height,
true, true,
removeKnownExtension(imp.getTitle())+"-DATI", "DATI-flat_fs"+four_sides+"_nr"+neib_radius+"_or"+neib_radius_other+"_md"+max_heib_diff+"_od"+max_neib_diff_other,
test_titles); test_titles);
System.out.println("extractFlatPixels(): found "+num_flat+" points.");
System.out.println(" four_sides= "+four_sides);
System.out.println(" neib_radius= "+neib_radius);
System.out.println(" neib_radius_other= "+neib_radius_other);
System.out.println(" max_heib_diff= "+max_heib_diff);
System.out.println(" max_neib_diff_other="+max_neib_diff_other);
} }
System.out.println("testDeconvolveSlices() Done"); return flat;
} }
private static double [] getDatiRegression(
final double [] temp,
final double [] dati,
final boolean [] mask) {
final Thread[] threads = ImageDtt.newThreadArray(QuadCLT.THREADS_MAX);
final AtomicInteger ai = new AtomicInteger(0);
final double [] as0 = new double[threads.length];
final double [] asx = new double[threads.length];
final double [] asx2 = new double[threads.length];
final double [] asy = new double[threads.length];
final double [] asxy= new double[threads.length];
final AtomicInteger ati = new AtomicInteger(0);
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
int thread_num = ati.getAndIncrement();
for (int ipix = ai.getAndIncrement(); ipix < temp.length; ipix = ai.getAndIncrement()) if (mask[ipix]){
double x = temp[ipix];
double y = dati[ipix];
if (!Double.isNaN(x+y)) {
as0 [thread_num] += 1;
asx [thread_num] += x;
asx2[thread_num] += x * x;
asy [thread_num] += y;
asxy[thread_num] += x * y;
}
} // for (int ipix
}
};
}
ImageDtt.startAndJoin(threads);
double s0 = 0.0;
double sx = 0.0;
double sx2 = 0.0;
double sy = 0.0;
double sxy= 0.0;
for (int i = 0; i < threads.length; i++) {
s0+= as0[i];
sx+= asx[i];
sx2+= asx2[i];
sy+= asy[i];
sxy+= asxy[i];
}
double dnm = s0 * sx2 - sx*sx;
double a = (sxy * s0 - sy * sx) / dnm;
double b = (sy * sx2 - sxy * sx) / dnm;
return new double[] {a,b};
}
private static int [][] getCirclePoints(
double radius) {
ArrayList<Point> point_list = new ArrayList<Point>();
double r2 =radius*radius;
int ir = (int) radius;
for (int y = -ir; y <= ir; y++) {
for (int x = -ir; x <= ir; x++) {
if ((y*y+x*x) < r2) point_list.add(new Point(x,y));
}
}
int [][] point_arr = new int[point_list.size()][2];
for (int i = 0; i < point_arr.length; i++) {
point_arr[i][0] = point_list.get(i).x;
point_arr[i][1] = point_list.get(i).y;
}
return point_arr;
}
public static double [] convolveWithKernel( public static double [] convolveWithKernel(
final double [] data, final double [] data,
final double [] kernel, final double [] kernel,
...@@ -976,7 +1713,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -976,7 +1713,7 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
int debugLevel) { int debugLevel) {
hor_sym |= all_sym; hor_sym |= all_sym;
vert_sym |= all_sym; vert_sym |= all_sym;
double [][] dbg_img = (debugLevel > 0)? new double [4][]: null; double [][] dbg_img = (debugLevel > -4)? new double [4][]: null;
int kernel_size = 2*kernel_radius + 1; int kernel_size = 2*kernel_radius + 1;
int fft_size = (int) Math.sqrt(data.length); int fft_size = (int) Math.sqrt(data.length);
double [] kernel = new double [kernel_size*kernel_size]; double [] kernel = new double [kernel_size*kernel_size];
...@@ -2081,6 +2818,653 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{ ...@@ -2081,6 +2818,653 @@ public class OrthoMap implements Comparable <OrthoMap>, Serializable{
return null; return null;
} }
public static double [] correlateWithPattern(
final double [] data,
final int width,
final int psize, // power of 2, such as 64
final double [] pattern, // [psize*psize]
final double phaseCoeff,
final int debugLevel) {
final int height = data.length/width;
final double [] dout = new double [data.length];
final double [] window = new double [psize*psize];
final double [] wnd1d = new double[psize/2];
for (int i = 0; i < psize/2; i++) {
wnd1d[i] = Math.sin((i+0.5)*Math.PI/psize);
}
for (int i = 0; i < psize/2; i++) {
for (int j = 0; j < psize/2; j++) {
double w = wnd1d[i]*wnd1d[j];
int k = i*psize + j;
window[k] = w;
window[window.length-1-k] = w;
k = (i+ 1)*psize - 1 -j;
window[k] = w;
window[window.length-1-k] = w;
}
}
final DoubleFHT doubleFHT0 = new DoubleFHT();
final double [] patternFD = pattern.clone();
doubleFHT0.transformPattern(patternFD);
final double [] filter = null; // probably not needed
if ((width== psize) && (height == psize)) {
testPhaseCorr(
data, // double [] data,
pattern, // double [] pattern,
window, // double [] wnd,
phaseCoeff); // double phaseCoeff)
return data;
}
final int tilesX = (int) Math.ceil(width/(psize/2)) + 1;
final int tilesY = (int) Math.ceil(height/(psize/2)) + 1;
final Thread[] threads = ImageDtt.newThreadArray();
final AtomicInteger ai = new AtomicInteger(0);
final int dbg_tileX = -20;
final int dbg_tileY = -1;
for (int passY = 0; passY < 2; passY++) { // to isolate threads results - no overlap
final int tileY0=passY;
final int ntilesY=(tilesY + 1 - passY)/ 2;
for (int passX = 0; passX < 2; passX++) {
final int tileX0=passX;
final int ntilesX=(tilesX + 1 - passX)/ 2;
final int ntiles = ntilesY*ntilesX;
ai.set(0);
for (int ithread = 0; ithread < threads.length; ithread++) {
threads[ithread] = new Thread() {
public void run() {
double [] dtile = new double [pattern.length];
TileNeibs tn = new TileNeibs(psize,psize);
DoubleFHT doubleFHT = new DoubleFHT();
for (int nTile = ai.getAndIncrement(); nTile < ntiles; nTile = ai.getAndIncrement()){
int tileX = tileX0 - 1 + 2 * (nTile % ntilesX); // nTile % ntilesX;
int tileY = tileY0 - 1 + 2 * (nTile / ntilesX); // nTile / ntilesX;
boolean dbg_tile = (Math.abs(tileX - dbg_tileX) < 1) && (Math.abs(tileY - dbg_tileY) < 1);
// int px0 = (psize/2) *(tileX0 - 1 + 2 * tileX); // absolute in the original/result image
// int py0 = (psize/2) *(tileY0 - 1 + 2 * tileY);
int px0 = (psize/2) * tileX; // absolute in the original/result image
int py0 = (psize/2) * tileY;
int x0 = Math.max(0, -px0);
int y0 = Math.max(0, -py0);
int x1 = Math.min(psize, width- px0);
int y1 = Math.min(psize, height-py0);
if (dbg_tile) {
System.out.println("correlateWithPattern(): tileX="+tileX+", tileY="+tileY);
System.out.println("correlateWithPattern(): px0="+px0+", py0="+py0);
}
int lwidth=x1-x0;
boolean has_NaN = false;
if ((x0>0) || (y0>0) || (x1 < psize) || (y1 < psize)) {
Arrays.fill(dtile,Double.NaN);
}
for (int y = y0; y < y1; y++) {
System.arraycopy(
data,
(py0 + y)*width+(px0+x0),
dtile,
y*psize+x0,
lwidth);
}
for (int i = 0; i < dtile.length; i++) {
if (Double.isNaN(dtile[i])) {
has_NaN = true;
break;
}
}
if (has_NaN) {
fillNaNs(dtile, tn, 3);
}
for (int i = 0; i < dtile.length; i++) {
dtile[i] *= window[i];
}
if (dbg_tile) {
String [] rslt_titles= {"windowed","pattern"};
ShowDoubleFloatArrays.showArrays(
new double[][] {dtile, pattern},
psize,
psize,
true,
"input_corr_tx"+tileX+"_ty"+tileY,
rslt_titles);
}
// now cross-correlate with fat zero dtile and
// important to assign: dtile is modified, but not to the result!
dtile = doubleFHT.phaseCorrelatePattern ( // new
dtile, // double [] first,
patternFD, // double [] secondFD,
phaseCoeff, // double phaseCoeff,
filter, // double [] filter, // high/low pass filtering
null); // double [] first_save ) //null-OK
if (dbg_tile) {
String [] rslt_titles= {"corr","pattern"};
ShowDoubleFloatArrays.showArrays(
new double[][] {dtile, pattern},
psize,
psize,
true,
"corr_tx"+tileX+"_ty"+tileY,
rslt_titles);
}
// multiply result by a window second time and add to the output array
for (int y = y0; y < y1; y++) {
for (int x = x0; x < x1; x++) {
int i = y * psize + x;
dout[(py0 + y) * width+(px0 + x)] += dtile[i]*window[i];
}
}
}
}
};
}
ImageDtt.startAndJoin(threads);
}
}
return dout;
}
public static void testPhaseCorr(
double [] data_in,
double [] pattern_in,
double [] wnd_in,
double phaseCoeff) {
phaseCoeff = 0.5;
double [] data = data_in.clone();
double [] pattern = pattern_in.clone();
double [] wnd = wnd_in.clone();
int size = (int) Math.sqrt(data.length);
double [] data_orig = data.clone();
for (int i = 0; i < data.length; i++) {
data[i] *= wnd[i];
}
{
String [] rslt_titles= {"original","window","windowed","pattern"};
ShowDoubleFloatArrays.showArrays(
new double[][] {data_orig, wnd, data, pattern},
size,
size,
true,
"input_corr",
rslt_titles);
}
DoubleFHT doubleFHT = new DoubleFHT();
double [] data_orig2 = data.clone();
double [] pattern_orig = pattern.clone();
boolean use_new = true; // false;
double [] corr_out = null;
if (use_new) {
final double [] patternFD = pattern.clone();
doubleFHT.transformPattern(patternFD);
corr_out = doubleFHT.phaseCorrelatePattern ( // new
data, // double [] first,
patternFD, // double [] secondFD,
phaseCoeff, // double phaseCoeff,
null, // double [] filter, // high/low pass filtering
null); // double [] first_save ) //null-OK
} else {
corr_out = doubleFHT. phaseCorrelate ( // new
data, // double [] first,
pattern, // double [] second,
phaseCoeff, // double phaseCoeff,
null, // double [] filter, // high/low pass filtering
null, // double [] first_save,
null); // double [] second_save ) //null-OK
}
{
String [] rslt_titles= {"corr", "original","window","windowed","pattern"};
ShowDoubleFloatArrays.showArrays(
new double[][] {corr_out, data_orig, wnd, data_orig2, pattern_orig},
size,
size,
true,
"output_corr"+phaseCoeff,
rslt_titles);
}
System.out.println("testPhaseCorr() done");
}
public static void fillNaNs(
double [] data,
TileNeibs tn,
int min_neibs) { // 3
boolean dbg=false;
double [] data_orig = dbg? data.clone():null;
boolean [] known = new boolean [data.length];
boolean [] wave = new boolean [data.length];
double wdiag=0.7;
double [] weights = {1.0,wdiag,1.0,wdiag,1.0,wdiag,1.0,wdiag};
ArrayList<Integer> front_list = new ArrayList<Integer>();
for (int i = 0; i < data.length; i++) {
if (Double.isNaN(data[i]) && !wave[i]) {
int num_neibs = 0;
for (int dir=0; dir < TileNeibs.DIRS; dir++) {
int i1 = tn.getNeibIndex(i, dir);
if ((i1 >= 0) && !Double.isNaN(data[i1])) {
num_neibs++;
}
}
if (num_neibs >= min_neibs) {
front_list.add(i);
wave[i]=true;
}
} else {
known[i] = true;
}
}
while (!front_list.isEmpty()) {
for (int i:front_list) {
double swd=0,sw=0;
for (int dir=0; dir < TileNeibs.DIRS; dir++) {
int i1 = tn.getNeibIndex(i, dir);
if ((i1 >= 0) && (known[i1])) {
double d = data[i1];
if (!Double.isNaN(d)) {
sw += weights[dir];
swd += weights[dir] * d;
}
}
}
data[i] = swd/sw;
}
int was_size = front_list.size();
for (int j = 0; j< was_size; j++) {
int i = front_list.remove(0);
known[i] = true;
for (int dir=0; dir < TileNeibs.DIRS; dir++) {
int i1 = tn.getNeibIndex(i, dir);
if ((i1 >= 0) && Double.isNaN(data[i1]) && !wave[i1]) {
// check number of known neighbors
int num_neibs = 0;
for (int dir1 = 0; dir1 < TileNeibs.DIRS; dir1++) {
int i2 = tn.getNeibIndex(i1, dir1);
if ((i2 >= 0) && known[i2]) {
num_neibs++;
}
}
if (num_neibs >= min_neibs) {
front_list.add(i1);
wave[i1]=true;
}
}
}
}
}
if (data_orig != null) {
String [] rslt_titles= {"original","filled"};
ShowDoubleFloatArrays.showArrays(
new double[][] {data_orig, data},
tn.getSizeX(),
tn.getSizeY(),
true,
"fill_NaN",
rslt_titles);
}
}
public static void testPatternCorrelate(
ImagePlus imp_src) {
String pattern_dir= "/media/elphel/SSD3-4GB/lwir16-proc/ortho_videos/debug/mines/pattern_25m_zoom1/synthetic/";
String pattern_file= "patterns_50m_zoom1_200x200.tiff";
int zoomout = 1;
int corr_size = 256;
double phaseCoeff = 0.5;
int debugLevel = 1;
GenericJTabbedDialog gd = new GenericJTabbedDialog("Correlate image with pattern",1200,500);
gd.addStringField ("Pattern directory", pattern_dir, 180, "Absolute path including trailing \"/\".");
gd.addStringField ("Pattern filename", pattern_file,180, "Pattern filename.");
gd.addNumericField("Zoom-out factor", zoomout, 0,4,"x", "Reduce pattern resolution to match image.");
gd.addNumericField("Phase correlation coefficient", phaseCoeff, 3,7,"","1.0 - pure phase correlation, 0.0 - regular correlation.");
gd.addNumericField("Debug level", debugLevel, 0,4,"", "Debug level.");
gd.showDialog();
if (gd.wasCanceled()) return;
pattern_dir= gd.getNextString();
pattern_file= gd.getNextString();
zoomout= (int) gd.getNextNumber();
phaseCoeff= gd.getNextNumber();
debugLevel= (int) gd.getNextNumber();
float [] fpixels = (float[]) imp_src.getProcessor().getPixels();
int width = imp_src.getWidth();
int height = imp_src.getHeight();
double [] data = new double [fpixels.length];
for (int i = 0; i < data.length; i++) {
data[i] = fpixels[i];
}
// get pattern(s)
String pattern_path=pattern_dir+pattern_file;
ImagePlus imp_pattern = new ImagePlus(pattern_path);
int pattern_size = imp_pattern.getWidth();
if (pattern_size == 0) {
System.out.println("testPatternCorrelate(): pattern \""+pattern_path+"\" is not found.");
return;
}
double [] kernel = getConvolutionKernel();
ImageStack stack_pattern = imp_pattern.getStack();
int nSlices = stack_pattern.getSize();
double [][] patterns = new double[nSlices][];
String [] pattern_labels = new String[nSlices];
for (int n = 0; n < patterns.length; n++) {
pattern_labels[n]=stack_pattern.getShortSliceLabel(n+1);
float [] fpixels_pattern = (float[]) stack_pattern.getPixels(n+1);
patterns[n]=new double[fpixels_pattern.length];
for (int i = 0; i < fpixels_pattern.length; i++) {
patterns[n][i] = fpixels_pattern[i];
}
}
/*
float [] fpixels_pattern = (float[]) imp_pattern.getProcessor().getPixels();
double [ ] pattern = new double[fpixels_pattern.length];
for (int i = 0; i < fpixels_pattern.length; i++) {
pattern[i] = fpixels_pattern[i];
}
*/
double [][] corrs_out = new double[patterns.length][];
for (int n = 0; n < patterns.length; n++) {
double [] corr_pattern = patternZoomCropPad(
patterns[n], // double [] pattern,
pattern_size, // int pattern_size,
corr_size, // int size,
zoomout, // int zoomout,
true); // out_normalize); // boolean normalize)
if (kernel != null) {
corr_pattern = convolveWithKernel(
corr_pattern, // final double [] data,
kernel, // final double [] kernel,
corr_size); // final int width)
}
corrs_out[n]= correlateWithPattern(
data, // final double [] data,
width, // final int width,
corr_size, // final int psize, // power of 2, such as 64
corr_pattern, //final double [] pattern, // [psize*psize]
phaseCoeff, // final double phaseCoeff,
debugLevel); // final int debugLevel) {
}
String [] rslt_titles =new String [patterns.length+1];
double [][] rslt = new double [patterns.length + 1][]; // ]{data, corr_out};
rslt[0] = data;
rslt_titles[0] = "original";
for (int i = 0; i < patterns.length; i++) {
rslt[i+1] = corrs_out[i];
rslt_titles[i+1] = "corr_"+pattern_labels[i];
}
ShowDoubleFloatArrays.showArrays(
rslt,
width,
height,
true,
removeKnownExtension(imp_src.getTitle())+"-PATTERN_CORR",
rslt_titles);
System.out.println("testPatternCorrelate(): correlation done");
}
public static void testPatternGenerate() {
int half_size = 100;
boolean half_pix = false; // center between pixels
double radius = 32; // 32;
double edge= 15; // 4;
double radius_in= 15; // if 0 - skip
double edge_in = 8;
double scale_in = -0.05;
double scale = -200; // black 200
boolean normalize = false;
int halves_number = 8; // hidden mines - number of diameter cuts
double cut_frac = 0.6; // half-cut width fraction of diameter
int corr_size = 128;
int zoomout = 2;
boolean out_normalize = true;
GenericJTabbedDialog gd = new GenericJTabbedDialog("Create circulat pattern",1200,600);
gd.addNumericField("Half pattern size", half_size, 0,4,"pix", "Pattern will be square, twice this side.");
gd.addCheckbox ("Pattern center between pixels", half_pix, "If false - center will be at integer pixel X/Y.");
gd.addNumericField("Radius at half height", radius, 3,7,"pix", ".");
gd.addNumericField("Edge width", edge, 3,7,"pix", "Transition from full on to zero.");
gd.addNumericField("Optional inner feature radius", radius_in, 3,7,"pix", "Put 0 if not needed.");
gd.addNumericField("Inner feature edge", edge_in, 3,7,"pix", "Transition for the inner feature.");
gd.addNumericField("Inner feature scale", scale_in, 3,7,"", "Negative OK.");
gd.addNumericField("Overall pattern scale", scale, 3,7,"", "Multiply pattern by this number.");
gd.addCheckbox ("Normalize", normalize, "Make sum of all pixel equal to 1.0");
gd.addMessage("--- obscured round objects, trying halves ---");
gd.addNumericField("Halves directions", halves_number, 0,4,"", "Number of directions for diameter cuts for obscured objects.");
gd.addNumericField("Cut relative widths", cut_frac, 3,7,"", "Width of diameter cuts as a fraction of diameter.");
gd.addMessage("--- testing scaling/cropping/padding ---");
gd.addNumericField("Output size (power of 2)", corr_size, 0,4,"pix", "Output size for correlation.");
gd.addNumericField("Zoom out factor", zoomout, 0,4,"", "Zoom out factor for correlation.");
gd.addCheckbox ("Normalize zoomed out", out_normalize, "Make sum of all pixel equal to 1.0");
gd.showDialog();
if (gd.wasCanceled()) return;
half_size= (int) gd.getNextNumber();
half_pix= gd.getNextBoolean();
radius= gd.getNextNumber();
edge= gd.getNextNumber();
radius_in= gd.getNextNumber();
edge_in= gd.getNextNumber();
scale_in= gd.getNextNumber();
scale= gd.getNextNumber();
normalize= gd.getNextBoolean();
halves_number= (int) gd.getNextNumber();
cut_frac= gd.getNextNumber();
corr_size= (int) gd.getNextNumber();
zoomout= (int) gd.getNextNumber();
out_normalize= gd.getNextBoolean();
double [] kernel = getConvolutionKernel();
int size = 2*half_size;
double [] pattern = generateDoubleCircularPattern (
half_size, // int half_size,
half_pix, // boolean half_pix, // center between pixels
radius, // double radius,
edge, // double edge,
radius_in, // double radius_in, // if 0 - skip
edge_in, // double edge_in,
scale_in); // double scale_in)
if (normalize) {
double sum_pix = 0;
for (int i = 0; i < pattern.length; i++) {
sum_pix+=pattern[i];
}
scale = 1.0/sum_pix;
}
for (int i = 0; i < pattern.length; i++) {
pattern[i]*=scale;
}
double [][] patterns = new double [halves_number+1][];
double xc= half_size + (half_pix? 0.5:0);
double yc= half_size + (half_pix? 0.5:0);
for (int n = 0; n < patterns.length; n++) {
patterns[n] = pattern.clone();
if (n > 0) {
double phi = 2*Math.PI*(n - 1)/halves_number;
double ck = Math.cos(phi)/(cut_frac * radius);
double sk = Math.sin(phi)/(cut_frac * radius);
for (int iy = 0; iy< size; iy++) {
double dy = iy - yc;
for (int ix = 0; ix < size; ix++) {
double dx = ix - xc;
double l = dx*ck+dy*sk;
if (l <= -1) {
patterns[n][iy* size + ix] = 0;
} else if (l < 1) {
patterns[n][iy* size + ix] *= 0.5*(1 + Math.sin(l* Math.PI/2));
}
}
}
}
}
String pattern_name = "patterns_r"+radius+"_e"+edge+"_ir"+radius_in+
"ie"+edge_in+"_is"+scale_in+"_h"+halves_number+"_w"+cut_frac+
"_s"+scale+ "_"+size+"x"+size;
String [] pattern_titles = new String [patterns.length];
pattern_titles[0] = "full";
for (int n = 1; n < pattern_titles.length; n++) {
pattern_titles[n] = String.format("%f5.1 deg", 360.0*(n-1)/halves_number);
}
// public static ImagePlus makeArrays(double[][] pixels, int width, int height, String title, String [] titles) {
ImagePlus imp = ShowDoubleFloatArrays.makeArrays(
patterns, // float[] pixels,
size,
size,
pattern_name,
pattern_titles);
imp.show();
double [][] corr_patterns = new double [patterns.length][];
for (int n = 0; n < corr_patterns.length; n++) {
corr_patterns[n] = patternZoomCropPad(
patterns[n], // double [] pattern,
size, // int pattern_size,
corr_size, // int size,
zoomout, // int zoomout,
out_normalize); // boolean normalize)
if (kernel != null) {
corr_patterns[n] = convolveWithKernel(
corr_patterns[n], // final double [] data,
kernel, // final double [] kernel,
corr_size); // final int width)
}
}
ShowDoubleFloatArrays.showArrays(
corr_patterns,
corr_size,
corr_size,
true,
"corr_patterns_h"+halves_number+"_w"+cut_frac+
"_"+corr_size+"x"+corr_size+"-zoomout"+zoomout,
pattern_titles);
System.out.println("testPatternGenerate(): Patterns created");
}
public static double [] patternZoomCropPad(
double [] pattern,
int pattern_size,
int size,
int zoomout,
boolean normalize) {
double [] pout = new double[size*size];
int xyc = size/2;
int pxyc = pattern_size/2;
double sum_pix = 0.0;
for (int y = 0; y < size; y++) {
int py = pxyc + (y-xyc)*zoomout;
if ((py >=0) && (py < pattern_size)) {
for (int x = 0; x < size; x++) {
int px = pxyc + (x-xyc)*zoomout;
if ((px >=0) && (px < pattern_size)) {
double d = pattern[py*pattern_size + px];
pout[y*size+x] = d;
sum_pix += d;
}
}
}
}
if (normalize) {
// double scale = 1.0/sum_pix;
double scale = -10000.0/sum_pix; //TODO: calculate better
for (int i = 0; i < pout.length; i++) {
pout[i] *= scale;
}
}
return pout;
}
public static double [] generateDoubleCircularPattern (
int half_size,
boolean half_pix, // center between pixels
double radius,
double edge,
double radius_in, // if 0 - skip
double edge_in,
double scale_in) {
double [] pattern = generateCircularPattern (
half_size, // int half_size,
half_pix, // boolean half_pix, // center between pixels
radius, // double radius,
edge); // double edge)
if (radius_in > 0) {
double [] pattern1 = generateCircularPattern (
half_size, // int half_size,
half_pix, // boolean half_pix, // center between pixels
radius_in, // double radius,
edge_in); // double edge)
for (int i = 0; i < pattern.length; i++) {
pattern[i] += scale_in*pattern1[i];
}
}
return pattern;
}
public static double [] generateCircularPattern (
int half_size,
boolean half_pix, // center between pixels
double radius,
double edge) {
double min_r = radius - edge/2;
double max_r = radius + edge/2;
double [] profile = new double [(int)Math.ceil(max_r)+1];
for (int i = 0; i < profile.length; i++) {
if (i < min_r) {
profile[i] = 1;
} else if (i < max_r) {
profile[i] = 0.5 * (1.0 +Math.cos(Math.PI *(i-min_r)/(max_r-min_r)));
}
}
return generateRotationalPattern (
half_size, // int half_size,
half_pix, // boolean half_pix, // center between pixels
profile); // double [] profile)
}
public static double [] generateRotationalPattern (
int half_size,
boolean half_pix, // center between pixels
double [] profile) {
int size = 2 * half_size;
double [] pattern = new double [size*size];
double xc= half_size + (half_pix? 0.5:0);
double yc= half_size + (half_pix? 0.5:0);
for (int iy = 0; iy< size; iy++) {
for (int ix = 0; ix < size; ix++) {
double r = Math.sqrt((ix-xc)*(ix-xc) + (iy-yc)*(iy-yc));
int ir = (int) Math.floor(r);
if (ir < profile.length) {
double d0 = profile[ir];
double d1 = (ir < (profile.length-1))?profile[ir+1] : 0;
double fr = r - ir;
pattern[iy* size + ix] = fr * d1 + (1-fr) * d0;
}
}
}
return pattern;
}
......
src/main/java/com/elphel/imagej/orthomosaic/OrthoMapsCollection.java
View file @ 81c431ad
...@@ -60,8 +60,8 @@ public class OrthoMapsCollection implements Serializable{ ...@@ -60,8 +60,8 @@ public class OrthoMapsCollection implements Serializable{
ortho_maps = new OrthoMap[paths.length]; ortho_maps = new OrthoMap[paths.length];
for (int n = 0; n < ortho_maps.length; n++) { for (int n = 0; n < ortho_maps.length; n++) {
String name = OrthoMap.getNameFromPath(paths[n]); String name = OrthoMap.getNameFromPath(paths[n]);
String scene_path = getSceneDir(scenes_path, name); // String scene_path = getSceneDir(scenes_path, name);
ortho_maps[n] = new OrthoMap(paths[n], scene_path); ortho_maps[n] = new OrthoMap(paths[n]); // , scene_path);
ortho_maps[n].setAffine(affine); ortho_maps[n].setAffine(affine);
} }
} }
...@@ -117,9 +117,15 @@ public class OrthoMapsCollection implements Serializable{ ...@@ -117,9 +117,15 @@ public class OrthoMapsCollection implements Serializable{
public void getAllTemperatures() { public void getAllTemperatures() {
for (int n = 0; n < ortho_maps.length; n++) { for (int n = 0; n < ortho_maps.length; n++) {
ortho_maps[n].getTemperature(); ortho_maps[n].getTemperature();
ortho_maps[n].getTemperatureAndTelemetry();
} }
} }
public void updateLLa() { // make automatic by saving/testing file modification stamp?
for (int n = 0; n < ortho_maps.length; n++) {
ortho_maps[n].updateLLA();
}
}
public static String[] getPathsFromSorceList( public static String[] getPathsFromSorceList(
String files_list) { String files_list) {
...@@ -683,6 +689,12 @@ public class OrthoMapsCollection implements Serializable{ ...@@ -683,6 +689,12 @@ public class OrthoMapsCollection implements Serializable{
} }
gpu_pair_img[n] = ortho_maps[gpu_pair[n]].getPaddedGPU (zoom_lev); // int zoom_level, gpu_pair_img[n] = ortho_maps[gpu_pair[n]].getPaddedGPU (zoom_lev); // int zoom_level,
} }
boolean invert_second = (debugLevel > 1000);
if (invert_second) {
for (int i = 0; i < gpu_pair_img[1].length; i++) {
gpu_pair_img[1][i]= -gpu_pair_img[1][i];
}
}
if (show_gpu_img) { if (show_gpu_img) {
String [] dbg_titles = {ortho_maps[gpu_pair[0]].getName(),ortho_maps[gpu_pair[1]].getName()}; String [] dbg_titles = {ortho_maps[gpu_pair[0]].getName(),ortho_maps[gpu_pair[1]].getName()};
ShowDoubleFloatArrays.showArrays( ShowDoubleFloatArrays.showArrays(
...@@ -1271,7 +1283,7 @@ public class OrthoMapsCollection implements Serializable{ ...@@ -1271,7 +1283,7 @@ public class OrthoMapsCollection implements Serializable{
num_passes, // int num_passes, num_passes, // int num_passes,
max_diff, // double max_diff, max_diff, // double max_diff,
ImageDtt.THREADS_MAX, // int threadsMax, ImageDtt.THREADS_MAX, // int threadsMax,
debugLevel); // int debug_level) debugLevel-1); // int debug_level)
} }
if (debugLevel > 0) { if (debugLevel > 0) {
String [] dbg_titles = {"x-raw","x_filled","y-raw","y_filled"}; String [] dbg_titles = {"x-raw","x_filled","y-raw","y_filled"};
......
src/main/java/com/elphel/imagej/orthomosaic/SensorTemperatureData.java 0 → 100644
View file @ 81c431ad
package com.elphel.imagej.orthomosaic;
import java.io.IOException;
import java.io.ObjectInputStream;
import java.io.ObjectOutputStream;
import java.io.Serializable;
public class SensorTemperatureData implements Serializable {
private static final long serialVersionUID = 1L;
double averagePixelValue;
double TLM_FPA_KELV;
double TLM_FFC_KELV;
double TLM_CORE_TEMP;
int TLM_FRAME_FFC;
int TLM_FRAME;
private void writeObject(ObjectOutputStream oos) throws IOException {
oos.defaultWriteObject();
}
private void readObject(ObjectInputStream ois) throws ClassNotFoundException, IOException {
ois.defaultReadObject();
}
}
src/main/java/com/elphel/imagej/tileprocessor/TileNeibs.java
View file @ 81c431ad
...@@ -100,10 +100,10 @@ public class TileNeibs{ ...@@ -100,10 +100,10 @@ public class TileNeibs{
int getY(int indx) {return indx / sizeX;}; int getY(int indx) {return indx / sizeX;};
int getSizeX() { public int getSizeX() {
return sizeX; return sizeX;
} }
int getSizeY() { public int getSizeY() {
return sizeY; return sizeY;
} }
......
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