多云时无法配准

cmd/log.go

@@ -49,7 +49,7 @@ func init() {
 		DisableColors:   true,
 	}
 
-	configureLogger(logrus.StandardLogger(), "log/SJY01-imgproc.log", logrus.InfoLevel)
+	configureLogger(logrus.StandardLogger(), "log/SJY01-imgproc.log", logrus.DebugLevel)
 }
 
 func NewLogger(logfile string) *logrus.Logger {

cmd/proc.go

@@ -6,13 +6,14 @@ import (
 	"strings"
 
 	"github.com/airbusgeo/godal"
-	log "github.com/sirupsen/logrus"
+	"github.com/sirupsen/logrus"
 	"github.com/spf13/cobra"
 	producer "starwiz.cn/sjy01/image-proc/producer"
 )
 
 var (
 	params    producer.Params
+	saveStrip bool
 )
 
 var procCmd = &cobra.Command{
@@ -52,12 +53,15 @@ var procCmd = &cobra.Command{
 
 		panScenes, mssScenes, err := reg.SubScenes()
 		if err != nil {
-			log.Error(err)
+			logrus.Error(err)
 		}
 
 		reg.SaveScenesToTiff(panScenes, mssScenes)
-		// reg.SaveOriginalPanToGDALGTiff(reg.Params.PanTiffFile)
+
-		// reg.SaveRegisteredMssToGDALGTiff(reg.Params.MssTiffFile)
+		if saveStrip {
+			reg.SaveOriginalPanToGDALGTiff(reg.Params.PanTiffFile)
+			reg.SaveRegisteredMssToGDALGTiff(reg.Params.MssTiffFile)
+		}
 
 		if params.PansharpenIHS {
 			reg.DoScenePansharpen(panScenes, mssScenes)
@@ -77,4 +81,5 @@ func init() {
 	procCmd.Flags().BoolVarP(&params.PansharpenIHS, "fus", "", false, "pansharpen using IHS")
 	procCmd.Flags().StringVarP(&params.OutputDir, "output-dir", "o", "data", "output directory")
 	procCmd.Flags().BoolVarP(&params.SubScenes, "sub-scenes", "", false, "process sub-scenes")
+	procCmd.Flags().BoolVarP(&saveStrip, "save-strip", "", false, "save original and registered images as GDAL GTiff")
 }

producer/gdal_pansharpen.go

@@ -86,7 +86,7 @@ func GDALPansharpen(panTiff, mssTiff string) (string, error) {
 	fusedTIff := strings.Replace(mssTiff, "MSS", "FUS", 1)
 	// pansharpenCmd := exec.Command("gdal_pansharpen.py", "-w 0.6 -w 0.1 -w 0.3 -w 0 -b 3 -b 2 -b 1", panTiff, mssTiff, fusedTIff, "-r cubic", "-of GTiff")
 	pansharpenCmd := exec.Command("gdal_pansharpen.py",
-		"-w", "0.6", "-w", "0.1", "-w", "0.3", "-w", "0",
+		"-w", "0.4", "-w", "0.2", "-w", "0.4", "-w", "0",
 		"-b", "3", "-b", "2", "-b", "1",
 		"-r", "cubic", "-of", "GTiff",
 		panTiff,

producer/image_proc_test.go

@@ -25,12 +25,12 @@ func TestGTiffToJPG(t *testing.T) {
 		t.Error(err)
 	}
 
-	tiff = "../data/051622/006/SJY01_PAN_20240516_101236_051622_096_006.tiff"
+	// tiff = "../data/051622/006/SJY01_PAN_20240516_101236_051622_096_006.tiff"
-	jpg = "../data/051622/006/SJY01_PAN_20240516_101236_051622_096_006.jpg"
+	// jpg = "../data/051622/006/SJY01_PAN_20240516_101236_051622_096_006.jpg"
-	err = GTiffToJPG(tiff, jpg, false)
+	// err = GTiffToJPG(tiff, jpg, false)
-	if err != nil {
+	// if err != nil {
-		t.Error(err)
+	// 	t.Error(err)
-	}
+	// }
 
 	tiff = "../data/051622/006/SJY01_FUS_20240516_101236_051622_096_006.tiff"
 	jpg = "../data/051622/006/SJY01_FUS_20240516_101236_051622_096_006.jpg"

producer/image_registration.go

@@ -21,7 +21,7 @@ const (
 	PanWidth                            = 9344 // 像素宽度
 	MssWidth                            = 2336
 	BlockNH                             = 4
-	BlockNW                             = 16
+	BlockNW                             = 8
 	OverlappedBlockLines                = 3000 // 重叠块的行数
 	DownSampled          ResampleMethod = "down_sample_pan"
 	UpSampled            ResampleMethod = "up_sample_mss"
@@ -221,7 +221,7 @@ func (r *Registrator) calcPhaseCorrelation(panImage gocv.Mat,
 
 				panBlock := panImage.Region(rect)
 				for band := 0; band < MssBands; band++ {
-					log.Debug("Band", band+1, ", processing block", bh, rect)
+					log.Debug("processing band:", band+1, ",block:", bh, rect)
 					mssBlock := mssImages[band].Region(rect)
 
 					// 处理每个分块
@@ -287,12 +287,17 @@ func (r *Registrator) calcDeltaCoeffs() error {
 
 			effectShift++
 			cx = append(cx, float64(shift.Block.coord.X+shift.Block.width/2)) // MSS 块在X方向没有分块
-			log.Debug("effective shift value:", effectShift, "cx:", shift.Block.coord.X, "dy:", shift.dy)
+			log.Debugf("effective shift value: %v, cx: %v, dx: %v, dy: %v",
+				effectShift, shift.Block.coord.X, shift.dx, shift.dy)
 			dx = append(dx, float64(shift.dx))
 			dy = append(dy, float64(shift.dy))
 
 		}
 
+		if len(cx) < 3 {
+			log.Errorf("No effective shift value found for band %d, skip delta coefficients calculation", i+1)
+			continue
+		} else {
 			var err error
 			if r.deltaXCoeffs[i], err = PolynomialFit(cx, dx, 1); err != nil {
 				log.Error("Error fitting deltaX coeffs: ", err)
@@ -304,8 +309,13 @@ func (r *Registrator) calcDeltaCoeffs() error {
 				return err
 			}
 		}
+	}
 
 	for i := 0; i < MssBands; i++ {
+		if len(r.deltaXCoeffs[i]) < 2 || len(r.deltaYCoeffs[i]) < 3 {
+			continue
+		}
+
 		log.Printf("Band %d:\n  delta_x = %.6f*x + %.6f, \n  delta_y = %.6f*x^2 + %.6f*x + %.6f\n",
 			i+1,
 			r.deltaXCoeffs[i][1], r.deltaXCoeffs[i][0],
@@ -317,6 +327,12 @@ func (r *Registrator) calcDeltaCoeffs() error {
 
 func (r *Registrator) DoCoRegestration() error {
 	for band := 0; band < MssBands; band++ {
+		if len(r.deltaXCoeffs[band]) < 2 || len(r.deltaYCoeffs[band]) < 3 {
+			log.Error("delta coefficients not calculated, skip co-registration")
+			r.registeredMssImages[band] = r.MssImages[band].Clone()
+			continue
+		}
+
 		mapX := gocv.NewMatWithSize(r.MssHeight, r.MssWidth, gocv.MatTypeCV32FC1)
 		mapY := gocv.NewMatWithSize(r.MssHeight, r.MssWidth, gocv.MatTypeCV32FC1)
 		for y := 0; y < r.MssHeight; y++ {

producer/jpg.go

@@ -19,7 +19,7 @@ func GTiffToJPG(ftiff, fjpg string, reversed bool) error {
 	defer ds.Close()
 
 	bands := ds.Bands()
-	log.Infof("TIFF file %s has %d bands", ftiff, len(bands))
+	log.Infof("creating JPG for TIFF file %s.", ftiff)
 
 	// 获取图像大小
 	width := bands[0].Structure().SizeX
@@ -62,9 +62,9 @@ func GTiffToJPG(ftiff, fjpg string, reversed bool) error {
 	img8 := gocv.NewMat()
 	defer img8.Close()
 	if bandsCnt == 1 {
-		img.ConvertToWithParams(&img8, gocv.MatTypeCV8UC1, 1.0/256.0*2, 0)
+		img.ConvertToWithParams(&img8, gocv.MatTypeCV8UC1, 1.0/256.0*1.5, 0)
 	} else {
-		img.ConvertToWithParams(&img8, gocv.MatTypeCV8UC3, 1.0/256.0*2, 0)
+		img.ConvertToWithParams(&img8, gocv.MatTypeCV8UC3, 1.0/256.0*1.5, 0)
 		if reversed {
 			channels := gocv.Split(img8)
 			gocv.Merge([]gocv.Mat{channels[2], channels[1], channels[0]}, &img8)

producer/output.go

@@ -13,7 +13,12 @@ import (
 )
 
 func (r *Registrator) SaveOriginalPanToGDALGTiff(tiffFile string) error {
-	return savePanToGDALGTiff(r.PanImage, tiffFile)
+	err := savePanToGDALGTiff(r.PanImage, tiffFile)
+	if err != nil {
+		return err
+	}
+
+	return nil
 }
 
 func (r *Registrator) SaveFilteredPanToGDALGTiff(tiffFile string) error {
@@ -68,10 +73,15 @@ func (r *Registrator) SaveRegisteredMssToGDALGTiff(tiffFile string) error {
 	r.rgbirImage = gocv.NewMat()
 	gocv.Merge(r.registeredMssImages[:], &r.rgbirImage)
 
-	return SaveBGRToGDALGTiff(r.rgbirImage,
+	err := SaveBGRToGDALGTiff(r.rgbirImage,
 		4, 5,
 		[]godal.ColorInterp{godal.CIBlue, godal.CIGreen, godal.CIRed, godal.CIUndefined},
 		tiffFile)
+	if err != nil {
+		return err
+	}
+
+	return nil
 }
 
 func (r *Registrator) SavePansharpenedToGDALGTiff(tiffFile string) error {

producer/phase_correlation.go

@@ -37,7 +37,7 @@ func (r *Registrator) calculateBlockPhaseShift(panBlock, mssBlock gocv.Mat) (goc
 
 	dx := shift.X
 	dy := shift.Y
-	log.Debugf("Block shift: dx = %f, dy = %f. response = %f \n", dx, dy, response)
+	log.Debugf("Block shift: dx = %f, dy = %f. response = %f", dx, dy, response)
 
 	return shift, response
 }

producer/scenes.go

@@ -100,6 +100,8 @@ func (r *Registrator) SaveScenesToTiff(panScenes []*Scene, mssScenes []*Scene) e
 			log.Errorf("save pan scene %d to tiff failed: %v", i+1, err)
 			return err
 		}
+
+		GTiffToJPG(scene.Tiff, strings.Replace(scene.Tiff, ".tiff", ".jpg", 1), false)
 	}
 
 	for i, scene := range mssScenes {
@@ -120,6 +122,8 @@ func (r *Registrator) SaveScenesToTiff(panScenes []*Scene, mssScenes []*Scene) e
 			log.Errorf("save mss scene %d to tiff failed: %v", i+1, err)
 			return err
 		}
+
+		GTiffToJPG(scene.Tiff, strings.Replace(scene.Tiff, ".tiff", ".jpg", 1), false)
 	}
 
 	return nil
@@ -133,7 +137,12 @@ func (r *Registrator) DoScenePansharpen(panScenes []*Scene, mssScenes []*Scene)
 		// filename := fmt.Sprintf("%s_%03d_FUS.tiff", path, i+1)
 
 		// SaveBGRToGDALGTiff(ihsImage, 3, 1.25, filename)
-		GDALPansharpen(panScenes[i].Tiff, mssScenes[i].Tiff)
+		fusedTiff, err := GDALPansharpen(panScenes[i].Tiff, mssScenes[i].Tiff)
+		if err != nil {
+			return err
+		}
+
+		GTiffToJPG(fusedTiff, strings.Replace(fusedTiff, ".tiff", ".jpg", 1), true)
 	}
 	return nil
 }

run.sh

@@ -5,3 +5,5 @@ go run cmd/*.go proc -p /Users/lan/workspace/sjy01/preprocessing/demo/output/051
 go run cmd/*.go proc -p /Users/lan/workspace/sjy01/preprocessing/demo/output/051922/SJY01_PAN_20240519_121433_051922_102.RAW -m /Users/lan/workspace/sjy01/preprocessing/demo/output/051922/SJY01_MSS_20240519_121433_051922_102.RAW -o data/051922 --fus
 
 go run cmd/*.go proc -p /Users/lan/workspace/sjy01/preprocessing/demo/output/051814/SJY01_PAN_20240517_111910_051814_098.RAW -m /Users/lan/workspace/sjy01/preprocessing/demo/output/051814/SJY01_MSS_20240517_111910_051814_098.RAW -o data/051814 --fus
+
+go run cmd/*.go proc -p /Users/lan/workspace/sjy01/preprocessing/demo/output/052723/SJY01_PAN_20240524_120820_052723_011.RAW -m /Users/lan/workspace/sjy01/preprocessing/demo/output/052723/SJY01_MSS_20240524_120820_052723_011.RAW -o data/052723 --fus