Scheduled Downtime
On Friday 21 April 2023 @ 5pm MT, this website will be down for maintenance and expected to return online the morning of 24 April 2023 at the latest

Corine Land Cover not taking into account by geogrid

Manuarii

Member
Hi,

I have a problem with my binary file created of corine land cover that seems to not be taking into account by geogrid. I have no warning or error on geogrid.log, but when plotting the geo_em.d02.nc for example, my Lu_INDEX have only 28 classes. It seems ok with default usage, but I have reclassify my CLC to USGS and I have some ID that goes up to 33. So I have to end up with 33 classes instead. The 31 to 33 representing the urban part. So the urban part is missing.To make my binary file I used a fortran code combine with the c code of write_geogrid.c. The output file seems ok, and I can see the data. Just in case, I have the same results, when putting empty binary file.

1- My reclassification is made on python with the reprojected tiff file, also reprojected on python. Here I reclassify with the following dictionnary :
reclassification_dict = {-128:0, 0:0, 48:0, 1:32, 2:31, 3:33, 4:33, 5:33, 6:33, 7:19, 8:19, 9:19, 10:7, 11:7, 12:2,
13:3, 14:3, 15:6, 16:6, 17:6, 18:2, 19:6, 20:4, 21:5, 22:6, 23:11, 24:14, 25:15, 26:7, 27:9,
28:9, 29:9, 30:19, 31:19, 32:19, 33:19, 34:24, 35:17, 36:17, 37:17, 38:17, 39:17, 40:16,
41:28, 42:28, 43:16, 44:16}

2- Next that new tiff file is send to the IDRIS server to convert the file to ASCII and into binary file. I run it on the server because of memory issues in my local computer. I used a slurm file containing the last step, where I do :

module load gdal/

gdal_translate -of AAIGrid cropped_CLC_wgs84_reclasssify.tif CLC_USGS.asc

gcc -D_UNDERSCORE -DBYTESWAP -DLINUX -DIO_NETCDF-DIO_BINARY-DIO_GRIB1 -D_GEOGRID -O -c write_geogrid.c

gfortran geogrid_clc.f90 write_geogrid.o

./a.out

-My index file contain :

type=categorical
row_order=top_bottom
category_min=1
category_max=33
projection=regular_ll
dx=0.00128354112745536
dy=0.00128354112745536
known_x=1.0
known_y=1.0
known_lat=68.3668663391942
known_lon=-31.8067559589092
wordsize=1
missing_value=0
tile_x=60602
tile_y=27641
tile_z=1
units="category"
description="CLC 2018 landuse"


-My namelist.wps contain : geog_data_res = 'srtm_3s+clc2006+usgs_lakes+default','srtm_3s+clc2006+usgs_lakes+default'

-And GEOGRID.TBL :

===============================
name=LANDUSEF
priority=2
dest_type=categorical
z_dim_name=land_cat
landmask_water = clc2006:16,28 # Calculate a landmask from this field
interp_option = clc2006:nearest_neighbor
rel_path= clc2006:specific/landuse_clc2018_wgs_test/
===============================

In case you ask, the problem is the same if I directly make a binary file from the tiff file. Also I have a reference to see what I have to obtain in my geo_em.d02.nc file for example. See figure attached. Don't pay attention to the levels of colorbar, I'm particularly interesting of the number of classes.


Is anyone have any idea ?
 

Attachments

  • geo_d02_ref.png
    geo_d02_ref.png
    262.8 KB · Views: 9
  • geo_d02_test.png
    geo_d02_test.png
    227.8 KB · Views: 9
Just to come back, now to be sure that my problem not come from binary convertion, I used the convert_geotiff function. But I still have the same issues and I'm not sure about where my problem come from and if I also use correctly convert_geotiff. To add more informations, I follow this full step :

1- Open CLC tiff file download of European project and crop it in a certain region with a python code.

2- Next I reproject on EPSG:4326 always with a custom python code.

3- Now I reclassify my data from this new tiff file reprojected, like said before, with custom python code.

4- Like I said before, I take the final tiff file that I put on IDRIS server and run the last step to convert in binary. And like I said, I try now with the convert_geotiff function to create my binary file and index file, to avoid any issues. Here I used convert_geotiff as follow :

./convert_geotiff -c 33 -t 10000 cropped_CLC_wgs84_reclassify.tif

I'm not sure if 10000 is correct, because the size of this tiff file is tile_x = 60602 and tile_y = 27641. Here 33 refer to the number of classes.
 
New Update! Now, with the functionality of the convert_geotiff function, I can accurately determine the number of classes required. However, a persisting problem remains. I have multiple binary files sharing the same window size, yet their land cover is consistently shifted compared to my reference plot, as illustrated in the accompanying figures. Addressing this spatial misalignment is the next focal point for resolution. I suspect that the mutiple binary file introduce this shifting. The command line that I'm currently using is : ./convert_geotiff -b 0 -c 33 -t 10000 OUT/cropped_CLC_wgs84_reclassify.tif.

The new index file generated is :

projection = regular_ll
known_x = 1
known_y = 27641
known_lat = 68.366867
known_lon = -31.806755
dx = 1.282648e-03
dy = 1.284267e-03
type = categorical
signed = yes
units = "NO UNITS"
description = "NO DESCRIPTION"
wordsize = 2
tile_x = 10000
tile_y = 10000
tile_z = 1
category_min = 1
category_max = 34
tile_bdr = 0
missing_value = 34.000000
scale_factor = 1.000000
row_order = bottom_top
endian = little
 

Attachments

  • geo_d02_ref_zoom.png
    geo_d02_ref_zoom.png
    333.2 KB · Views: 8
  • geo_d02_test_zoom.png
    geo_d02_test_zoom.png
    298.8 KB · Views: 8
Hi,
I would like to apologize for the long delay in response. Between holidays, unplanned time out of the office, the AMS conference, and planning for our WRF tutorial taking place this week, we got pretty behind on forum questions. Thank you so much for your patience.

Are you still experiencing issues with this?
 
Hi,

I'm pleased to inform you that everything has been successfully rectified. The issue stemmed from inaccuracies within my index file, where certain information was misaligned. Now, through the utilization of gdal_translate to convert my TIFF file to binary format, followed by renaming it to the specified format (00001-xxxxx.00001-xxxxx), the process works seamlessly.

Furthermore, I've developed a Python script to streamline the generation of my index file, ensuring its accuracy and conformity. Consequently, my file now accurately reflects the entire area of reclassified Corine Land Cover (CLC) data to United States Geological Survey (USGS) standards.

In summary, my current workflow involves the utilization of the following tools and processes:

gdal_translate -of ENVI OUT/CLC_wgs84_reclassify.tif output.bin

mv output.bin 00001-99999.00001-37778

python py_code/get_index_file.py -CLC "path_out"

Additional details: Upon obtaining the CLC file, I perform a reprojection to EPSG:4326 and apply a slight cropping to adjust for the x-axis size limitation, capping it at 99,999. Additionally, I adjust the starting point for latitude, known as know_lat, due to the original CLC file's known_y not being 1, but rather the height of the file, which is 37,778. Consequently, the resulting index file comprises the following specifications:

projection = regular_ll
known_x = 1
known_y= 1
known_lat= 24.2846829
known_lon = -56.5051419
dx = 0.001280602116034
dy = 0.001280602116034
type = categorical
signed = yes
units = "NO UNITS"
description = "NO DESCRIPTION"
wordsize = 1
tile_x = 99999
tile_y = 37778
tile_z = 1
category_min = 1
category_max = 41
tile_bdr = 0
missing_value = 0
scale_factor = 1.000000
row_order = top_bottom
endian = little
 
This is excellent news, and thank you so much for providing the script information. It may very well help another user in the future!
 
I got the Corine Land Cover working 2 years ago using a set ready to download. For more info see:
 
Top