Download this notebook from github.
Example to compute image cross spectrum between different looks within a WV Sentinel-1 SLC product
in WV images (20 km X 20 km) there is only one burst (thus no overlapping area).
[1]:
import xsar
import warnings
import logging
from importlib import reload
reload(logging)
logging.basicConfig(level=logging.INFO)
logging.info('go')
# warnings.simplefilter(action='ignore', category=FutureWarning)
warnings.filterwarnings("ignore")
#tmp = xsar.get_test_file('S1A_WV_SLC__1SSV_20160510T101603_20160510T102143_011195_010EA1_Z010.SAFE')
tmp = xsar.get_test_file('S1A_WV_SLC__1SSV_20191005T100804_20191005T101314_029322_035533_0597.SAFE') # very small product
# print(tmp)
# tmp = '/home/datawork-cersat-public/cache/project/mpc-sentinel1/data/esa/sentinel-1a/L2/WV/S1A_WV_OCN__2S/2016/131/S1A_WV_OCN__2SSV_20160510T101603_20160510T102143_011195_010EA1_F972.SAFE'
# tmp = '/home/datawork-cersat-public/cache/project/mpc-sentinel1/data/esa/sentinel-1a/L1/WV/S1A_WV_SLC__1S/2016/131/S1A_WV_SLC__1SSV_20160510T101603_20160510T102143_011195_010EA1_7038.SAFE'
# tmp = '/home/datawork-cersat-public/cache/project/mpc-sentinel1/data/esa/sentinel-1a/L1/WV/S1A_WV_SLC__1S/2023/014/S1A_WV_SLC__1SSV_20230114T130811_20230114T133351_046780_059BC5_FED5.SAFE'
strwv = 'SENTINEL1_DS:'+tmp+':WV_017'
print(strwv)
xsar_instance = xsar.Sentinel1Dataset(strwv)
xsar_instance
INFO:root:go
SENTINEL1_DS:/home1/scratch/agrouaze/xsardatasync/xsardata/S1A_WV_SLC__1SSV_20191005T100804_20191005T101314_029322_035533_0597.SAFE:WV_017
[1]:
<Sentinel1Dataset full coverage object>
[2]:
xsar_instance.dataset
[2]:
<xarray.Dataset>
Dimensions: (line: 4923, sample: 5681, pol: 1)
Coordinates:
* line (line) int64 0 1 2 3 4 5 6 ... 4917 4918 4919 4920 4921 4922
* sample (sample) int64 0 1 2 3 4 5 ... 5675 5676 5677 5678 5679 5680
* pol (pol) object 'VV'
Data variables:
digital_number (pol, line, sample) complex64 dask.array<chunksize=(1, 4923, 5000), meta=np.ndarray>
time (line) datetime64[ns] 2019-10-05T10:11:58.450274048 ... 2...
sampleSpacing float64 1.498
lineSpacing float64 4.133
Attributes: (12/15)
name: SENTINEL1_DS:/home1/scratch/agrouaze/xsardatasync/xsar...
short_name: SENTINEL1_DS:S1A_WV_SLC__1SSV_20191005T100804_20191005...
product: SLC
safe: S1A_WV_SLC__1SSV_20191005T100804_20191005T101314_02932...
swath: WV
multidataset: False
... ...
start_date: 2019-10-05 10:11:58.450260
stop_date: 2019-10-05 10:12:01.432986
footprint: POLYGON ((-61.07570128124852 28.34141558287481, -61.28...
coverage: 20km * 20km (line * sample )
orbit_pass: Descending
platform_heading: -167.4153352100693[3]:
xsar_instance.dataset.digital_number.values[0:10,0:10]
[3]:
array([[[ 0. +0.j, 0. +0.j, 0. +0.j, ..., 0. +0.j,
0. +0.j, 0. +0.j],
[ 0. +0.j, 0. +0.j, 0. +0.j, ..., 0. +0.j,
0. +0.j, 0. +0.j],
[ 0. +0.j, 0. +0.j, 0. +0.j, ..., 14. +15.j,
64. -30.j, 43. -46.j],
...,
[ -12.+243.j, -69. -81.j, 13. -53.j, ..., 82. +53.j,
190. +4.j, 108.-126.j],
[ 87. -29.j, -120. -67.j, 83. +82.j, ..., 123.+140.j,
232. -8.j, -123.-134.j],
[ -29. +21.j, -34. -6.j, -44. +99.j, ..., 213. +40.j,
325.+204.j, 99. -20.j]]], dtype=complex64)
display digital numbers
[4]:
import holoviews as hv
import numpy as np
from holoviews.operation.datashader import datashade,rasterize
hv.extension('bokeh')
[5]:
tmp = abs(xsar_instance.dataset['digital_number'].sel(pol='VV'))
DN = tmp.rolling({"sample":10,"line":10}).mean()
DN
[5]:
<xarray.DataArray 'digital_number' (line: 4923, sample: 5681)>
dask.array<truediv, shape=(4923, 5681), dtype=float32, chunksize=(4923, 5009), chunktype=numpy.ndarray>
Coordinates:
* line (line) int64 0 1 2 3 4 5 6 7 ... 4916 4917 4918 4919 4920 4921 4922
* sample (sample) int64 0 1 2 3 4 5 6 ... 5674 5675 5676 5677 5678 5679 5680
pol <U2 'VV'
Attributes:
comment: denoised digital number, read at full resolution
history: digital_number: measurement/s1a-wv1-slc-vv-20191005t101158-2019...[6]:
rasterize(hv.Image(DN.values).opts(width=800,height=700,cmap='Greys',colorbar=True))# clim=(0,40)
[6]:
compute the periodograms and the looks inside the periodograms to get a cross spectra
[7]:
from xsarslc.processing.xspectra import compute_WV_intraburst_xspectra
import time
t0 = time.time()
xs0 = compute_WV_intraburst_xspectra(dt=xsar_instance.datatree,
polarization='VV',periodo_width={"line":2000,"sample":2000},
periodo_overlap={"line":1000,"sample":1000})
print('time to compute x-spectra on WV :%1.1f sec'%(time.time()-t0))
xs0
time to compute x-spectra on WV :236.2 sec
[7]:
<xarray.Dataset>
Dimensions: (bt_thresh: 5, 0tau: 3, freq_line: 97,
freq_sample: 274, 1tau: 2, 2tau: 1, c_sample: 2,
c_line: 2, k_gp: 4, phi_hf: 5,
lambda_range_max_macs: 10)
Coordinates:
pol <U2 'VV'
* bt_thresh (bt_thresh) int64 5 50 100 150 200
k_rg (freq_sample) float64 0.0 0.003143 ... 0.858
k_az (freq_line) float64 -0.1505 -0.1473 ... 0.1505
line int64 2461
sample int64 2840
longitude float64 -61.2
latitude float64 28.27
* k_gp (k_gp) int64 1 2 3 4
* phi_hf (phi_hf) int64 1 2 3 4 5
* lambda_range_max_macs (lambda_range_max_macs) float64 50.0 ... 275.0
Dimensions without coordinates: 0tau, freq_line, freq_sample, 1tau, 2tau,
c_sample, c_line
Data variables: (12/26)
incidence float64 24.19
ground_heading float64 -168.5
sensing_time datetime64[ns] 2019-10-05T10:11:59.941637632
sigma0 float64 0.1942
nesz float64 0.001677
bright_targets_histogram (bt_thresh) int64 54 0 0 0 0
... ...
corner_line (c_line) int64 25 4897
corner_sample (c_sample) int64 25 5655
burst_corner_longitude (c_sample, c_line) float64 -61.08 ... -61.32
burst_corner_latitude (c_sample, c_line) float64 28.34 28.16 28.38 28.2
cwave_params (k_gp, phi_hf, 2tau) float64 4.855 ... -0.334
macs (lambda_range_max_macs, 2tau) complex128 (0.167...
Attributes: (12/21)
name: SENTINEL1_DS:/home1/scratch/agrouaze/xsardatasync...
short_name: SENTINEL1_DS:S1A_WV_SLC__1SSV_20191005T100804_201...
product: SLC
safe: S1A_WV_SLC__1SSV_20191005T100804_20191005T101314_...
swath: WV
multidataset: False
... ...
radar_frequency: 5405000704.0
azimuth_time_interval: 0.0006059988896297211
tile_width_line: 20138.934607
tile_width_sample: <xarray.DataArray 'cumulative ground length' ()>\...
tile_overlap_sample: 0.0
tile_overlap_line: 0.0[8]:
import numpy as np
from xsarslc.processing import xspectra
xs = xs0.swap_dims({'freq_line': 'k_az', 'freq_sample': 'k_rg'})
xs = xspectra.symmetrize_xspectrum(xs['xspectra_2tau'], dim_range='k_rg', dim_azimuth='k_az')
############################################ real part ############################
coS=np.abs(xs.mean(dim='2tau').real)
coS_reduced = coS.where(np.logical_and(np.abs(coS.k_rg)<=0.14, np.abs(coS.k_az)<=0.14), drop=True)
from matplotlib import pyplot as plt
%matplotlib inline
from matplotlib import colors as mcolors
cmap = mcolors.LinearSegmentedColormap.from_list("", ["white","violet","mediumpurple","cyan","springgreen","yellow","red"])
PuOr = plt.get_cmap('PuOr')
plt.figure(figsize=(8,6))
#coS_reduced_av=coS_reduced.rolling(k_rg=3, center=True).mean().rolling(k_az=3, center=True).mean()
coS_reduced.plot(cmap=cmap, levels=20, vmin=0)
plt.grid()
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.axis('scaled')
############################################ imag.part ############################
ims = xs.mean(dim='2tau').imag.squeeze()
xS_red = ims.where(np.logical_and(np.abs(ims.k_rg)<=0.14, np.abs(ims.k_az)<=0.14), drop=True)
#xS_av=xS_red.rolling(k_rg=3, center=True).mean().rolling(k_az=3, center=True).mean()
plt.figure(figsize=(8,6))
xS_red.plot(x='k_rg',cmap=PuOr)
plt.grid()
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.axis('scaled')
[8]:
(-0.13984955506418678,
0.13984955506418678,
-0.13949481344309983,
0.13949481344309983)
