[GIS] Displaying netCDF data with correct CRS

gdalnetcdfqgisraster

I successfully extracted some data from a netCDF file with a single variable E. (the original file was taking ages to load).
I used the following code to extract the data:

ncks -d lon,24.,36. -d lat,21.,32. E_2015.nc extracted_E.nc

The resulting file I open using the popular python plugin "NetCDF Browser" by Mr Etienne Tourigny.
Usually, the file opens and displays fantastically. However, I have come across some that have rather peculiar behaviour. (Please see below)

netCDF data in B/W Vector layer of where the data should be in purple

The purple layer is a vector layer for correct orientation reference, the netCDF data appears to be mirrored (with longitudes/latitudes switched and transversed).

gdalinfo gives me :

Driver: netCDF/Network Common Data Format
Files: extracted_E.nc
Size is 44, 48
Coordinate System is `'
Origin = (32.000000000000000,36.000000000000000)
Pixel Size = (-0.250000000000000,-0.250000000000000)
Metadata:
  E#long_name=Some Variable
  E#standard_name=Some Variable
  E#units=mm/day
  E#_FillValue=-999
  lat#long_name=latitude
  lat#standard_name=latitude
  lat#units=degrees_north
  lat#_FillValue=-999
  lon#long_name=longitude
  lon#standard_name=longitude
  lon#units=degrees_east
  lon#_FillValue=-999
  NC_GLOBAL#history=Sun Nov 20 16:14:38 2016: ncks -d lon,24.,36. -d lat,21.,32. E_2015.nc extracted_E.nc
  NC_GLOBAL#NCO="4.5.4"
  NETCDF_DIM_EXTRA={time}
  NETCDF_DIM_time_DEF={365,6}
  NETCDF_DIM_time_VALUES={16436,16437,16438,16439,16440,16441,16442,16443,16444,16445,16446,16447,16448,16449,16450,16451,16452,16453,16454,16455,16456,16457,16458,16459,16460,16461,16462,16463,16464,16465,16466,16467,16468,16469,16470,16471,16472,16473,16474,16475,16476,16477,16478,16479,16480,16481,16482,16483,16484,16485,16486,16487,16488,16489,16490,16491,16492,16493,16494,16495,16496,16497,16498,16499,16500,16501,16502,16503,16504,16505,16506,16507,16508,16509,16510,16511,16512,16513,16514,16515,16516,16517,16518,16519,16520,16521,16522,16523,16524,16525,16526,16527,16528,16529,16530,16531,16532,16533,16534,16535,16536,16537,16538,16539,16540,16541,16542,16543,16544,16545,16546,16547,16548,16549,16550,16551,16552,16553,16554,16555,16556,16557,16558,16559,16560,16561,16562,16563,16564,16565,16566,16567,16568,16569,16570,16571,16572,16573,16574,16575,16576,16577,16578,16579,16580,16581,16582,16583,16584,16585,16586,16587,16588,16589,16590,16591,16592,16593,16594,16595,16596,16597,16598,16599,16600,16601,16602,16603,16604,16605,16606,16607,16608,16609,16610,16611,16612,16613,16614,16615,16616,16617,16618,16619,16620,16621,16622,16623,16624,16625,16626,16627,16628,16629,16630,16631,16632,16633,16634,16635,16636,16637,16638,16639,16640,16641,16642,16643,16644,16645,16646,16647,16648,16649,16650,16651,16652,16653,16654,16655,16656,16657,16658,16659,16660,16661,16662,16663,16664,16665,16666,16667,16668,16669,16670,16671,16672,16673,16674,16675,16676,16677,16678,16679,16680,16681,16682,16683,16684,16685,16686,16687,16688,16689,16690,16691,16692,16693,16694,16695,16696,16697,16698,16699,16700,16701,16702,16703,16704,16705,16706,16707,16708,16709,16710,16711,16712,16713,16714,16715,16716,16717,16718,16719,16720,16721,16722,16723,16724,16725,16726,16727,16728,16729,16730,16731,16732,16733,16734,16735,16736,16737,16738,16739,16740,16741,16742,16743,16744,16745,16746,16747,16748,16749,16750,16751,16752,16753,16754,16755,16756,16757,16758,16759,16760,16761,16762,16763,16764,16765,16766,16767,16768,16769,16770,16771,16772,16773,16774,16775,16776,16777,16778,16779,16780,16781,16782,16783,16784,16785,16786,16787,16788,16789,16790,16791,16792,16793,16794,16795,16796,16797,16798,16799,16800}
  time#calender=standard
  time#long_name=time
  time#standard_name=time
  time#units=days since 1970-01-01 00:00:00 UTC
  time#_FillValue=-999
Corner Coordinates:
Upper Left  (  32.0000000,  36.0000000) 
Lower Left  (  32.0000000,  24.0000000) 
Upper Right (  21.0000000,  36.0000000) 
Lower Right (  21.0000000,  24.0000000)
(followed by band data)

Both the original netCDF and the extracted display correctly in netCDF browsers/viewers such as PanoplyJ. Other extracted netCDF files also display correctly in QGIS. Is there something special about this file? How can I fix this?

I have tried gdalwarp:

gdalwarp -t_srs '+proj=lonlat +datum=WGS84 +no_defs t_srs EPSG:4326' -of E_2015.nc x.nc

I get:

ERROR 1: Translating source or target SRS failed:

Furthurmore, 

gdal_translate -a_ullr 24.0 32.0 36.0 21.0 -a_srs 'EPSG:4326' -of netCDF E_2015.nc x.nc

Gives:

Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute

With the resulting image as shown below.

enter image description here

B/W file is the new data. This needs to be at right angles to what it is.

The resulting file when treated with the same code as above returns the original rotation (see first image above), but bound within the area of the purple vector reference. It is almost as if the Image is rotated around those two points the (ulx,uly) and (lrx,lry) overridden by gdal_translate (similar to flipping a card around it's diagonals). Also no amount of swapping x and y, or entering the opposing diagonal points in the gdal_translate command results in the correct orientation(on either the original data, or the result of gdal_translate).

I've also attempted to convert to Geotiff and reproject with no luck.

Best Answer

You can use the netCDF operators (NCO) to add spatial reference information that will be understood by QGIS. Essentially, we need to:

  1. Create a crs variable that has no associated dimensions
  2. Attach attributes to the crs variable that describe our coordinate system
  3. Use a grid_mapping attribute to associate the variables that contain spatial data to attributes defined under crs.

Here's an example in WGS84:

ncap -h -O -s 'crs=-9999' my_data.nc my_data.nc
ncatted -h -O \
    -a spatial_ref,crs,c,c,'GEOGCS[\"GCS_WGS_1984\",DATUM[\"WGS_1984\",SPHEROID[\"WGS_84\",6378137.0,298.257223563]],PRIMEM[\"Greenwich\",0.0],UNIT[\"Degree\",0.017453292519943295]]' \
    -a grid_mapping_name,crs,c,c,'latitude_longitude' \
    -a longitude_of_prime_meridian,crs,c,d,0 \
    -a semi_major_axis,crs,c,d,6378137 \
    -a inverse_flattening,crs,c,d,298.257223563 \
    -a grid_mapping,my_var_1,c,c,'crs' \
    -a grid_mapping,my_var_2,c,c,'crs' \
    my_data.nc
Related Question