Difficult 4 days of Big Problem of Noise Removal of Sentinel-1 imagery

First Question:
Sentinel-1 original HH or HV tiff image convention is still failing, including in ENVI  and SNAP software. Question on two results by using these methods convention contribute to the extents of these two images are different from the original.

Finally, I using gdal_translate solved the problem.But, I think this convention step is superfluous.

However, there is still GCPs information in the result, then direct adding of GCS_WGS_1984 projection information is invalid and will get a tiff file that pixel coordinate is over 5000+ degree northern latitude. This is obviously WRONG! So, I think Geo Points i.e. GCPs information leads to the error.

Therefore, there is second question.

Second Question:
How to add correct projection of a Geotiff image to a tiff without projection information?
Most relevant link is [1], however, there is no answer to the question. Then, I went on searching in Stack Overflow and Google!

Maybe link [3] can answer to the Q. Let me try!
I think my Question is very similar to link [3], 

Solving method is like below:

you can use gdal_edit.py for this:

gdal_edit.py -gcp 0.0 0.0 68.7535734331 81.6383694455 two.tiff

However, I don't how he extracted all of these points into one string. That's too bad!

GeoTiff official website! 3 formats composite GTiff files!

The following website[6] introduces how to georeference. As @Manimalis said, My Sentinel-1 image is similar to the website. It meas 'At this point, the image is not georeferenced yet.' Then I tested 'gdalwarp' method and only to find it is used to georeference. Therefore, it cannot be used. I just need to add these GCPs to Tiff file rather to georeference.

down vote

So, I succeed to georeference a tif file with GDAL with the use of 4 gcps (ground control points). To do this reprojection, I use gdal in command line.

First, use gdal translate like this :

gdal_translate -of GTiff -gcp 0 0 -6.848326  45.501053 -gcp 6862 0 -6.490975 45.501503 -gcp 0 1379 -6.762872 45.377363 -gcp 6862 1379 -6.545354 45.382523 "Inputimage.tif" "OutputImage.tif"

The ground control points are build like this : pixel coordinate in the image (x, y) and then geographical location (longitude, latitude).

I suggest to check with gdalinfo if the output file after gdaltranslate has been correctly fed with the gcps. With gdalinfo, you have to see in dos command something like that :

Coordinate System is `'
GCP Projection =
GCP[  0]: Id=1, Info=
          (0,0) -> (45.501053,-6.848326,0)
GCP[  1]: Id=2, Info=
          (6862,0) -> (45.501503,-6.490975,0)
GCP[  2]: Id=3, Info=
          (0,1379) -> (45.377363,-6.762872,0)
GCP[  3]: Id=4, Info=
          (6862,1379) -> (45.382523,-6.545354,0)

At this point, the image is not georeferenced yet. You have to use gdalwarp like this (there is several way to reproject, here I use "near" :

gdalwarp -r near -order 1 -co COMPRESS=NONE -dstalpha 
"OutputImageF_From_gdal_translate.tif" 
"FinalImage.tif"

shareeditflag                                                                                                     

answered Aug 24 '15 at 12:44

Manimalis

117113

Finally, I changed language(English-->Chinese) and website(Google-->Chinese) for searching. When I type in ' tiff 提取 gcps' only to find the function extracting GCPs in Geotiff files, that is 'GetGCPCount', 'GetGCPProjection', 'GetGCPs' in website[7](as following picture), all this function I also find using Google in English in link [8]. However, I think those functions can only use in C language and there is no practical examples in the website.
The following is the example in website[7].

I think I have found the answer for the question. Those functions of gdal is just answer.(Oct. 12, 2018)

However, I don't how to using SetGCPs......

Saturday(10.13), I attended practice test of civil servant in the morning, I still about 30 questions not to answer. I began to try to find the wrong reason for my code in the afternoon. In the evening, I found that change SetGCPs function to 'ds.SetGCPs(  gcp )' can get this error, 'TypeError: Dataset_SetGCPs() takes exactly 3 arguments (2 given)',  error of 'TypeError: Dataset_SetGCPs() takes exactly 3 arguments (4 given)' when I use 'ds.SetGCPs(  gcpcount , gcp, gcpproj )'. Therefore, I guess gcp should change format before set it into SetGCPs function.

Then, I found maybe the answer of website[3] is useful. I can try it! Besides, website[9] is GCP() function introduction. It can be following two format. I want to use first type.

GCP​(double x, double y, double pixel, double line)
GCP​(double x, double y, double z, double pixel, double line)

Or if you want to do it in a more pythonic way (adapted from gdal_edit.py):

from osgeo import gdal

gcp_items = filter(None, gcp_string.split("-gcp"))

gcp_list = []
for item in gcp_items:
    pixel, line, x, y = map(float, item.split())
    z = 0
    gcp = gdal.GCP(x, y, z, pixel, line)
    gcp_list.append(gcp)

ds = gdal.Open("two.tiff", gdal.GA_Update)
wkt = ds.GetProjection()
ds.SetGCPs(gcp_list, wkt)
ds = None

Note: this is assuming that the projection of one.tiff and two.tiff are the same. Otherwise you will have to get the wkt from one.tiff.

shareedit

edited Aug 12 '16 at 14:39

answered Aug 12 '16 at 14:03

John Truckenbrodt

14616

However, when I using this method and modified my code, after runing well, I got same error as before. gdal.GCP() function is used for converting the format as same to function of GetGCPs(). I do again some idle work. Sad! 

However, this is no error when I just put gcp_list and gcpproj as above picture, but noise tif cannot get GCPs and projection. WHY?????????!!!!!!!!!!!!

I finally solved the error !!!!!! Missing the word 'ds = None', the answer is on website[3]. 
The following picture is result of noise tiff adding GCPs and projections from Sentinel-1 images. My God! Relevant question I posted on Stack Overflow is website[10]

4 days have run off! 

Until now(10.14), Second question has successfully solved by myself!

However, I found I have do the Thermal Noise Removal in re-processing Python program(10.15). Therefore, I don't know Tingting Zhu(WHU) in her PhD graduation thesis to do this work by herself. The following picture is SNAP TNR step introduction

My pre-processing p Python program( website[11]) includes 4 steps.

Eventually, I have to again resume my GLCM program modification work.

Last week, my largest harvest is improve my ability of writing and convention from point data  to raster interpolation based on XML. This experience is very inspiring! ( 17:09 in Oct. 15, 2018)

Reference Link:
[1] https://stackoverflow.com/questions/48770002/project-raster-file-using-gcps
[2] Gdal_translate Webpage: https://www.gdal.org/gdal_translate.html
[3] Add GCP to GeoTIFF https://gis.stackexchange.com/questions/206544/in-python-gdal-how-to-add-ground-control-point-to-geotiff-image-by-using-gdal-tr 
[4]  Use GDAL command line to copy projections  https://gis.stackexchange.com/questions/108673/use-gdal-command-line-to-copy-projections
[5] https://gis.stackexchange.com/questions/257904/converting-geotiff-with-gcps-into-projcs-with-python-gdal
[6] https://gis.stackexchange.com/questions/159004/creating-geotiff-from-tiff-using-gdal/159436#159436
[7] https://blog.csdn.net/bcbobo21cn/article/details/51014246?locationNum=12&fps=1
[8] https://www.gdal.org/classGDALDataset.html
[9] https://gdal.org/java/org/gdal/gdal/GCP.html
[10] https://gis.stackexchange.com/questions/298806/error-of-copying-gcps-information-from-geotiff-to-another-tiff
[11] https://www.dropbox.com/s/bk25khnha7w85wg/procSentinelRTC_recipe.py?dl=0

Import GDAL and arcpy ERROR in 64-bit Spyder

Question:

Environment: 64-bit andconda(Spyder) and IPython console Information is：

Two modules cannot are imported to Spyder under this environment.

1. import GDAL Error:

Solving Reference Link:
 https://stackoverflow.com/questions/52565731/how-to-solve-import-gdal-error-in-64-bit-anaconda 

Solving Steps:

To type the following command sentence in Anaconda Prompt.

conda install -c conda-forge georasters

2. import arcpy Error:

Two Solving Reference Link:

① https://gis.stackexchange.com/questions/81066/installing-64-bit-geoprocessing-of-arcgis-for-desktop

②https://gis.stackexchange.com/questions/85408/debugging-importerror-no-module-named-arcpy

Solving Steps:

Firstly, as @crafty762 said in hyperlink ①, I have to install ArcGIS for Desktop Background Geoprocessing(64-bit).

Then look for hyperlink②, there are many methods to solve the problem. In the end, I solved the problem using @PolyGeo 's answer.

In the bottom of  Importing ArcPy is:

In my PC, I investigated to C:\Python27\ArcGISx6410.2\Lib\site-packages, which is the installing directory of ArcGIS for Desktop Background Geoprocessing(64-bit), there is a similar file named DTBGGP64.pth then opened it in Notepad++ :

Until then, I just realized the error is wrong first directory importing into PYTHONPATH manager, it should be the bin64 folder rather than bin folder:

Finally, I added these three directories to PYTHONPATH manager in Spyder, there is no error:

Happy Night! 

Oct. 2,2018 at midnight

NU 微信公众号挂靠在学校认证下流程

        Wechat 对个人公众号需管理员与主体一致，管理员之更改益难，遂去年偷懒，将此trouble遗留至今，无奈W师姐迁恼于Z，Z求助于我，念吾之过，盖将此任耽于己。求助网络，wechat以及其他，得二法，其一，山大白杨社之奇策，其二乃“账号迁移”，然该法费事费钱费力。于Z商讨，因法二为根本，唯此方可根本解此难题。Z欲吾问他人，余先网络搜索，叹终归无力，乃求D学部之人F，F实友善人，下图为其关键信息。遂废弃法二，其实为给定法三。

        翌日，问认证之来龙去脉，得其部负责之C老师为重要人员，亦可询我院之负责老师，吾问W老师，建议W老师往C老师方向刨根问底，W老师亦推荐问我院W同学，因W同学与C老师为熟，后证为假，这是后话。W老师人脉甚广，一刻钟内即回复响应策略，如下图，终得解，已一周飞逝。

附录：

介绍信

    为了凝聚校内外师生、传播学科知识、加强思想政治宣传教育，XX学院院特申请认证微信公众号，拟认证的微信公众号信息如下：

公众号类型：订阅号

公众号名称：

公众号ID：

公众号管理员：

工作证号:

身份证号:

联系方式：

邮箱：

特此申请，为盼！

                          单位公章

                         单位负责人签字：

                       年   月   日

最后，这种方法还是废掉了，因为微信公众号的名称比较正式与官方，学院不允许使用如此正式的名字发布研会的消息。

最终以院行政老师自己申请个人账号作为公众号，以其可能常驻此工作为媒介，解决学生因为毕业而流动的问题。只不过，当其也离职时，必然再次面临这个问题，不过短期内的比较容易地解决办法也只有如此了！

Python Read XML Files

This code is used to extract noise value and rows and columns of sample points from auxiliary XML files of Sentinel-1 image source folders.The file the program reads is attached to the last.

Source Code:

Auxiliary XML Files Links：
https://mega.nz/#!dOIhyajZ!ymdjwhzSgoV1QOubBZpyuXo6WtpMP5t9yG3ydElFZdA

About Sentinel-1 Subswath Noise Removal Problem using FFT

I don't why I cannot to find the correct solution to achieve my goal, removing sub-swath noise in Sentinel-1 HV polarization imagery. I’ve tried as my best in StackOverflow website and CSDN in Chinese. But, those efforts are meaningless. I still don't know how to using FFT in python to get my succeed.

I've referred to these website(https://stackoverflow.com/questions/38476359/fft-on-image-with-python , https://stackoverflow.com/questions/32951036/opencv-python-mask-on-fftimage-why-do-we-need-two-channels ) and some PPT on introducing Fourier Transformation. Besides, I get five steps to make FFT, as following:

1.read image

2.get fft of image --> f

3.crate mask

4.multiply f with mask --> g

5.get inverse of g

My puzzle encountered in coding is how to put FFT apply to tiff format images, websites in above is for processing .jpg or .png format. Another question is when I found a sample demo about processing tif( Coding is given in the end of the article), I unexpectedly to find the result after FFT existed a fatal bug, data value of tiff image changed. Until now, I've wasted a week days in the begining of  new semester of my master degree in the last year. Oh, my god!

I surrender and give up for this way. Maybe I won't remove boundary noise in S-1 HV images.

Coding:

# -*- coding: utf-8 -*-

"""

Created on Tue Sep 04 14:12:09 2018

@author: Baikal

"""

import cv2

import scipy

import numpy as np

tif = cv2.imread( r'C:\Users\Baikal\Desktop\delete\Aim_Beta0_HV.tif', 2 )

rows, cols = tif.shape

## Get FFT imgage.

fft = np.fft.fft2( tif )

## Move to center

fftshift = np.fft.fftshift( fft )

## Create Mask tif file.

rowl = rows / 2 - 100

rowr = rows / 2 + 100

columnl = cols / 2 - 100

columnr = cols / 2 + 100

mask = np.ones( tif.shape, np.float32 )

for row in range( rowl, rowr ):

    for col in range( columnl, columnr ):

        mask[ row, col ] = 0

        

##  Filter.    

fftshift_mask = fftshift * mask 

## Get inverse of fft.

ifft = np.fft.ifft2( np.fft.ifftshift( fftshift_mask ) )

# To show normally , to get abs of complex number(fft).

ifft_img = np.abs( ifft ) 

# Save FFT image to disk.

scipy.misc.imsave( 'ifft.tif', ifft_img )

• Original TIFF image(Data range below 0)

• FFT result(Data range from 0 to 255)

LibSVM 参数说明 && 参数寻优

参考链接：
https://blog.csdn.net/lanchunhui/article/details/72729815
https://blog.csdn.net/garfield2005/article/details/9111995



参数寻优参考链接：
http://xiaosheng.me/2016/06/09/article17/