create grid and reproject

[moh713...@gmail.com]

Dear List,

as you can see the image attached, the study area(black dash line : it is an island) is covered by 4 sentinel-2 tiles. I have two inquires:

1- how can I create the ARD with less number of tiles (design the grid)

2- how can I edit the projection in the parameter file for the data cube creation?

 I hope this is very clear.

I use the latest version of FORCE on ubuntu 18.04

thanks in advance

[Stefan]

Dear Mohammed,

this is well documented in the tutorials.

https://force-eo.readthedocs.io/en/latest/howto/datacube.html

Stefan

[moh713...@gmail.com]

Hello Stefan,

this is what really happens 

The data cube parameters :

# DATA CUBES

# ------------------------------------------------------------------------

DO_REPROJ = TRUE

DO_TILE = TRUE

FILE_TILE = /home/localuser/force/mis/allow-list.txt

TILE_SIZE = 30000

BLOCK_SIZE = 3000

RESOLUTION_LANDSAT = 30

RESOLUTION_SENTINEL2 = 10

ORIGIN_LON =  12     

ORIGIN_LAT =  53 

PROJECTION = PROJCRS["ProjWiz_Custom_Cylindrical_Equal_Area",BASEGEOGCRS["WGS 84",DATUM["World Geodetic System 1984",ELLIPSOID["WGS 84",6378137,298.257223563, LENGTHUNIT["metre",1]],

            ID["EPSG",6326]],

        PRIMEM["Greenwich",0,

            ANGLEUNIT["Degree",0.0174532925199433]]],

    CONVERSION["unnamed",

        METHOD["Lambert Cylindrical Equal Area (Spherical)",

            ID["EPSG",9834]],

        PARAMETER["Latitude of 1st standard parallel",6.3630421,

            ANGLEUNIT["Degree",0.0174532925199433],

            ID["EPSG",8823]],

        PARAMETER["Longitude of natural origin",53.9222717,

            ANGLEUNIT["Degree",0.0174532925199433],

            ID["EPSG",8802]],

        PARAMETER["False easting",0,

            LENGTHUNIT["metre",1],

            ID["EPSG",8806]],

        PARAMETER["False northing",0,

            LENGTHUNIT["metre",1],

            ID["EPSG",8807]]],

    CS[Cartesian,2],

        AXIS["(E)",east,

            ORDER[1],

            LENGTHUNIT["metre",1,

                ID["EPSG",9001]]],

        AXIS["(N)",north,

            ORDER[2],

            LENGTHUNIT["metre",1,

                ID["EPSG",9001]]]]

RESAMPLING = NN

and this is the log file after running

/home/localuser/force/level1/sentinel/1-2017/T39PZQ/S2A_MSIL1C_20170101T064252_N0204_R120_T39PZQ_20170101T064359.SAFE

: Transformation failed.

Computing tile origin in dst_srs failed. Starting datacube failed.

Starting datacube(s) failed.

so what is going wrong?

thanks in advance

[david frantz]

Hi Mohammed, 

please check your settings for 

ORIGIN_LON =  12     

ORIGIN_LAT =  53 

I suspect that you have swapped X and Y.

Best,
David

[moh713...@gmail.com]

Hi DAVID,

I have gone through the original code :

cube->origin_geo.x = pl2->orig_lon;

    cube->origin_geo.y = pl2->orig_lat;

    if ((warp_geo_to_any(cube->origin_geo.x,  cube->origin_geo.y,

                        &cube->origin_map.x, &cube->origin_map.y,

                         cube->proj)) == FAILURE){

      printf("Computing tile origin in dst_srs failed. ");

      free_multicube(multicube);

      return NULL;

    }

 trying to explore how I swapped X and Y but without getting progress to solve the problem.could you clarify it please?

[moh713...@gmail.com]

hopefully with this image , the problem will be more clear

[david frantz]

Dear Mohammed, 

this has nothing to do with the code. Each projection is only valid within defined limits.

Your area of interest is at 12° Latitude and 53° Longitude.

In the parameter file, you defined 12°Longitude and 53° latitude (this is somewhere in Germany)

X = Longitude, Y = Latitude.

[moh713...@gmail.com]

Dear David,

no things has been changed !!

[david frantz]

1) Did you validate the projection?

2) if you cannot make your projection work, use one of the built-in projections, e.g. GLANCE7-AF

[moh713...@gmail.com]

1- it is validated well

Validate Succeeds

PROJ.4 : +proj=cea +lat_ts=6.3630421 +lon_0=53.9222717 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs

OGC WKT2:2018 :

PROJCRS["ProjWiz_Custom_Cylindrical_Equal_Area",

    BASEGEOGCRS["WGS 84",

        DATUM["World Geodetic System 1984",

2-   that already worked well but I want to use this projection to keep minimum distortion in area as  I plan to create land cover for this area .

hopefully this issue will be solved

thanks

[moh713...@gmail.com]

Hi David,

After I validated the projection, I rewrote it in parameter file as below and it worked :

PROJECTION = PROJCRS["ProjWiz_Custom_Cylindrical_Equal_Area",BASEGEOGCRS["WGS 84",DATUM["World Geodetic System 1984",ELLIPSOID["WGS 84",6378137,298.257223563, LENGTHUNIT["metre",1]],ID["EPSG",6326]],PRIMEM["Greenwich",0,ANGLEUNIT["Degree",0.0174532925199433]]],CONVERSION["unnamed",METHOD["Lambert Cylindrical Equal Area (Spherical)",ID["EPSG",9834]],PARAMETER["Latitude of 1st standard parallel",6.3630421,ANGLEUNIT["Degree",0.0174532925199433],ID["EPSG",8823]],PARAMETER["Longitude of natural origin",53.9222717,ANGLEUNIT["Degree",0.0174532925199433],ID["EPSG",8802]],PARAMETER["False easting",0,LENGTHUNIT["metre",1],ID["EPSG",8806]],PARAMETER["False northing",0,LENGTHUNIT["metre",1],ID["EPSG",8807]]],CS[Cartesian,2],AXIS["(E)",east,ORDER[1],LENGTHUNIT["metre",1,ID["EPSG",9001]]],AXIS["(N)",north,ORDER[2],LENGTHUNIT["metre",1,ID["EPSG",9001]]]]

a new issue appeared when I checked the tiles as the most eastern part of the study area is not covered by the tiles as shown below :

I tried to check if the problem related to the projection choosing one from  the built-in projections but the problem is the same

thanks in advance.

[david frantz]

Hi Mohammed,

I am not really sure what problem you are referring to.. Please elaborate more. From your image, I do not see any problem. Your study area is at the edge of a swath.. This means, the entire island will never be observed on the same day.

Best,

David

[moh713...@gmail.com]

Hi David,

I used multiple images with different dates with a cloud cover less than 75%  covering all the study area(four tiles covering the whole  island as shown at the first message in this conversation).

the problem is that all tiles should  have to be produced to cover the whole study area as logically expected  and based on the input data, unfortunately the most right part of the study area has not been covered by the produced tiles as expected despite changing the projection, the cloud cover percentage of the processed images.

hopefully this is now very clear

[Stefan]

Dear Mohammed,

see how the Sentinel-2 tiling grid relates to the Sentinel-2 orbits:

https://sentinel.esa.int/web/sentinel/user-guides/sentinel-2-msi/revisit-coverage and https://sentinel.esa.int/web/sentinel/missions/sentinel-2/data-products

While FORCE provides a lot of easy to use tools, it does rely on its users having a solid base knowledge of remote sensing concepts and the particularities of the data used. In this free support group here we are trying to do our best to answer questions that are directly related to the functionality of FORCE. Providing support outside of this scope is unfortunately beyond our possibilites.

Stefan

[moh713...@gmail.com]

Dear David, 

thank you for all the efforts you do for FORCE users.

I will do my best to see to which kind of problems this issue is related .

I will also try in another area to confirm from where this issue arise?

thanks again for  all your efforts  for the scientific community. 

[moh713...@gmail.com]

Dear Stefan,

I think there is a misunderstanding to what I meant in my issue reporting.

I used in my process  sentinel-2 images that cover all the study area with characteristics that have already mentioned but after processing No tiles covering the right corner of the study area. I used all the multidate produced tiles to discover that issue.

I checked the cloud cover of the images to see if that is related to the parameters we set in the parameter file, I change the projection to the one in-built to see again but the issue is still .

hopefully this is very clear now

[Stefan]

If you are 100% sure that you had L1 images from the orbit to the east, it's most likely that the L2 processing of these images failed due to high cloud cover / low registration targets. You can check the log files to see which scenes were processed successfully and which failed.

[moh713...@gmail.com]

I really took into considerations all what mentioned and the log file did not report any issue regarding failing due to high cloud cover / low registration targets instead it reports this :

dc:  66.02%. wc:  66.14%. sc:   0.02%. cc:   0.36%. AOD: 0.1991. # of targets: 170/2.  1 product(s) written. Success! Processing time: 05 mins 57 secs

I will do the process again and see what is the issue. and report here.

thanks again

[david frantz]

Dear Mohammed,

please indicate how many images you processed, and please post the original scene IDs.

Thanks,

David

[moh713...@gmail.com]

Dear David,

Those are the original scene IDs : 

T39PYP

 T39PYQ

T39PZP

T39PZQ

Regarding the numbers of images I processed, I used four images covering the whole study area >

images with less 75 % cloud cover and  I tried with less than 5 % cloud cover when  I did not see the tiles covering the eastern part of the study area and the parameter file setting allow producing tiles to 75 % cloud cover.

when failed because the coregistration, I disable the coregistration in the parameter file but the problem is still.

Thanks

[Stefan]

Mohammed,

please download the shapefile and kml from the links I have posted above. As mentioned before: it is not possible to cover your study area with S2 images from one single day.

What you posted are sentinel tile identifiers, not scene IDs. Again: check how orbits and S2 tiles overlay.

[moh713...@gmail.com]

hi

so sorry  for this confusion.

some thing has gone wrong with me !!

I am reprocess the images again.

so sorry again for this confusion