diff --git a/cmd/proc.go b/cmd/proc.go index ef681a1..34b0e4d 100644 --- a/cmd/proc.go +++ b/cmd/proc.go @@ -13,10 +13,11 @@ import ( ) var ( - params producer.Params - saveStrip bool - doRRC bool - lutDir string + params producer.Params + saveStrip bool + doRRC bool + doMomentMatching bool + lutDir string ) var procCmd = &cobra.Command{ @@ -50,6 +51,10 @@ var procCmd = &cobra.Command{ reg.DoRRC(lutDir) } + if doMomentMatching { + reg.DoMomentMatching() + } + if err := reg.DoPhaseCorrelation(); err != nil { logrus.Fatal(err) } @@ -99,6 +104,7 @@ func init() { procCmd.Flags().StringVarP(¶msXML, "params", "x", "", "params xml file path") procCmd.Flags().StringVarP(&lutDir, "lut", "l", "data/lut", "LUT directory") procCmd.Flags().BoolVarP(&doRRC, "rrc", "", false, "do RRC") + procCmd.Flags().BoolVarP(&doMomentMatching, "mm", "", false, "do moment matching") } func initParams() producer.Params { diff --git a/cmd/rrc.go b/cmd/rrc.go index bf4135a..48f1eeb 100644 --- a/cmd/rrc.go +++ b/cmd/rrc.go @@ -6,8 +6,9 @@ import ( ) var ( - panDS string - mssDS string + panDS string + mssDS string + lutOut string ) var rrcCmd = &cobra.Command{ @@ -15,7 +16,7 @@ var rrcCmd = &cobra.Command{ Short: "Run RRC algorithm on an image", Long: `Run RRC algorithm on an image`, Run: func(cmd *cobra.Command, args []string) { - rrc := rrc.NewRRC() + rrc := rrc.NewRRC(lutOut) rrc.StatisticalPAN(panDS) rrc.StatisticalMSS(mssDS) }, @@ -24,6 +25,7 @@ var rrcCmd = &cobra.Command{ func init() { rrcCmd.Flags().StringVarP(&panDS, "pan-dataset", "p", "data/RAW/pan.txt", "path to pan dataset") rrcCmd.Flags().StringVarP(&mssDS, "mss-dataset", "m", "data/RAW/mss-none.txt", "path to mss dataset") + rrcCmd.Flags().StringVarP(&lutOut, "lut-out-dir", "o", "data/rrc", "path to output lut file") rootCmd.AddCommand(rrcCmd) } diff --git a/pkg/producer/rrc.go b/pkg/producer/rrc.go index f128b1d..36408b4 100644 --- a/pkg/producer/rrc.go +++ b/pkg/producer/rrc.go @@ -40,3 +40,12 @@ func (r *Registrator) DoRRC(lutDir string) error { return nil } + +func (r *Registrator) DoMomentMatching() error { + rrc.DoMomentMatching(r.PanImage) + for i := 0; i < 4; i++ { + rrc.DoMomentMatching(r.MssImages[i]) + } + + return nil +} diff --git a/pkg/rrc/moment_matching.go b/pkg/rrc/moment_matching.go new file mode 100644 index 0000000..1955173 --- /dev/null +++ b/pkg/rrc/moment_matching.go @@ -0,0 +1,56 @@ +package rrc + +import ( + "math" + + log "github.com/sirupsen/logrus" + "gocv.io/x/gocv" +) + +// Destriping multisensor imagery with moment matching + +func DoMomentMatching(originalImg gocv.Mat) { + probes := originalImg.Cols() + log.Printf("do moment matching for %d probes", probes) + // 第i个探元的像元均值和标准差 + means := make([]float64, probes) + stds := make([]float64, probes) + // 计算每个探元的均值和标准差 + for x := 0; x < originalImg.Cols(); x++ { + var total int64 + n := originalImg.Rows() + var dn uint16 + for y := 0; y < n; y++ { + dn = uint16(originalImg.GetShortAt(y, x)) + total += int64(dn) + } + means[x] = float64(total) / float64(n) + + var a float64 + for y := 0; y < n; y++ { + dn = uint16(originalImg.GetShortAt(y, x)) + a += math.Pow(float64(dn)-means[x], 2) + } + stds[x] = math.Sqrt(a / float64(n)) + } + + // 列参考值和列参考标准差 + var mu float64 + var sig float64 + for x := 0; x < probes; x++ { + mu += means[x] + sig += stds[x] + } + mu = mu / float64(probes) + sig = sig / float64(probes) + + log.Printf("mean reference value: %f, std reference value: %f", mu, sig) + // 修正 DN_adjusted[i] = (DN[i] - means[i]) *sig/stds[i]+mu + for x := 0; x < originalImg.Cols(); x++ { + for y := 0; y < originalImg.Rows(); y++ { + dn := uint16(originalImg.GetShortAt(y, x)) + dn_adjusted := uint16((float64(dn)-means[x])*sig/stds[x] + mu) + originalImg.SetShortAt(y, x, int16(dn_adjusted)) + } + } +} diff --git a/pkg/rrc/rrc.go b/pkg/rrc/rrc.go index 3ea3865..14b4d4d 100644 --- a/pkg/rrc/rrc.go +++ b/pkg/rrc/rrc.go @@ -4,6 +4,7 @@ import ( "bufio" "fmt" "os" + "path/filepath" "runtime" "sync" @@ -21,17 +22,18 @@ const ( ) type RRC struct { + LUTOutDir string Histograms [5]ProbeHistogram } -func NewRRC() *RRC { - r := RRC{} +func NewRRC(dir string) *RRC { + r := RRC{LUTOutDir: dir} r.Histograms[0].init(PANCameraProbeNum) for i := 1; i < 5; i++ { r.Histograms[i].init(MSSCameraProbeNum) } - os.MkdirAll("data/rrc", 0755) + os.MkdirAll(r.LUTOutDir, 0755) return &r } @@ -110,7 +112,7 @@ func (rrc *RRC) StatisticalPAN(dsfile string) error { log.Println("compute PAN histogram...") rrc.Histograms[0].compute() log.Println("save PAN gray table matrix.") - rrc.Histograms[0].saveLUT("data/rrc/B0.LUT") + rrc.Histograms[0].saveLUT(filepath.Join(rrc.LUTOutDir, "B0.LUT")) return nil } @@ -189,7 +191,7 @@ func (rrc *RRC) StatisticalMSS(dsfile string) error { log.Println("compute MSS histogram...") rrc.Histograms[i].compute() log.Println("save MSS gray table matrix.") - rrc.Histograms[i].saveLUT(fmt.Sprintf("data/rrc/B%d.LUT", i)) + rrc.Histograms[i].saveLUT(filepath.Join(rrc.LUTOutDir, fmt.Sprintf("B%d.LUT", i))) } return nil