diff --git a/cxx/isce3/cuda/matchtemplate/pycuampcor/cuAmpcorParameter.cu b/cxx/isce3/cuda/matchtemplate/pycuampcor/cuAmpcorParameter.cu index 6e058df0c..4ac88bf73 100644 --- a/cxx/isce3/cuda/matchtemplate/pycuampcor/cuAmpcorParameter.cu +++ b/cxx/isce3/cuda/matchtemplate/pycuampcor/cuAmpcorParameter.cu @@ -88,8 +88,8 @@ void cuAmpcorParameter::setupParameters() windowSizeWidth = windowSizeWidthRaw*rawDataOversamplingFactor; // windowSizeHeight = windowSizeHeightRaw*rawDataOversamplingFactor; - searchWindowSizeWidthRaw = windowSizeWidthRaw + 2*halfSearchRangeDownRaw; - searchWindowSizeHeightRaw = windowSizeHeightRaw + 2*halfSearchRangeAcrossRaw; + searchWindowSizeWidthRaw = windowSizeWidthRaw + 2*halfSearchRangeAcrossRaw; + searchWindowSizeHeightRaw = windowSizeHeightRaw + 2*halfSearchRangeDownRaw; searchWindowSizeWidthRawZoomIn = windowSizeWidthRaw + 2*halfZoomWindowSizeRaw; searchWindowSizeHeightRawZoomIn = windowSizeHeightRaw + 2*halfZoomWindowSizeRaw; @@ -237,14 +237,28 @@ void cuAmpcorParameter::setChunkStartPixels() if(sChunkEA < vpixel) sChunkEA = vpixel; } } + // Clamp chunk start pixels to valid image range + if (mChunkSD < 0) mChunkSD = 0; + if (mChunkSA < 0) mChunkSA = 0; + if (sChunkSD < 0) sChunkSD = 0; + if (sChunkSA < 0) sChunkSA = 0; + // Compute end pixels and clamp to image range + int mEndD = mChunkED + windowSizeHeightRaw; + if (mEndD > referenceImageHeight) mEndD = referenceImageHeight; + int mEndA = mChunkEA + windowSizeWidthRaw; + if (mEndA > referenceImageWidth) mEndA = referenceImageWidth; + int sEndD = sChunkED + searchWindowSizeHeightRaw; + if (sEndD > secondaryImageHeight) sEndD = secondaryImageHeight; + int sEndA = sChunkEA + searchWindowSizeWidthRaw; + if (sEndA > secondaryImageWidth) sEndA = secondaryImageWidth; referenceChunkStartPixelDown[idxChunk] = mChunkSD; referenceChunkStartPixelAcross[idxChunk] = mChunkSA; secondaryChunkStartPixelDown[idxChunk] = sChunkSD; secondaryChunkStartPixelAcross[idxChunk] = sChunkSA; - referenceChunkHeight[idxChunk] = mChunkED - mChunkSD + windowSizeHeightRaw; - referenceChunkWidth[idxChunk] = mChunkEA - mChunkSA + windowSizeWidthRaw; - secondaryChunkHeight[idxChunk] = sChunkED - sChunkSD + searchWindowSizeHeightRaw; - secondaryChunkWidth[idxChunk] = sChunkEA - sChunkSA + searchWindowSizeWidthRaw; + referenceChunkHeight[idxChunk] = mEndD - mChunkSD; + referenceChunkWidth[idxChunk] = mEndA - mChunkSA; + secondaryChunkHeight[idxChunk] = sEndD - sChunkSD; + secondaryChunkWidth[idxChunk] = sEndA - sChunkSA; if(maxReferenceChunkHeight < referenceChunkHeight[idxChunk]) maxReferenceChunkHeight = referenceChunkHeight[idxChunk]; if(maxReferenceChunkWidth < referenceChunkWidth[idxChunk] ) maxReferenceChunkWidth = referenceChunkWidth[idxChunk]; if(maxSecondaryChunkHeight < secondaryChunkHeight[idxChunk]) maxSecondaryChunkHeight = secondaryChunkHeight[idxChunk]; @@ -256,60 +270,37 @@ void cuAmpcorParameter::setChunkStartPixels() /// check whether reference and secondary windows are within the image range void cuAmpcorParameter::checkPixelInImageRange() { +#ifdef CUAMPCOR_DEBUG int endPixel; for(int row=0; row= referenceImageHeight) - { - fprintf(stderr, "Reference Window end pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, endPixel); - exit(EXIT_FAILURE); - } + fprintf(stderr, "Warning: Reference Window end pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, endPixel); endPixel = referenceStartPixelAcross[i] + windowSizeWidthRaw; if(endPixel >= referenceImageWidth) - { - fprintf(stderr, "Reference Window end pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, endPixel); - exit(EXIT_FAILURE); - } + fprintf(stderr, "Warning: Reference Window end pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, endPixel); //secondary - if(secondaryStartPixelDown[i] <0) - { - fprintf(stderr, "Secondary Window start pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, secondaryStartPixelDown[i]); - exit(EXIT_FAILURE); - } - if(secondaryStartPixelAcross[i] <0) - { - fprintf(stderr, "Secondary Window start pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, secondaryStartPixelAcross[i]); - exit(EXIT_FAILURE); - } + if(secondaryStartPixelDown[i] < 0) + fprintf(stderr, "Warning: Secondary Window start pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, secondaryStartPixelDown[i]); + if(secondaryStartPixelAcross[i] < 0) + fprintf(stderr, "Warning: Secondary Window start pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, secondaryStartPixelAcross[i]); endPixel = secondaryStartPixelDown[i] + searchWindowSizeHeightRaw; if(endPixel >= secondaryImageHeight) - { - fprintf(stderr, "Secondary Window end pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, endPixel); - exit(EXIT_FAILURE); - } + fprintf(stderr, "Warning: Secondary Window end pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, endPixel); endPixel = secondaryStartPixelAcross[i] + searchWindowSizeWidthRaw; if(endPixel >= secondaryImageWidth) - { - fprintf(stderr, "Secondary Window end pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, endPixel); - exit(EXIT_FAILURE); - } - + fprintf(stderr, "Warning: Secondary Window end pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, endPixel); } } +#endif } diff --git a/cxx/isce3/cuda/matchtemplate/pycuampcor/cuArraysCopy.cu b/cxx/isce3/cuda/matchtemplate/pycuampcor/cuArraysCopy.cu index d9e36461b..3e80f2808 100644 --- a/cxx/isce3/cuda/matchtemplate/pycuampcor/cuArraysCopy.cu +++ b/cxx/isce3/cuda/matchtemplate/pycuampcor/cuArraysCopy.cu @@ -68,7 +68,7 @@ void cuArraysCopyToBatch(cuArrays *image1, cuArrays *image2, // kernel for cuArraysCopyToBatchWithOffset template -__global__ void cuArraysCopyToBatchWithOffset_kernel(const T_in *imageIn, const int inNY, +__global__ void cuArraysCopyToBatchWithOffset_kernel(const T_in *imageIn, const int inNX, const int inNY, T_out *imageOut, const int outNX, const int outNY, const int nImages, const int *offsetX, const int *offsetY) { @@ -77,8 +77,12 @@ __global__ void cuArraysCopyToBatchWithOffset_kernel(const T_in *imageIn, const int outy = threadIdx.y + blockDim.y*blockIdx.y; if(idxImage>=nImages || outx >= outNX || outy >= outNY) return; int idxOut = idxImage*outNX*outNY + outx*outNY + outy; - int idxIn = (offsetX[idxImage]+outx)*inNY + offsetY[idxImage] + outy; - imageOut[idxOut] = T_out{imageIn[idxIn]}; + int inx = offsetX[idxImage] + outx; + int iny = offsetY[idxImage] + outy; + if (inx >= 0 && inx < inNX && iny >= 0 && iny < inNY) + imageOut[idxOut] = T_out{imageIn[inx*inNY + iny]}; + else + imageOut[idxOut] = T_out{0.0f}; } /** @@ -98,14 +102,14 @@ void cuArraysCopyToBatchWithOffset(cuArrays *image1, const int lda1, cuA dim3 blockSize(nthreads, nthreads, 1); dim3 gridSize(IDIVUP(image2->height,nthreads), IDIVUP(image2->width,nthreads), image2->count); cuArraysCopyToBatchWithOffset_kernel<<>> ( - image1->devData, lda1, + image1->devData, image1->height, lda1, image2->devData, image2->height, image2->width, image2->count, offsetH, offsetW); getLastCudaError("cuArraysCopyToBatchAbsWithOffset_kernel"); } // same as above, but from complex to real(take amplitudes) -__global__ void cuArraysCopyToBatchAbsWithOffset_kernel(const float2 *imageIn, const int inNY, +__global__ void cuArraysCopyToBatchAbsWithOffset_kernel(const float2 *imageIn, const int inNX, const int inNY, float2 *imageOut, const int outNX, const int outNY, const int nImages, const int *offsetX, const int *offsetY) { @@ -114,8 +118,12 @@ __global__ void cuArraysCopyToBatchAbsWithOffset_kernel(const float2 *imageIn, c int outy = threadIdx.y + blockDim.y*blockIdx.y; if(idxImage>=nImages || outx >= outNX || outy >= outNY) return; int idxOut = idxImage*outNX*outNY + outx*outNY + outy; - int idxIn = (offsetX[idxImage]+outx)*inNY + offsetY[idxImage] + outy; - imageOut[idxOut] = make_float2(complexAbs(imageIn[idxIn]), 0.0); + int inx = offsetX[idxImage] + outx; + int iny = offsetY[idxImage] + outy; + if (inx >= 0 && inx < inNX && iny >= 0 && iny < inNY) + imageOut[idxOut] = make_float2(complexAbs(imageIn[inx*inNY + iny]), 0.0f); + else + imageOut[idxOut] = make_float2(0.0f, 0.0f); } /** @@ -135,7 +143,7 @@ void cuArraysCopyToBatchAbsWithOffset(cuArrays *image1, const int lda1, dim3 blockSize(nthreads, nthreads, 1); dim3 gridSize(IDIVUP(image2->height,nthreads), IDIVUP(image2->width,nthreads), image2->count); cuArraysCopyToBatchAbsWithOffset_kernel<<>> ( - image1->devData, lda1, + image1->devData, image1->height, lda1, image2->devData, image2->height, image2->width, image2->count, offsetH, offsetW); getLastCudaError("cuArraysCopyToBatchAbsWithOffset_kernel"); @@ -158,7 +166,7 @@ void cuArraysCopyToBatchWithOffsetR2C(cuArrays *image1, const int lda1, c dim3 blockSize(nthreads, nthreads, 1); dim3 gridSize(IDIVUP(image2->height,nthreads), IDIVUP(image2->width,nthreads), image2->count); cuArraysCopyToBatchWithOffset_kernel<<>> ( - image1->devData, lda1, + image1->devData, image1->height, lda1, image2->devData, image2->height, image2->width, image2->count, offsetH, offsetW); getLastCudaError("cuArraysCopyToBatchWithOffsetR2C_kernel"); diff --git a/cxx/isce3/matchtemplate/pycuampcor/cuAmpcorParameter.cpp b/cxx/isce3/matchtemplate/pycuampcor/cuAmpcorParameter.cpp index b82313fc9..c28b39880 100644 --- a/cxx/isce3/matchtemplate/pycuampcor/cuAmpcorParameter.cpp +++ b/cxx/isce3/matchtemplate/pycuampcor/cuAmpcorParameter.cpp @@ -90,8 +90,8 @@ void cuAmpcorParameter::setupParameters() windowSizeWidth = windowSizeWidthRaw*rawDataOversamplingFactor; // windowSizeHeight = windowSizeHeightRaw*rawDataOversamplingFactor; - searchWindowSizeWidthRaw = windowSizeWidthRaw + 2*halfSearchRangeDownRaw; - searchWindowSizeHeightRaw = windowSizeHeightRaw + 2*halfSearchRangeAcrossRaw; + searchWindowSizeWidthRaw = windowSizeWidthRaw + 2*halfSearchRangeAcrossRaw; + searchWindowSizeHeightRaw = windowSizeHeightRaw + 2*halfSearchRangeDownRaw; searchWindowSizeWidthRawZoomIn = windowSizeWidthRaw + 2*halfZoomWindowSizeRaw; searchWindowSizeHeightRawZoomIn = windowSizeHeightRaw + 2*halfZoomWindowSizeRaw; @@ -239,14 +239,28 @@ void cuAmpcorParameter::setChunkStartPixels() if(sChunkEA < vpixel) sChunkEA = vpixel; } } + // Clamp chunk start pixels to valid image range + if (mChunkSD < 0) mChunkSD = 0; + if (mChunkSA < 0) mChunkSA = 0; + if (sChunkSD < 0) sChunkSD = 0; + if (sChunkSA < 0) sChunkSA = 0; + // Compute end pixels and clamp to image range + int mEndD = mChunkED + windowSizeHeightRaw; + if (mEndD > referenceImageHeight) mEndD = referenceImageHeight; + int mEndA = mChunkEA + windowSizeWidthRaw; + if (mEndA > referenceImageWidth) mEndA = referenceImageWidth; + int sEndD = sChunkED + searchWindowSizeHeightRaw; + if (sEndD > secondaryImageHeight) sEndD = secondaryImageHeight; + int sEndA = sChunkEA + searchWindowSizeWidthRaw; + if (sEndA > secondaryImageWidth) sEndA = secondaryImageWidth; referenceChunkStartPixelDown[idxChunk] = mChunkSD; referenceChunkStartPixelAcross[idxChunk] = mChunkSA; secondaryChunkStartPixelDown[idxChunk] = sChunkSD; secondaryChunkStartPixelAcross[idxChunk] = sChunkSA; - referenceChunkHeight[idxChunk] = mChunkED - mChunkSD + windowSizeHeightRaw; - referenceChunkWidth[idxChunk] = mChunkEA - mChunkSA + windowSizeWidthRaw; - secondaryChunkHeight[idxChunk] = sChunkED - sChunkSD + searchWindowSizeHeightRaw; - secondaryChunkWidth[idxChunk] = sChunkEA - sChunkSA + searchWindowSizeWidthRaw; + referenceChunkHeight[idxChunk] = mEndD - mChunkSD; + referenceChunkWidth[idxChunk] = mEndA - mChunkSA; + secondaryChunkHeight[idxChunk] = sEndD - sChunkSD; + secondaryChunkWidth[idxChunk] = sEndA - sChunkSA; if(maxReferenceChunkHeight < referenceChunkHeight[idxChunk]) maxReferenceChunkHeight = referenceChunkHeight[idxChunk]; if(maxReferenceChunkWidth < referenceChunkWidth[idxChunk] ) maxReferenceChunkWidth = referenceChunkWidth[idxChunk]; if(maxSecondaryChunkHeight < secondaryChunkHeight[idxChunk]) maxSecondaryChunkHeight = secondaryChunkHeight[idxChunk]; @@ -258,60 +272,37 @@ void cuAmpcorParameter::setChunkStartPixels() /// check whether reference and secondary windows are within the image range void cuAmpcorParameter::checkPixelInImageRange() { +#ifdef CUAMPCOR_DEBUG int endPixel; for(int row=0; row= referenceImageHeight) - { - fprintf(stderr, "Reference Window end pixel out ot range in Down, window (%d,%d), pixel %d\n", row, col, endPixel); - exit(EXIT_FAILURE); - } + fprintf(stderr, "Warning: Reference Window end pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, endPixel); endPixel = referenceStartPixelAcross[i] + windowSizeWidthRaw; if(endPixel >= referenceImageWidth) - { - fprintf(stderr, "Reference Window end pixel out ot range in Across, window (%d,%d), pixel %d\n", row, col, endPixel); - exit(EXIT_FAILURE); - } + fprintf(stderr, "Warning: Reference Window end pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, endPixel); //secondary - if(secondaryStartPixelDown[i] <0) - { - fprintf(stderr, "Secondary Window start pixel out ot range in Down, window (%d,%d), pixel %d\n", row, col, secondaryStartPixelDown[i]); - exit(EXIT_FAILURE); - } - if(secondaryStartPixelAcross[i] <0) - { - fprintf(stderr, "Secondary Window start pixel out ot range in Across, window (%d,%d), pixel %d\n", row, col, secondaryStartPixelAcross[i]); - exit(EXIT_FAILURE); - } + if(secondaryStartPixelDown[i] < 0) + fprintf(stderr, "Warning: Secondary Window start pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, secondaryStartPixelDown[i]); + if(secondaryStartPixelAcross[i] < 0) + fprintf(stderr, "Warning: Secondary Window start pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, secondaryStartPixelAcross[i]); endPixel = secondaryStartPixelDown[i] + searchWindowSizeHeightRaw; if(endPixel >= secondaryImageHeight) - { - fprintf(stderr, "Secondary Window end pixel out ot range in Down, window (%d,%d), pixel %d\n", row, col, endPixel); - exit(EXIT_FAILURE); - } + fprintf(stderr, "Warning: Secondary Window end pixel out of range in Down, window (%d,%d), pixel %d\n", row, col, endPixel); endPixel = secondaryStartPixelAcross[i] + searchWindowSizeWidthRaw; if(endPixel >= secondaryImageWidth) - { - fprintf(stderr, "Secondary Window end pixel out ot range in Across, window (%d,%d), pixel %d\n", row, col, endPixel); - exit(EXIT_FAILURE); - } - + fprintf(stderr, "Warning: Secondary Window end pixel out of range in Across, window (%d,%d), pixel %d\n", row, col, endPixel); } } +#endif } diff --git a/cxx/isce3/matchtemplate/pycuampcor/cuArraysCopy.cpp b/cxx/isce3/matchtemplate/pycuampcor/cuArraysCopy.cpp index 0a7949b71..daebb2ac3 100644 --- a/cxx/isce3/matchtemplate/pycuampcor/cuArraysCopy.cpp +++ b/cxx/isce3/matchtemplate/pycuampcor/cuArraysCopy.cpp @@ -47,12 +47,17 @@ void cuArraysCopyToBatchWithOffset(cuArrays *image1, const int lda1, cuA const int outNX = image2->height; const int outNY = image2->width; + const int inNX = image1->height; for (int idxImage = 0; idxImage < image2->count; idxImage++) { for (int outx = 0; outx < outNX; outx++) { for (int outy = 0; outy < outNY; outy++) { int idxOut = idxImage*outNX*outNY + outx*outNY + outy; - int idxIn = (offsetX[idxImage]+outx)*inNY + offsetY[idxImage] + outy; - imageOut[idxOut] = imageIn[idxIn]; + int inx = offsetX[idxImage] + outx; + int iny = offsetY[idxImage] + outy; + if (inx >= 0 && inx < inNX && iny >= 0 && iny < inNY) + imageOut[idxOut] = imageIn[inx*inNY + iny]; + else + imageOut[idxOut] = make_float2(0.0f, 0.0f); } } } @@ -77,12 +82,17 @@ void cuArraysCopyToBatchAbsWithOffset(cuArrays *image1, const int lda1, const float2* imageIn = image1->devData; const int inNY = lda1; float2* imageOut = image2->devData; + const int inNX = image1->height; for (int idxImage = 0; idxImage < image2->count; idxImage++) { for (int outx = 0; outx < outNX; outx++) { for (int outy = 0; outy < outNY; outy++) { int idxOut = idxImage*outNX*outNY + outx*outNY + outy; - int idxIn = (offsetX[idxImage]+outx)*inNY + offsetY[idxImage] + outy; - imageOut[idxOut] = make_float2(complexAbs(imageIn[idxIn]), 0.0); + int inx = offsetX[idxImage] + outx; + int iny = offsetY[idxImage] + outy; + if (inx >= 0 && inx < inNX && iny >= 0 && iny < inNY) + imageOut[idxOut] = make_float2(complexAbs(imageIn[inx*inNY + iny]), 0.0f); + else + imageOut[idxOut] = make_float2(0.0f, 0.0f); } } } diff --git a/cxx/isce3/matchtemplate/pycuampcor/cuCorrNormalizer.h b/cxx/isce3/matchtemplate/pycuampcor/cuCorrNormalizer.h index 94f9ff6ff..b67e1f09a 100644 --- a/cxx/isce3/matchtemplate/pycuampcor/cuCorrNormalizer.h +++ b/cxx/isce3/matchtemplate/pycuampcor/cuCorrNormalizer.h @@ -22,9 +22,8 @@ namespace isce3::matchtemplate::pycuampcor { */ class cuNormalizeProcessor { public: - // default constructor and destructor - cuNormalizeProcessor() {} - ~cuNormalizeProcessor() {} + cuNormalizeProcessor() = default; + virtual ~cuNormalizeProcessor() = default; // execute interface virtual void execute(cuArrays * correlation, cuArrays *reference, cuArrays *secondary) = 0; }; diff --git a/python/packages/nisar/workflows/dense_offsets.py b/python/packages/nisar/workflows/dense_offsets.py index 4a1feddd3..f10438e1a 100644 --- a/python/packages/nisar/workflows/dense_offsets.py +++ b/python/packages/nisar/workflows/dense_offsets.py @@ -178,12 +178,22 @@ def set_optional_attributes(ampcor_obj, cfg, length, width): ampcor_obj.numberWindowAcross = cfg['offset_width'] else: offset_width = (width - margin_rg) // ampcor_obj.skipSampleAcross + if offset_width <= 0: + err_str = (f"Image width ({width}) is too small for the configured parameters " + f"(margin_rg={margin_rg}): no valid ampcor windows fit across range.") + error_channel.log(err_str) + raise ValueError(err_str) ampcor_obj.numberWindowAcross = offset_width if cfg['offset_length'] is not None: ampcor_obj.numberWindowDown = cfg['offset_length'] else: offset_length = (length - margin_az) // ampcor_obj.skipSampleDown + if offset_length <= 0: + err_str = (f"Image length ({length}) is too small for the configured parameters " + f"(margin_az={margin_az}): no valid ampcor windows fit along azimuth.") + error_channel.log(err_str) + raise ValueError(err_str) ampcor_obj.numberWindowDown = offset_length if cfg['cross_correlation_domain'] is not None: @@ -266,7 +276,7 @@ def set_optional_attributes(ampcor_obj, cfg, length, width): if cfg['merge_gross_offset'] is not None: ampcor_obj.mergeGrossOffset = 1 if cfg['merge_gross_offset'] else 0 - # Check pixel in image range + # Check pixel in image range; warns to stderr when out of range but does not abort ampcor_obj.checkPixelInImageRange() return ampcor_obj diff --git a/share/nisar/schemas/insar.yaml b/share/nisar/schemas/insar.yaml index bbeab4faa..8006d324a 100644 --- a/share/nisar/schemas/insar.yaml +++ b/share/nisar/schemas/insar.yaml @@ -312,10 +312,10 @@ dense_offsets_options: window_azimuth: int(min=16, max=512, required=False) # Number of columns for search chip/template in secondary image - half_search_range: int(min=4, max=128, required=False) + half_search_range: int(min=4, max=256, required=False) # Number of lines for search chip/template in secondary image - half_search_azimuth: int(min=4, max=128, required=False) + half_search_azimuth: int(min=4, max=256, required=False) # Number of columns to skip in reference image for next offset estimate skip_range: int(min=1, required=False) @@ -487,10 +487,10 @@ offset_layer_options: window_azimuth: int(min=16, max=512, required=False) # Number of columns for search chip/template in secondary image - half_search_range: int(min=4, max=64, required=False) + half_search_range: int(min=4, max=256, required=False) # Number of lines for search chip/template in secondary image - half_search_azimuth: int(min=4, max=64, required=False) + half_search_azimuth: int(min=4, max=256, required=False) rubbersheet_options: # Path to dense offsets outputs (offsets, snr, covariance).