Commit c432b7da authored by Andrey Filippov's avatar Andrey Filippov

CODEX: Add resident LWIR conditioning kernels

Co-authored-by: 's avatarCodex <codex@elphel.com>
parent fcd254f6
...@@ -20,6 +20,7 @@ ...@@ -20,6 +20,7 @@
#include "tp_paths.h" #include "tp_paths.h"
#include "tp_utils.h" // copyalloc_kernel_gpu / copyalloc_pointers_gpu / copyalloc_image_gpu / alloc_kernel_gpu #include "tp_utils.h" // copyalloc_kernel_gpu / copyalloc_pointers_gpu / copyalloc_image_gpu / alloc_kernel_gpu
#include "tp_files.h" // readFloatsFromFile / readAllFloatsFromFile #include "tp_files.h" // readFloatsFromFile / readAllFloatsFromFile
#include "tp_condition.h" // LWIR scale/offset conditioning // By Codex on 07/15/2026
#include "tp_consolidate.h" // TD_TILE_FLOATS (host use in tp_proc_exec_consolidate) // By Claude on 07/12/2026 #include "tp_consolidate.h" // TD_TILE_FLOATS (host use in tp_proc_exec_consolidate) // By Claude on 07/12/2026
#include "tp_peak.h" // PEAK_ROW_FLOATS/PEAK_MAX_PIX (host use in tp_proc_exec_corr2d_peak) // By Claude on 07/13/2026 #include "tp_peak.h" // PEAK_ROW_FLOATS/PEAK_MAX_PIX (host use in tp_proc_exec_corr2d_peak) // By Claude on 07/13/2026
#include "tp_lma.h" // lma_normal_products prototype (roadmap rung 3) // By Codex on 07/14/2026 #include "tp_lma.h" // lma_normal_products prototype (roadmap rung 3) // By Codex on 07/14/2026
...@@ -40,6 +41,7 @@ static const char* TP_DEFINES = ...@@ -40,6 +41,7 @@ static const char* TP_DEFINES =
"#define LIST_TEXTURE_BIT 8\n#define TEXT_NTILE_SHIFT 9\n#define FAT_ZERO_WEIGHT 0.0001\n" "#define LIST_TEXTURE_BIT 8\n#define TEXT_NTILE_SHIFT 9\n#define FAT_ZERO_WEIGHT 0.0001\n"
"#define THREADS_DYNAMIC_BITS 5\n#define RBYRDIST_LEN 5001\n#define RBYRDIST_STEP 0.0004\n#define TILES_PER_BLOCK_GEOM 2\n"; "#define THREADS_DYNAMIC_BITS 5\n#define RBYRDIST_LEN 5001\n#define RBYRDIST_STEP 0.0004\n#define TILES_PER_BLOCK_GEOM 2\n";
static const char* SRC_FILES[] = {"dtt8x8.h","dtt8x8.cu","geometry_correction.h","geometry_correction.cu","TileProcessor.h","TileProcessor.cu", static const char* SRC_FILES[] = {"dtt8x8.h","dtt8x8.cu","geometry_correction.h","geometry_correction.cu","TileProcessor.h","TileProcessor.cu",
"tp_condition.h","tp_condition.cu", // LWIR two-map conditioning + double validation tier // By Codex on 07/15/2026
"tp_consolidate.h","tp_consolidate.cu", // TD consolidation chain (roadmap 2c) // By Claude on 07/12/2026 "tp_consolidate.h","tp_consolidate.cu", // TD consolidation chain (roadmap 2c) // By Claude on 07/12/2026
"tp_peak.h","tp_peak.cu", // peak measurement kernel (roadmap 3-C rung 2) // By Claude on 07/13/2026 "tp_peak.h","tp_peak.cu", // peak measurement kernel (roadmap 3-C rung 2) // By Claude on 07/13/2026
"tp_lma.h","tp_lma.cu"}; // normal-equation products (roadmap rung 3) // By Codex on 07/14/2026 "tp_lma.h","tp_lma.cu"}; // normal-equation products (roadmap rung 3) // By Codex on 07/14/2026
...@@ -48,6 +50,7 @@ static const char* KERNELS[] = { ...@@ -48,6 +50,7 @@ static const char* KERNELS[] = {
"textures_nonoverlap","generate_RBGA","clear_texture_list","mark_texture_tiles","mark_texture_neighbor_tiles", "textures_nonoverlap","generate_RBGA","clear_texture_list","mark_texture_tiles","mark_texture_neighbor_tiles",
"gen_texture_list","clear_texture_rbga","textures_accumulate","create_nonoverlap_list","erase_clt_tiles", "gen_texture_list","clear_texture_rbga","textures_accumulate","create_nonoverlap_list","erase_clt_tiles",
"calculate_tiles_offsets","calc_rot_deriv","calcReverseDistortionTable", "calculate_tiles_offsets","calc_rot_deriv","calcReverseDistortionTable",
"condition_lwir_f","condition_lwir_d", // production float maps + Java-matched validation // By Codex on 07/15/2026
"clt_average_sensors","index_consolidate","consolidate_oob","clt_average_sensors_list", // By Claude on 07/12/2026 "clt_average_sensors","index_consolidate","consolidate_oob","clt_average_sensors_list", // By Claude on 07/12/2026
"corr2D_peak_eig","corr2D_peak_eig_f", // double + float instantiations (precision two-step) // By Claude on 07/13-14/2026 "corr2D_peak_eig","corr2D_peak_eig_f", // double + float instantiations (precision two-step) // By Claude on 07/13-14/2026
"lma_normal_products"}; // deterministic double normal-equation products // By Codex on 07/14/2026 "lma_normal_products"}; // deterministic double normal-equation products // By Codex on 07/14/2026
...@@ -592,6 +595,12 @@ struct TpProc { ...@@ -592,6 +595,12 @@ struct TpProc {
bool rbr_ready; // reverse-distortion table is valid for the uploaded gc bool rbr_ready; // reverse-distortion table is valid for the uploaded gc
std::vector<float*> kernels_h, offs_h, images_h, clt_h, clt_ref_h; std::vector<float*> kernels_h, offs_h, images_h, clt_h, clt_ref_h;
float **gpu_kernels, **gpu_kernel_offsets, **gpu_images, **gpu_clt, **gpu_clt_ref; float **gpu_kernels, **gpu_kernel_offsets, **gpu_images, **gpu_clt, **gpu_clt_ref;
// LWIR conditioning calibration is resident and reused across scene uploads.
// The float pair is the production two-map model; the doubles preserve the
// current Java expression as a validation oracle. By Codex on 07/15/2026.
float *cond_scale_maps, *cond_offset_maps;
double *cond_scales, *cond_offsets, *cond_scales2, *cond_fpn;
bool cond_maps_ready, cond_oracle_ready;
// The live kernel task pointer aliases one of two persistent slots. Slot 0 is // The live kernel task pointer aliases one of two persistent slots. Slot 0 is
// also the backward-compatible tp_proc_set_tasks target; slot 1 preserves the // also the backward-compatible tp_proc_set_tasks target; slot 1 preserves the
// motion-blur partner so consolidation can consume both post-offset streams // motion-blur partner so consolidation can consume both post-offset streams
...@@ -650,6 +659,9 @@ TpProc* tp_proc_create(TpModule* m){ ...@@ -650,6 +659,9 @@ TpProc* tp_proc_create(TpModule* m){
p->task_slot_ntiles[0]=p->task_slot_ntiles[1]=0; p->task_slot_ntiles[0]=p->task_slot_ntiles[1]=0;
p->active_task_slot=-1; // By Codex on 07/15/2026 p->active_task_slot=-1; // By Codex on 07/15/2026
p->gpu_kernels=p->gpu_kernel_offsets=p->gpu_images=p->gpu_clt=p->gpu_clt_ref=nullptr; p->gpu_kernels=p->gpu_kernel_offsets=p->gpu_images=p->gpu_clt=p->gpu_clt_ref=nullptr;
p->cond_scale_maps=p->cond_offset_maps=nullptr;
p->cond_scales=p->cond_offsets=p->cond_scales2=p->cond_fpn=nullptr;
p->cond_maps_ready=p->cond_oracle_ready=false; // By Codex on 07/15/2026
p->gc=p->cv=p->rbr=p->rot=nullptr; p->gc=p->cv=p->rbr=p->rot=nullptr;
p->rbr_ready=false; p->rbr_ready=false;
p->have_rbg=false; p->gpu_corr_images=nullptr; p->have_rbg=false; p->gpu_corr_images=nullptr;
...@@ -711,6 +723,59 @@ int tp_proc_set_kernel_offsets(TpProc* p, int cam, const float* d, int n){ if(!p ...@@ -711,6 +723,59 @@ int tp_proc_set_kernel_offsets(TpProc* p, int cam, const float* d, int n){ if(!p
// pitched image upload (reuses harness update_image_gpu: dstride in floats) // pitched image upload (reuses harness update_image_gpu: dstride in floats)
int tp_proc_set_image(TpProc* p, int cam, const float* d){ if(!p||cam<0||cam>=p->num_cams)return -1; cuCtxSetCurrent(p->mod->ctx); int tp_proc_set_image(TpProc* p, int cam, const float* d){ if(!p||cam<0||cam>=p->num_cams)return -1; cuCtxSetCurrent(p->mod->ctx);
update_image_gpu((float*)d, p->images_h[cam], p->dstride_img, p->img_w, p->img_h); return 0; } // dstride in BYTES (cudaMemcpy2D dpitch) update_image_gpu((float*)d, p->images_h[cam], p->dstride_img, p->img_w, p->img_h); return 0; } // dstride in BYTES (cudaMemcpy2D dpitch)
// De-pitch one resident source image. Primarily a correctness gate for
// conditioning; also useful for explicit ingest diagnostics.
int tp_proc_get_image(TpProc* p, int cam, float* out){ if(!p||!out||cam<0||cam>=p->num_cams)return -1; cuCtxSetCurrent(p->mod->ctx);
return cudaMemcpy2D(out,(size_t)p->img_w*sizeof(float), p->images_h[cam],p->dstride_img,
(size_t)p->img_w*sizeof(float),p->img_h,cudaMemcpyDeviceToHost)==cudaSuccess?0:-2; }
// Upload resident production maps, contiguous [camera][height][width]. The
// scale map is deliberately per-pixel even while current calibration is
// sensor-global: this freezes the interface needed for multiplicative dust.
int tp_proc_set_conditioning_maps(TpProc* p, const float* scales, const float* offsets, int num_pixels){
if(!p||!scales||!offsets||(num_pixels!=p->img_w*p->img_h)) return -1;
cuCtxSetCurrent(p->mod->ctx); p->cond_maps_ready=false;
const size_t n=(size_t)p->num_cams*num_pixels;
if(!p->cond_scale_maps && cudaMalloc((void**)&p->cond_scale_maps,n*sizeof(float))!=cudaSuccess){ seterr("conditioning scale-map alloc failed"); return -2; }
if(!p->cond_offset_maps && cudaMalloc((void**)&p->cond_offset_maps,n*sizeof(float))!=cudaSuccess){ seterr("conditioning offset-map alloc failed"); return -3; }
if(cudaMemcpy(p->cond_scale_maps,scales,n*sizeof(float),cudaMemcpyHostToDevice)!=cudaSuccess){ seterr("conditioning scale-map upload failed"); return -4; }
if(cudaMemcpy(p->cond_offset_maps,offsets,n*sizeof(float),cudaMemcpyHostToDevice)!=cudaSuccess){ seterr("conditioning offset-map upload failed"); return -5; }
p->cond_maps_ready=true; return 0;
}
// Upload the unfused Java-model calibration. fpn is contiguous
// [camera][height][width]; dc_level is passed per launch so a later sequence
// switch does not require rebuilding the per-pixel array.
int tp_proc_set_conditioning_oracle(TpProc* p, const double* scales, const double* offsets,
const double* scales2, const double* fpn, int num_pixels){
if(!p||!scales||!offsets||!scales2||!fpn||(num_pixels!=p->img_w*p->img_h)) return -1;
cuCtxSetCurrent(p->mod->ctx); p->cond_oracle_ready=false;
const size_t nc=(size_t)p->num_cams, n=nc*num_pixels;
if(!p->cond_scales && cudaMalloc((void**)&p->cond_scales, nc*sizeof(double))!=cudaSuccess) return -2;
if(!p->cond_offsets && cudaMalloc((void**)&p->cond_offsets,nc*sizeof(double))!=cudaSuccess) return -3;
if(!p->cond_scales2 && cudaMalloc((void**)&p->cond_scales2,nc*sizeof(double))!=cudaSuccess) return -4;
if(!p->cond_fpn && cudaMalloc((void**)&p->cond_fpn, n*sizeof(double))!=cudaSuccess) return -5;
if(cudaMemcpy(p->cond_scales, scales, nc*sizeof(double), cudaMemcpyHostToDevice)!=cudaSuccess) return -6;
if(cudaMemcpy(p->cond_offsets,offsets,nc*sizeof(double), cudaMemcpyHostToDevice)!=cudaSuccess) return -7;
if(cudaMemcpy(p->cond_scales2,scales2,nc*sizeof(double),cudaMemcpyHostToDevice)!=cudaSuccess) return -8;
if(cudaMemcpy(p->cond_fpn,fpn,n*sizeof(double),cudaMemcpyHostToDevice)!=cudaSuccess) return -9;
p->cond_oracle_ready=true; return 0;
}
// In-place conditioning of all camera image buffers. use_double selects the
// correctness oracle; production uses the float two-map kernel.
int tp_proc_exec_conditioning(TpProc* p, int use_double, double dc_level){
if(!p){ seterr("exec_conditioning: null"); return -1; }
if(use_double ? !p->cond_oracle_ready : !p->cond_maps_ready){ seterr("exec_conditioning: calibration not uploaded"); return -2; }
cuCtxSetCurrent(p->mod->ctx);
const char* name=use_double?"condition_lwir_d":"condition_lwir_f";
CUfunction f=getfun(p->mod,name); if(!f){ seterr("%s missing",name); return -3; }
int nc=p->num_cams,w=p->img_w,h=p->img_h; size_t stride=p->dstride_img/sizeof(float);
void* a_float[]={&nc,&w,&h,&stride,&p->gpu_images,&p->cond_scale_maps,&p->cond_offset_maps};
void* a_double[]={&nc,&w,&h,&stride,&p->gpu_images,&p->cond_scales,&p->cond_offsets,
&p->cond_scales2,&p->cond_fpn,&dc_level};
return launch1(f,(w+31)/32,(h+7)/8,nc,32,8,1,use_double?a_double:a_float,name)?-4:0;
}
// use_center_image: broadcast one center image to all sensors (FPN back-prop mode) // use_center_image: broadcast one center image to all sensors (FPN back-prop mode)
int tp_proc_set_center_image(TpProc* p, const float* d){ if(!p)return -1; cuCtxSetCurrent(p->mod->ctx); int tp_proc_set_center_image(TpProc* p, const float* d){ if(!p)return -1; cuCtxSetCurrent(p->mod->ctx);
for(int c=0;c<p->num_cams;c++) update_image_gpu((float*)d, p->images_h[c], p->dstride_img, p->img_w, p->img_h); return 0; } for(int c=0;c<p->num_cams;c++) update_image_gpu((float*)d, p->images_h[c], p->dstride_img, p->img_w, p->img_h); return 0; }
...@@ -1341,6 +1406,8 @@ void tp_proc_destroy(TpProc* p){ ...@@ -1341,6 +1406,8 @@ void tp_proc_destroy(TpProc* p){
for(int c=0;c<p->num_cams;c++){ if(c<(int)p->kernels_h.size())cudaFree(p->kernels_h[c]); if(c<(int)p->offs_h.size())cudaFree(p->offs_h[c]); for(int c=0;c<p->num_cams;c++){ if(c<(int)p->kernels_h.size())cudaFree(p->kernels_h[c]); if(c<(int)p->offs_h.size())cudaFree(p->offs_h[c]);
if(c<(int)p->images_h.size())cudaFree(p->images_h[c]); if(c<(int)p->clt_h.size())cudaFree(p->clt_h[c]); if(c<(int)p->clt_ref_h.size())cudaFree(p->clt_ref_h[c]); } if(c<(int)p->images_h.size())cudaFree(p->images_h[c]); if(c<(int)p->clt_h.size())cudaFree(p->clt_h[c]); if(c<(int)p->clt_ref_h.size())cudaFree(p->clt_ref_h[c]); }
cudaFree(p->gpu_kernels); cudaFree(p->gpu_kernel_offsets); cudaFree(p->gpu_images); cudaFree(p->gpu_clt); cudaFree(p->gpu_clt_ref); cudaFree(p->gpu_kernels); cudaFree(p->gpu_kernel_offsets); cudaFree(p->gpu_images); cudaFree(p->gpu_clt); cudaFree(p->gpu_clt_ref);
cudaFree(p->cond_scale_maps); cudaFree(p->cond_offset_maps);
cudaFree(p->cond_scales); cudaFree(p->cond_offsets); cudaFree(p->cond_scales2); cudaFree(p->cond_fpn); // By Codex on 07/15/2026
cudaFree(p->active); cudaFree(p->num_active); cudaFree(p->active); cudaFree(p->num_active);
cudaFree(p->task_slots[0]); cudaFree(p->task_slots[1]); // p->ftasks aliases one of these // By Codex on 07/15/2026 cudaFree(p->task_slots[0]); cudaFree(p->task_slots[1]); // p->ftasks aliases one of these // By Codex on 07/15/2026
cudaFree(p->gc); cudaFree(p->cv); cudaFree(p->rbr); cudaFree(p->rot); cudaFree(p->gc); cudaFree(p->cv); cudaFree(p->rbr); cudaFree(p->rot);
......
...@@ -25,7 +25,7 @@ fi ...@@ -25,7 +25,7 @@ fi
# Kernel + host sources shared by every test ("load all kernels"): # Kernel + host sources shared by every test ("load all kernels"):
COMMON="$SRC/TileProcessor.cu $SRC/dtt8x8.cu $SRC/geometry_correction.cu \ COMMON="$SRC/TileProcessor.cu $SRC/dtt8x8.cu $SRC/geometry_correction.cu \
$SRC/tp_consolidate.cu $SRC/tp_peak.cu $SRC/tp_lma.cu \ $SRC/tp_condition.cu $SRC/tp_consolidate.cu $SRC/tp_peak.cu $SRC/tp_lma.cu \
$SRC/tp_utils.cu $SRC/tp_files.cu tp_test_data.cu" $SRC/tp_utils.cu $SRC/tp_files.cu tp_test_data.cu"
TESTS="${@:-$(ls test_*.cu | sed 's/\.cu$//')}" TESTS="${@:-$(ls test_*.cu | sed 's/\.cu$//')}"
......
/*
* test_condition.cu
*
* Standalone synthetic gate for the LWIR conditioning kernels. It checks the
* double validation tier bit-for-bit against the Java expression order and the
* production two-map tier bit-for-bit against host fmaf(), including pitched
* images and dimensions that are not block multiples.
*
* Created on: Jul 15, 2026
* Author: Codex (OpenAI), for Elphel
*/
#include <cstdio>
#include <cstdlib>
#include <cmath>
#include <cstdint>
#include <cstring>
#include <vector>
#include <cuda_runtime.h>
#include <helper_cuda.h>
#include "tp_condition.h"
static uint32_t bits(float v)
{
uint32_t u;
memcpy(&u, &v, sizeof(u));
return u;
}
int main(int argc, char ** argv)
{
findCudaDevice(argc, (const char **) argv);
const int cams = 3;
const int width = 37;
const int height = 19;
const int pixels = width * height;
const double dc = -3512.375;
std::vector<std::vector<float> > raw(cams, std::vector<float>(pixels));
std::vector<double> scales(cams), offsets(cams), scales2(cams, 0.0);
std::vector<double> fpn((size_t) cams * pixels);
std::vector<float> scale_maps((size_t) cams * pixels);
std::vector<float> offset_maps((size_t) cams * pixels);
std::vector<std::vector<float> > expected_d = raw;
std::vector<std::vector<float> > expected_f = raw;
for (int c = 0; c < cams; c++) {
scales[c] = 0.985 + 0.013 * c;
offsets[c] = 20117.25 + 9.75 * c;
for (int i = 0; i < pixels; i++) {
const size_t k = (size_t) c * pixels + i;
raw[c][i] = (float) (20750.0 + 0.125 * (i % 97) - 3.5 * c);
fpn[k] = 7.0 * sin(0.017 * i + 0.3 * c) + 0.03125 * (i % 11);
const double d = (double) raw[c][i] - offsets[c];
expected_d[c][i] = (float) (scales2[c] * d * d + scales[c] * d - fpn[k] - dc);
scale_maps[k] = (float) scales[c];
offset_maps[k] = (float) (scales[c] * offsets[c] + fpn[k] + dc);
expected_f[c][i] = std::fmaf(raw[c][i], scale_maps[k], -offset_maps[k]);
}
}
std::vector<float *> images(cams, nullptr);
size_t pitch = 0;
for (int c = 0; c < cams; c++) {
size_t this_pitch = 0;
checkCudaErrors(cudaMallocPitch((void **) &images[c], &this_pitch,
(size_t) width * sizeof(float), height));
if (c == 0) pitch = this_pitch;
if (this_pitch != pitch) {
fprintf(stderr, "inconsistent pitches\n");
return EXIT_FAILURE;
}
}
float ** gpu_images = nullptr;
float * gpu_scale_maps = nullptr, * gpu_offset_maps = nullptr;
double * gpu_scales = nullptr, * gpu_offsets = nullptr, * gpu_scales2 = nullptr, * gpu_fpn = nullptr;
checkCudaErrors(cudaMalloc((void **) &gpu_images, cams * sizeof(float *)));
checkCudaErrors(cudaMemcpy(gpu_images, images.data(), cams * sizeof(float *), cudaMemcpyHostToDevice));
checkCudaErrors(cudaMalloc((void **) &gpu_scale_maps, scale_maps.size() * sizeof(float)));
checkCudaErrors(cudaMalloc((void **) &gpu_offset_maps, offset_maps.size() * sizeof(float)));
checkCudaErrors(cudaMalloc((void **) &gpu_scales, scales.size() * sizeof(double)));
checkCudaErrors(cudaMalloc((void **) &gpu_offsets, offsets.size() * sizeof(double)));
checkCudaErrors(cudaMalloc((void **) &gpu_scales2, scales2.size() * sizeof(double)));
checkCudaErrors(cudaMalloc((void **) &gpu_fpn, fpn.size() * sizeof(double)));
checkCudaErrors(cudaMemcpy(gpu_scale_maps, scale_maps.data(), scale_maps.size() * sizeof(float), cudaMemcpyHostToDevice));
checkCudaErrors(cudaMemcpy(gpu_offset_maps, offset_maps.data(), offset_maps.size() * sizeof(float), cudaMemcpyHostToDevice));
checkCudaErrors(cudaMemcpy(gpu_scales, scales.data(), scales.size() * sizeof(double), cudaMemcpyHostToDevice));
checkCudaErrors(cudaMemcpy(gpu_offsets, offsets.data(), offsets.size() * sizeof(double), cudaMemcpyHostToDevice));
checkCudaErrors(cudaMemcpy(gpu_scales2, scales2.data(), scales2.size() * sizeof(double), cudaMemcpyHostToDevice));
checkCudaErrors(cudaMemcpy(gpu_fpn, fpn.data(), fpn.size() * sizeof(double), cudaMemcpyHostToDevice));
const dim3 block(32, 8, 1);
const dim3 grid((width + block.x - 1) / block.x,
(height + block.y - 1) / block.y, cams);
auto upload = [&]() {
for (int c = 0; c < cams; c++) checkCudaErrors(cudaMemcpy2D(
images[c], pitch, raw[c].data(), (size_t) width * sizeof(float),
(size_t) width * sizeof(float), height, cudaMemcpyHostToDevice));
};
auto compare = [&](const std::vector<std::vector<float> > & expected, const char * name) {
long bad = 0;
double max_diff = 0.0;
std::vector<float> got(pixels);
for (int c = 0; c < cams; c++) {
checkCudaErrors(cudaMemcpy2D(got.data(), (size_t) width * sizeof(float),
images[c], pitch, (size_t) width * sizeof(float), height, cudaMemcpyDeviceToHost));
for (int i = 0; i < pixels; i++) {
if (bits(got[i]) != bits(expected[c][i])) bad++;
const double d = fabs((double) got[i] - expected[c][i]);
if (d > max_diff) max_diff = d;
}
}
printf("%s: %d pixels, bit mismatches=%ld max|diff|=%g -> %s\n",
name, cams * pixels, bad, max_diff, bad ? "FAIL" : "PASS");
return bad;
};
upload();
condition_lwir_d<<<grid, block>>>(cams, width, height, pitch / sizeof(float), gpu_images,
gpu_scales, gpu_offsets, gpu_scales2, gpu_fpn, dc);
checkCudaErrors(cudaGetLastError());
checkCudaErrors(cudaDeviceSynchronize());
long bad = compare(expected_d, "condition_lwir_d");
upload();
condition_lwir_f<<<grid, block>>>(cams, width, height, pitch / sizeof(float), gpu_images,
gpu_scale_maps, gpu_offset_maps);
checkCudaErrors(cudaGetLastError());
checkCudaErrors(cudaDeviceSynchronize());
bad += compare(expected_f, "condition_lwir_f");
double max_model_diff = 0.0;
for (int c = 0; c < cams; c++) for (int i = 0; i < pixels; i++) {
const double d = fabs((double) expected_f[c][i] - expected_d[c][i]);
if (d > max_model_diff) max_model_diff = d;
}
printf("float two-map vs double Java model: max|diff|=%g counts\n", max_model_diff);
for (float * p : images) cudaFree(p);
cudaFree(gpu_images); cudaFree(gpu_scale_maps); cudaFree(gpu_offset_maps);
cudaFree(gpu_scales); cudaFree(gpu_offsets); cudaFree(gpu_scales2); cudaFree(gpu_fpn);
return bad ? EXIT_FAILURE : EXIT_SUCCESS;
}
/*
* tp_condition.cu - LWIR scale/offset conditioning kernels.
* See tp_condition.h for the frozen contracts.
*
* Created on: Jul 15, 2026
* Author: Codex (OpenAI), for Elphel
*/
#ifndef JCUDA // NVRTC concatenation already placed tp_condition.h before this file
#include "tp_condition.h"
#endif
extern "C" __global__ void condition_lwir_f(
int num_cams,
int width,
int height,
size_t image_stride,
float ** gpu_images,
const float * scale_maps,
const float * offset_maps)
{
const int x = (int) (blockIdx.x * blockDim.x + threadIdx.x);
const int y = (int) (blockIdx.y * blockDim.y + threadIdx.y);
const int c = (int) blockIdx.z;
if ((c >= num_cams) || (x >= width) || (y >= height)) return;
const size_t pixel = (size_t) y * width + x;
const size_t map_index = (size_t) c * width * height + pixel;
float * image = gpu_images[c];
image[(size_t) y * image_stride + x] = fmaf(
image[(size_t) y * image_stride + x],
scale_maps[map_index],
-offset_maps[map_index]);
}
extern "C" __global__ void condition_lwir_d(
int num_cams,
int width,
int height,
size_t image_stride,
float ** gpu_images,
const double * scales,
const double * offsets,
const double * scales2,
const double * fpn,
double dc_level)
{
const int x = (int) (blockIdx.x * blockDim.x + threadIdx.x);
const int y = (int) (blockIdx.y * blockDim.y + threadIdx.y);
const int c = (int) blockIdx.z;
if ((c >= num_cams) || (x >= width) || (y >= height)) return;
const size_t pixel = (size_t) y * width + x;
const size_t map_index = (size_t) c * width * height + pixel;
float * image = gpu_images[c];
// Explicit round-to-nearest operations prevent NVRTC from fusing across the
// Java statement boundaries. This kernel is a validation oracle, not the
// production arithmetic path.
const double raw = (double) image[(size_t) y * image_stride + x];
const double d = __dsub_rn(raw, offsets[c]);
double value;
if (scales2[c] == 0.0) {
value = __dmul_rn(scales[c], d);
} else {
const double quadratic = __dmul_rn(__dmul_rn(scales2[c], d), d);
value = __dadd_rn(quadratic, __dmul_rn(scales[c], d));
}
value = __dsub_rn(value, fpn[map_index]);
value = __dsub_rn(value, dc_level);
image[(size_t) y * image_stride + x] = (float) value;
}
/*
* tp_condition.h
*
* LWIR input conditioning on the resident TpProc image buffers. Row/column
* denoising is intentionally absent: production Boson cameras perform that
* correction internally. The production contract is the two-map model
*
* output = raw * scale_map - offset_map
*
* where offset_map includes scale*sensor_offset, per-pixel FPN, and the fixed
* per-sequence DC level. condition_lwir_d is the validation tier: it mirrors
* CuasConditioning.condition()'s double expression and rounds only when the
* result is stored back to the float image.
*
* Created on: Jul 15, 2026
* Author: Codex (OpenAI), for Elphel
*/
#ifndef SRC_TP_CONDITION_H_
#define SRC_TP_CONDITION_H_
#ifndef JCUDA
#include <cuda_runtime.h>
#endif
extern "C" __global__ void condition_lwir_f(
int num_cams,
int width,
int height,
size_t image_stride,
float ** gpu_images,
const float * scale_maps,
const float * offset_maps);
extern "C" __global__ void condition_lwir_d(
int num_cams,
int width,
int height,
size_t image_stride,
float ** gpu_images,
const double * scales,
const double * offsets,
const double * scales2,
const double * fpn,
double dc_level);
#endif /* SRC_TP_CONDITION_H_ */
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