Thursday, April 19, 2018

Pemetaan Tumpahan Minyak di Pantai Balikpapan Menggunakan Citra Radar Sentinel 1A dan Google Earth Engine

Kasus tumpahan minyak di Pantai Balikpapan awal April 2018 menjadi isu hangat bulan ini. Berbagai media memberitakan dampak negatif dari bencana lingkungan ini. Sebagai geoscientist, isu ini menarik buat saya pribadi karena sebaran tumpahan minyak di laut saat ini sudah bisa dideteksi dengan baik dan dipetakan sebarannya menggunakan data radar SAR.

Dengan bermodal layanan cloud image processing dari Google Earth Engine (GEE) saya mencoba bereksperimen memetakan sebaran area tumpahan minyak di Pantai Balikpapan. GEE menyediakan akses ke repositori data Sentinel 1 SAR yang saat ini operasional, dan kebetulan selang beberapa hari dari kejadian tumpahan minyak, satelit Sentinel 1A lewat di atas Balikpapan.

Saya menggunakan kode awal dari Simon Ilyuschenko dari GEE Team yang digunakan untuk pemetaan area banjir dari citra Sentinel, yang saya modifikasi untuk lebih peka terhadap tumpahan minyak. Dengan menggunakan teknik pemrosesan image tresholding dan image difference dari citra yang merekam pada kondisi sebelum ( minggu ke 3 maret 2018) dan sesudah kejadian tumpahan minyak (minggu pertama April 2018), area sebaran tumpahan minyak dapat diidentifikasi. Tumpahan minyak dapat diidentifikasi dengan mudah dari citra SAR karena  minyak dan air mempunyai konstanta dielektrik yang berbeda, dan spektrum gelombang mikro pada panjang gelombang C Band sangat sensitif dengan perbedaan konstanta dielektrik, terutama pada polarisasi VV (Vertical send Vertical return)

Meskipun demikian, kode ini belum sempurna karena masih terdapat area - area yang misklasifikasi, sehingga revisi dan pemrosesan citra lebih lanjut perlu dilakukan. Sebenarnya identifikasi area tumpahan minyak ini lebih gampang dilakukan kalau kita menggunakan software desktop yang bisa mengolah data radar seperti ESA SNAP, ENVI Sarscape atau ERDAS dan PCI Geomatica, tetapi GEE menawarkan sesuatu yang berbeda, yaitu kepraktisan, karena kita tidak perlu mendownload data yang ukurannya cukup besar, dan menyiapkan hardware komputer yang mumpuni. Semua bisa dikerjakan di Cloud hanya dengan mengetikkan beberapa baris kode javascript API. Berikut ini screenshot hasilnya.


Rekan - rekan yang tertarik mencoba bisa menggunakan kode yang saya gunakan di bawah ini. Sukur -sukur bisa memperbaiki dan mengembangkan lebih lanjut. Jangan lupa colek saya ya kalau ada revisi.

var geometry = geometry;

// Load Sentinel-1 images to map Oil Spill in Balikpapan Indonesia April 2018.
// This script was originally written by Simon Ilyushchenko (GEE team)
// and adapted by Simon Gascoin (CNRS/CESBIO) and Michel Le Page (IRD/CESBIO)
// and also adapted by Bramantiyo Marjuki (Indonesia Ministry of Public Works/Regional Development Departement of Diponegoro University)
// Default location
var pt = geometry; //Coulommiers

// Load Sentinel-1 C-band SAR Ground Range collection (log scaling, VV co-polar)
var collection =  ee.ImageCollection('COPERNICUS/S1_GRD').filterBounds(pt)
.filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VV'))
.select('VV');

// Filter by date 
var before = collection.filterDate('2018-03-10', '2018-03-29').mosaic();
var after = collection.filterDate('2018-03-31', '2018-04-06').mosaic();
var diff = after.subtract(before);
// Threshold smoothed radar intensities to identify "spill" areas.
var SMOOTHING_RADIUS = 200;
var DIFF_UPPER_THRESHOLD = -0.9;

var diff_smoothed = after.focal_median(SMOOTHING_RADIUS, 'circle', 'meters').subtract(before.focal_median(SMOOTHING_RADIUS, 'circle', 'meters'));

var diff_thresholded = diff_smoothed.lt(DIFF_UPPER_THRESHOLD);

// Display map
Map.centerObject(pt, 13);
Map.addLayer(before, {min:-30,max:0}, 'Before oil spill');
Map.addLayer(after, {min:-30,max:0}, 'After oil spill');
Map.addLayer(before.addBands(after).addBands(diff), {min:-10,max:10},'composite', 0);
Map.addLayer(after.subtract(before), {min:-10,max:10}, 'After - before', 0);
Map.addLayer(diff_smoothed, {min:-10,max:10}, 'diff smoothed', 0);
Map.addLayer(diff_thresholded.updateMask(diff_thresholded),
  {palette:"0000FF"},'oil spill areas - blue',1);

Saturday, February 17, 2018

Using Mapproxy to cache and host both online and offline basemaps (Basemap series Part 5)

Here we go again

In this post now I will continue about how to simultaneous host and cache basemaps, either from online or offline sources.

First, If you want to just use online basemaps to be available or hosted on your favorite protocols (like serving google Street Map to WMS or ArcGIS REST service) you could use Portable Basemap Server. I have previously covered it, check this LINK.

If you want to only cache the tiles to be used for offline needs, you can use MOBAC or SasPlanet.

Now if you want to have capabilities like both serving basemaps directly to your user using client server services protocols (either online or local online or even pure offline situation) while maintain online access, plus adding caching capabilities (useful if you want your application to be as fast as possible), you can use MAPPROXY.

It can do what PBS can do plus it has caching ability, so you wont have to worry if your internet networks get troubled, because the generated caches will serve as backup renderer. Basically it is combine PBS, and offline cacher like MOBAC or Sasplanet.

Now about the procedures:

1. Install mapproxy using python PIP from their official page.
2. You are python blind, dont worry this SITE provide you offline installer (install it as windows services for your convenient). Because the program will have their own webserver, beware of Port Conflict, make sure you choose unused HTTP PORT on your server/computer. 
3. Now deep down to settings, make YAML configuration file. The structure of the code will looks like this
services:
  demo:
  tms:
  wms:
    md: 
      title: Google MapProxy WMS Server
  wmts:
    kvp: false
    restful: true

sources:
  google_image:
    type: tile
    url: http://mt1.google.com/vt/lyrs=s&x=%(x)s&y=%(y)s&z=%(z)s&s=Galileo
    grid: grid_webmerc

grids:
  grid_webmerc:
    base: GLOBAL_WEBMERCATOR
  grid_merc:
    base: GLOBAL_MERCATOR

caches:
  cache_google_image:
    grids: [grid_merc]
    sources: [google_image]
    cache:
      type: mbtiles

layers:
  - name: layer_google_image
    title: Google Satellite Map
    sources: [cache_google_image]

globals:
  image:
    resampling_method: bilinear
    jpeg_quality: 90

4. The above YAML config file is dedicated to serve google image to your clients. Refer to this LINK for documentation. The most important things is, you must know the exact tile query URL of your target basemap providers.

5. Now save the YAML to directory C:\ProgramData\Mapproxy

6. Next, open the mapproxy_app.py or mapproxy_srv.py in C:\ProgramFiles(x86)\Mapproxy, and change the app_config line to your created YAML file.



7. Load up the mapproxy by launch it from CMD from the program files directory, like this (dont forget to run as administrator to avoid denied access error).


8. or from Windows Services (when it installed as services). Load the windows services from Control Panel or windows search tool at the taskbar and type "services",


9. Open browser and check if the webserver has been running, the typical URL is localhost:port/mapproxy/demo. In this page, you will see all the capabilities and protocols and endpoints of the map services. From this point, you can starts to seed the services to your GIS desktop, WebGIS, mobileGIS whatever.


10. A sneak peek.




Basemap is not so myths after all,

Stay Tuned


Friday, January 26, 2018

Modifying Vector Tile Style (Basemap series Part 4)

Please read my previous blog post about how to author vector tile.

Okay now I will continue about how to modify vector tile styling to suit custom needs.

we start from the vector tiles layer hosted on ArcGIS online.



what you must do first is copy the style using clik on layer name, and choose copy, because I use indonesian languange, it is called SALIN.


then followed by Save layer from the copied layer.


In order to edit the style of the copied layer, you can use two built in application from ArcGIS online, which is

1. Simple Editor . This is code based, so you can freely modify as long as you understand the syntax. The syntax is similar with CSS for your information.


2. Vector Basemap Style Editor. This tool is more GUI based and perform better in many cases.


Give the two apps permission to works on your arcGIS online account

Using these two applications, you can modifying the style of the tiled layer like changing colors, thicken the outline, make the layer transparent, changing font or any other modification permitted by Vector based tiled layer specifications. Dont forget to save your works.

Besides these two built in apps, you can edit the style using any other JSON supported text editor and validator (i.e notepad++, adobe brackets , psPad, or just plain NOTEPAD *lol)

Here are two screenshots of original and modified vector layer style coming from same vector tile package source.




Cheers!

Vector Tile (Basemap Series Part 3)

My two previous posts concerned about tile layers served as basemap in raster format. Tile layer in raster format is superb because it is deliver fast reading, great compatibility between systems and standards, Tile layer in raster format is also works for any kind of data and independent from the source data format. Either you have thematic maps in vector format or satellite imagery/DEM in raster format, Tile layer in raster could handle them well.

But technology always changing. In recent years, development of tile layer in vector format got its fruition. It is somewhat better than raster, but it is not dedicated to replace the raster based tiles. It has specific function which is mainly to serve basemap in VECTOR format.

By using vector format, you can get :
1. Smaller file size to serve, so higher performance is guaranteed.
2. Dynamic and good resolution/display in every scale or zoom level (vector ftw man).
3. Faster generation and processing time compared to raster.
4. No artifacts or jagged features due to resolution limitation like we commonly found on raster.
5. And the most important is, points are stored in points, line stored in line, polygon stored in polygon, and label stored in font (not pixels). In this way, vector tile enable the user to change the symbolization or label format of the tiles according to their needs. This is the feature that raster based tiles couldn't do.

The limitation of vector tiles is of course they can't server data sources stored in raster format like satellite imagery or Digital Elevation Model.

You can make vector tile using online cloud GIS provider like Mapbox, or ArcGIS pro.

Here is an example to authorize vector tile using ArcGIS pro and hosted under ArcGIS online.

1. Load up your map and its final symbolization and labeling in ArcGIS pro.


2. Open Geoprocessing tools > Data Management Tools > Packages , To make proper vector tile, you must make first the tile index, so open the Create Vector Tile Index Tool. Set it up according your needs, and go.



3. Now the tiles creation, load up Create Vector Tile Tool, and make your tiles, dont forget to set the tile index created from previous step as the index polygons. The output will be saved in VTPK format.


4. After it is done, upload the created data to ArcGIS Online under organization account (public account is not working). You can apply trial to ESRI if you dont have access to ArcGIS online subscription. Here is how the map looks in ArcGIS online. 





I will continue the tutorial about how to change the symbolization of vector tile in next post, stay tuned.


Btw, here is a LINK which all tools and theoretical background about vector tiles currently has been developed
https://github.com/mapbox/awesome-vector-tiles/



Thursday, January 25, 2018

Making Mbtiles in QGIS (Basemap series part 2)

I will continue my previous blog post.

Okay, now how to make the mbtiles of our data directly from our own maps.

In ArcGIS you can make Tile Package (TPK) which can be converted into other tiling scheme using certain tool (one of them has been covered in part 1)

In QGIS (or Global Mapper) you can actually generate directly tile packages to be consumed in WebGIS, either in compressed tile ZIP or Mbtiles format.

For this task, we going to need QTILES Plugin which can be downloaded and installed from QGIS Plugins menu.


And then, just open your QGIS map project (usually stored in qgs file) or just make it from the scratch.

QTILES can be accessed from Plugins drop down menu, once it is loaded, just set all the parameters according your needs.





Here is an example of mbtiles hosting using TILESERVER.PHP running on my local XAMPP machine.





stay tuned for the part 3, part 4, part 5. lols

Wednesday, January 24, 2018

ArcGIS Tile Package Conversion to Mbtiles using tpkutils (Basemap series part 1)

Imagine if all the published basemaps over the internet is not suit to your needs, what you can do?. The answer is simply, just make your own basemaps.

And how to make them ? by convert the created map into raster tile package of course.
From the created tile package, we can feed it off to map service provider like ArcGIS.com, Mapbox, etc, or just put it on your own server and webgis application/geoportal.

Basemap can be authored using ArcGIS or QGIS (just to name the BIG TWO software these days). I am going to talk about ArcGIS first for this post.

in ESRI ArcGIS environment, the tile package is called ArcGIS Tile Layer Package, and has *.tpk file extension.

How to make them?, simply
1. Make your own map in ArcMap (symbolization, labelling, etc and save it on MXD.
2. Boot ArcToolbox, go to Data Management Tools > Packages > Create Map Tile Package
3. Set the tiling scheme (just use arcgis online/google tiling scheme), tiling format, level or details, tag, summary, extent and output location.

In order to make the maxium benefit of the tpk, load it on your arcgis online organization, or arcgis server. The question is what if we dont have access to them? and what if we want to use the package to our own custom webgis application?. We could go to MBTILES. MBTILES can be made directly on Global Mapper or may be other software, but not in the case if we use ArcGIS (there is only TPK man).

There are few ways to convert TPK to MBTILES, one of them is using combination of PBS and MOBAC, but it is tedious. As alternative we can use TPKUTILS available on github in this LINK.

How to do it?

1. Clone the tool on Github using, PyPI, Github desktop or any other tool you familiarize with. Here is an example using PyPI by my friend Muhammad Iqnaul Siregar


2. Download Phyton 3 or newer release, and install it.
3. Check if PIP has been installed together with other python lib by typing PIP command in scripts directory of Python installation. Here is an example.



3. If it succeeded, Install the TPKUTILS using PIP INSTALL in phyton command line interface. Here is an example.
Note that I placed the downloaded TPKUTILS on my python installation directory, you can place it anywhere, just make sure windows has access r/w to the directory,

4. To keep it simple, copy your TPK file to your python installation directory, the run the python CLI. 
5. Call the tpkutils using command:
from tpkutils import TPK
6. continue with this command:
tpk = TPK('yourtpk.tpk')
7. followed by this command and you are done.
tpk.to_mbtiles('yourtpk.mbtiles')
8. Here is an example of the complete code and the result




9. And here is how it looks when I put the converted MBTILES to TILESERVER running on my local XAMPP machine.



Credit to the GIS Jedi master Marhensa Aditya Hadi and Muhammad Iqnaul Siregar
Also thanks to Benny Istanto which introduce me to Vector and Raster Tile Package back then in 2011.

Adios.



Sunday, January 21, 2018

download extract data from non hosted arcgis online public account

how to download arcgis online data from non hosted or public account



For Your Information



This extraction method is purely dedicated only to extract data stored in public

Account or non hosted ArcGIS online service



Non hosted data is pretty much different from hosted data, because it doesnt

have REST endpoint (Map Services) like commonly we found at ArcGIS

Server/online hosted data.



To be successfully extract the data, you must be the account owner (i.e the data

belongs to you), or the webmap has been set to public access by the account

owner.



if you still confused about non hosted data, try to make arcgis online public

account (not organizational account), and upload some data to it, then make a

webmap. The data stored in those webmap that we called non hosted data.