Skip to content
Projects
Groups
Snippets
Help
Loading...
Help
Submit feedback
Contribute to GitLab
Sign in / Register
Toggle navigation
GeoSentinel
Project
Project
Details
Activity
Releases
Cycle Analytics
Repository
Repository
Files
Commits
Branches
Tags
Contributors
Graph
Compare
Charts
Issues
0
Issues
0
List
Board
Labels
Milestones
Merge Requests
0
Merge Requests
0
CI / CD
CI / CD
Pipelines
Jobs
Schedules
Charts
Wiki
Wiki
Snippets
Snippets
Members
Members
Collapse sidebar
Close sidebar
Activity
Graph
Charts
Create a new issue
Jobs
Commits
Issue Boards
Open sidebar
Mario Chirinos Colunga
GeoSentinel
Commits
9c250220
Commit
9c250220
authored
May 24, 2018
by
Mario Chirinos Colunga
💬
Browse files
Options
Browse Files
Download
Email Patches
Plain Diff
new names
parent
e8a76029
Changes
3
Hide whitespace changes
Inline
Side-by-side
Showing
3 changed files
with
121 additions
and
0 deletions
+121
-0
normalizedIndex.py
geosentinel/normalizedIndex.py
+55
-0
productListToIndex.sh
tools/productListToIndex.sh
+66
-0
productToRGBJPG.sh
tools/productToRGBJPG.sh
+0
-0
No files found.
geosentinel/n
dvi
.py
→
geosentinel/n
ormalizedIndex
.py
View file @
9c250220
...
...
@@ -7,45 +7,40 @@ import numpy as np
import
gdal
from
osgeo
import
osr
def
geoint_pr_n
dvi
(
nir_fn
,
vir
_fn
,
outfile
):
def
geoint_pr_n
ormalizedIndex
(
b1_fn
,
b2
_fn
,
outfile
):
'''!
Normalized Difference
Vegetation Index (NDV
I)
@param
nir_fr String: Near Infra Red band (S2_B08
) Image file name
@param
vir String: Red Band (S2_B04
) Image file name
@param outfile String:
NDVI
image file name
Normalized Difference
Index (ND*
I)
@param
b1_fn String: band #1 (S2_B*
) Image file name
@param
b2_fn String: Band #2 (S2_B*
) Image file name
@param outfile String:
processed
image file name
'''
nir
=
gdal
.
Open
(
nir
_fn
)
red
=
gdal
.
Open
(
vir
_fn
)
b1
=
gdal
.
Open
(
b1
_fn
)
b2
=
gdal
.
Open
(
b2
_fn
)
rows
,
cols
,
geotransform
=
nir
.
RasterYSize
,
nir
.
RasterXSize
,
nir
.
GetGeoTransform
()
rows
,
cols
,
geotransform
=
b1
.
RasterYSize
,
b1
.
RasterXSize
,
b1
.
GetGeoTransform
()
#
# Read the input bands as numpy arrays.
np_
nir
=
nir
.
GetRasterBand
(
1
)
.
ReadAsArray
(
0
,
0
,
nir
.
RasterXSize
,
nir
.
RasterYSize
)
np_
red
=
red
.
GetRasterBand
(
1
)
.
ReadAsArray
(
0
,
0
,
red
.
RasterXSize
,
red
.
RasterYSize
)
# Read the input bands as numpy arrays.
np_
b1
=
b1
.
GetRasterBand
(
1
)
.
ReadAsArray
(
0
,
0
,
b1
.
RasterXSize
,
b1
.
RasterYSize
)
np_
b2
=
b2
.
GetRasterBand
(
1
)
.
ReadAsArray
(
0
,
0
,
b2
.
RasterXSize
,
b2
.
RasterYSize
)
# Convert the np arrays to 32-bit floating point to make sure division will occur properly.
np_
nir_as32
=
np_nir
.
astype
(
np
.
float32
)
np_
red_as32
=
np_red
.
astype
(
np
.
float32
)
np_
b1_as32
=
np_b1
.
astype
(
np
.
float32
)
np_
b2_as32
=
np_b2
.
astype
(
np
.
float32
)
# Calculate the NDVI formula.
np
.
seterr
(
divide
=
'ignore'
,
invalid
=
'ignore'
)
numerator
=
np
.
subtract
(
np_
nir_as32
,
np_red
_as32
)
denominator
=
np
.
add
(
np_
nir_as32
,
np_red
_as32
)
numerator
=
np
.
subtract
(
np_
b1_as32
,
np_b2
_as32
)
denominator
=
np
.
add
(
np_
b1_as32
,
np_b2
_as32
)
result
=
np
.
divide
(
numerator
,
denominator
)
result
=
np
.
multiply
((
result
+
1
),
(
2
**
7
-
1
))
driverTiff
=
gdal
.
GetDriverByName
(
'GTiff'
)
# driverJP2 = gdal.GetDriverByName('JP2OpenJPEG')
# driverJP2ECW = gdal.GetDriverByName('JP2ECW')
# raster = np.zeros((rows, cols), dtype=np.uint8)
output
=
driverTiff
.
Create
(
outfile
,
cols
,
rows
,
1
,
gdal
.
GDT_Byte
)
# output = driverTiff.CreateCopy(outfile+".tiff", nir, strict=0)
# output = driverJP2.CreateCopy(outfile, nir, strict=0)
output
.
GetRasterBand
(
1
)
.
SetNoDataValue
(
-
99
)
output
.
GetRasterBand
(
1
)
.
WriteArray
(
result
)
output
.
GetRasterBand
(
1
)
.
FlushCache
()
output
.
SetGeoTransform
(
geotransform
)
wkt
=
nir
.
GetProjection
()
wkt
=
b1
.
GetProjection
()
sr
=
osr
.
SpatialReference
()
sr
.
ImportFromWkt
(
wkt
)
...
...
tools/productListTo
NDVI
.sh
→
tools/productListTo
Index
.sh
View file @
9c250220
#!/bin/sh
BASEDIR
=
$1
OUTDIR
=
$2
WINDOW
=
$3
# Get Normalized index (B1-B2)/(B1+B2) for all products in PRODCUTSDIR
PRODCUTSDIR
=
$1
#Products directory
JP2DIR
=
$2
#JP2 Directory
WINDOW
=
$3
#Crop Window
OUTPUTDIR
=
$4
#Process output
$B1
=
$5
$B2
=
$6
count
=
1
ext
=
.jpg
BOX
=
$(
python3
"
$HOME
/git/GeoSentinel/geosentinel/polygonToBox.py"
"
$
3
"
)
BOX
=
$(
python3
"
$HOME
/git/GeoSentinel/geosentinel/polygonToBox.py"
"
$
WINDOW
"
)
echo
$BOX
cd
$
BASE
DIR
cd
$
PRODCUTS
DIR
for
f
in
$(
find
.
-type
f
-name
'*.zip'
)
do
filepattern
=
$(
echo
$f
|
cut
-d
"_"
-f6
)
"_"
$(
echo
$f
|
cut
-d
"_"
-f3
)
"_B"
filename
=
$
OUT
DIR$filepattern
filename
=
$
JP2
DIR$filepattern
echo
$filename
if
[
!
-e
$filename
"04
.jp2"
]
;
then
if
[
!
-e
$filename
$B1
"
.jp2"
]
;
then
unzip
-n
-j
$f
*
B04.jp2
-d
$OUT
DIR
unzip
" -n -j
$f
*"
$B2
".jp2 -d "
$JP2
DIR
fi
if
[
!
-e
$filename
"08
.jp2"
]
;
then
if
[
!
-e
$filename
$B2
"
.jp2"
]
;
then
unzip
-n
-j
$f
*
B08.jp2
-d
$OUT
DIR
unzip
" -n -j
$f
*"
$B1
".jp2 -d "
$JP2
DIR
fi
cd
$
OUT
DIR
NIR
=
$(
ls
$filepattern
"08
.jp2"
-t
|head
-1
)
VIR
=
$(
ls
$filepattern
"04
.jp2"
-t
|head
-1
)
cd
$
JP2
DIR
$B1FN
=
$(
ls
$filepattern$B1
"
.jp2"
-t
|head
-1
)
$B2FN
=
$(
ls
$filepattern$B1
"
.jp2"
-t
|head
-1
)
echo
"
NIR-"
$NIR
echo
"
VIR-"
$VIR
echo
"
B1-"
$B1FN
echo
"
B2-"
$B2FN
if
[
!
-d
"ndvi"
]
;
then
mkdir
ndvi
if
[
!
-d
$OUTPUTDIR
]
;
then
mkdir
$OUTPUTDIR
fi
fileout1
=
$(
echo
$filepattern
|
cut
-d
"_"
-f1
)
fileout2
=
$(
echo
$filepattern
|
cut
-d
"_"
-f2
)
if
[
!
-e
ndvi/
$fileout2
"_"
$fileout1
".tiff"
]
;
then
python3
$HOME
/git/GeoSentinel/geosentinel/n
dvi.py
$NIR
$VIR
ndvi/
$fileout2
"_"
$fileout1
".tiff"
if
[
!
-e
$OUTPUTDIR
"/"
$fileout2
"_"
$fileout1
".tiff"
]
;
then
python3
$HOME
/git/GeoSentinel/geosentinel/n
ormalizedIndex.py
$B1FN
$B2FN
$OUTPUTDIR
"/"
$fileout2
"_"
$fileout1
".tiff"
fi
if
[
!
-d
"video
_ndvi
"
]
;
then
mkdir
video
_ndvi
if
[
!
-d
"video"
]
;
then
mkdir
video
fi
# set GDAL_PAM_ENABLED=NO
gdal_translate
-projwin
$BOX
-projwin_srs
WGS84
-ot
Byte
-scale
0 255 0 255
-of
JPEG
ndvi/
$fileout2
"_"
$fileout1
".tiff"
video_ndvi
/
$fileout2
"_"
$fileout1
".jpg"
gdal_translate
-projwin
$BOX
-projwin_srs
WGS84
-ot
Byte
-scale
0 255 0 255
-of
JPEG
$OUTPUTDIR
"/"
$fileout2
"_"
$fileout1
".tiff"
video
/
$fileout2
"_"
$fileout1
".jpg"
# set GDAL_PAM_ENABLED=YES
rm
*
.xml
cd
..
cd
$PRODCUTSDIR
done
cd
$
OUT
DIR
"/video_ndvi"
cd
$
JP2
DIR
"/video_ndvi"
ffmpeg
-i
%
*
.jpg
-c
:v libx264
-vf
fps
=
10
-pix_fmt
yuvj422p ndvi.mp4
tools/
sentinel
ToRGBJPG.sh
→
tools/
product
ToRGBJPG.sh
View file @
9c250220
File moved
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment