--- netcdfdataset.cpp.orig	2011-04-15 14:49:17.000000000 -0300
+++ netcdfdataset.cpp	2011-04-15 19:44:13.000000000 -0300
@@ -485,6 +485,9 @@
             case NC_BYTE:
                 /* don't do default fill-values for bytes, too risky */
                 dfNoData = 0.0;
+                /* from etiennesky: why is this risky? use default value, 
+                   if not values that are 0.0 will be treated as missing */
+                /* dfNoData = NC_FILL_BYTE; */
                 break;
             case NC_CHAR:
                 dfNoData = NC_FILL_CHAR;
@@ -2402,10 +2405,10 @@
 
     GDALDataType eDT;
 
-    status = nc_def_dim( fpImage, "x", nXSize, &nXDimID );
+    status = nc_def_dim( fpImage, "lon", nXSize, &nXDimID );
     CPLDebug( "GDAL_netCDF", "status nc_def_dim X = %d\n", status );
 
-    status = nc_def_dim( fpImage, "y", nYSize, &nYDimID );
+    status = nc_def_dim( fpImage, "lat", nYSize, &nYDimID );
     CPLDebug( "GDAL_netCDF", "status nc_def_dim Y = %d\n", status );
 
     CPLDebug( "GDAL_netCDF", "nYDimID = %d\n", nXDimID );
@@ -2440,13 +2443,22 @@
         double adfGeoTransform[6];
         char   szGeoTransform[ MAX_STR_LEN ];
         char   szTemp[ MAX_STR_LEN ];
+        char   szGdalGeo[ MAX_STR_LEN ];
 
-
-/* -------------------------------------------------------------------- */
-/*      Copy GeoTransform array from source                             */
-/* -------------------------------------------------------------------- */
+        double dfX0=0.0;
+        double dfDX=0.0;
+        double dfY0=0.0;
+        double dfDY=0.0;
+        double *pafLonLat  = NULL;
+        int    anLatLonDims[1];
+        size_t startLonLat[1];
+        size_t countLonLat[1];
+
+	/* -------------------------------------------------------------------- */
+	/*      Copy GeoTransform array from source                             */
+	/* -------------------------------------------------------------------- */
         poSrcDS->GetGeoTransform( adfGeoTransform );
-        bWriteGeoTransform = TRUE;
+        bWriteGeoTransform = FALSE;
 
         *szGeoTransform = '\0';
         for( i=0; i<6; i++ ) {
@@ -2455,69 +2467,212 @@
             strcat( szGeoTransform, szTemp );
         }
         CPLDebug( "GDAL_netCDF", "szGeoTranform = %s", szGeoTransform );
-        strcpy( pszNetcdfProjection, "GDAL_Geographics" );
 
-        status = nc_def_var( fpImage, 
-                             pszNetcdfProjection, 
-                             NC_CHAR, 
-                             0, NULL, &NCDFVarID );
-	
         dfNN = adfGeoTransform[3];
         dfSN = ( adfGeoTransform[5] * nYSize ) + dfNN;
         dfWE = adfGeoTransform[0];
         dfEE = ( adfGeoTransform[1] * nXSize ) + dfWE;
-
-        status = nc_put_att_double( fpImage,
-                                    NCDFVarID, 
-                                    "Northernmost_Northing",
-                                    NC_DOUBLE,
-                                    1,
-                                    &dfNN );
-        status = nc_put_att_double( fpImage,
-                                    NCDFVarID, 
-                                    "Southernmost_Northing",
-                                    NC_DOUBLE,
-                                    1,
-                                    &dfSN );
-        status = nc_put_att_double( fpImage,
-                                    NCDFVarID,
-                                    "Easternmost_Easting",
-                                    NC_DOUBLE,
-                                    1,
-                                    &dfEE );
-        status = nc_put_att_double( fpImage,
-                                    NCDFVarID,
-                                    "Westernmost_Easting",
-                                    NC_DOUBLE,
-                                    1,
-                                    &dfWE );
         pszWKT = (char *) poSrcDS->GetProjectionRef() ;
 
-        nc_put_att_text( fpImage, 
-                         NCDFVarID, 
-                         "spatial_ref",
-                         strlen( pszWKT ),
-                         pszWKT );
+	/* -------------------------------------------------------------------- */
+	/*      Write GDAL_Geographics attributes                               */
+	/* -------------------------------------------------------------------- */
 
-        nc_put_att_text( fpImage, 
-                         NCDFVarID, 
-                         "GeoTransform",
-                         strlen( szGeoTransform ),
-                         szGeoTransform );
+	if( bWriteGeoTransform == TRUE ) {
+	    
+            strcpy( pszNetcdfProjection, "GDAL_Geographics" );
+            nc_def_var( fpImage, 
+                        pszNetcdfProjection, 
+                        NC_CHAR, 
+                        0, NULL, &NCDFVarID );
+	    
+            nc_put_att_double( fpImage,
+                               NCDFVarID, 
+                               "Northernmost_Northing",
+                               NC_DOUBLE,
+                               1,
+                               &dfNN );
+            nc_put_att_double( fpImage,
+                               NCDFVarID, 
+                               "Southernmost_Northing",
+                               NC_DOUBLE,
+                               1,
+                               &dfSN );
+            nc_put_att_double( fpImage,
+                               NCDFVarID,
+                               "Easternmost_Easting",
+                               NC_DOUBLE,
+                               1,
+                               &dfEE );
+            nc_put_att_double( fpImage,
+                               NCDFVarID,
+                               "Westernmost_Easting",
+                               NC_DOUBLE,
+                               1,
+                               &dfWE );
+	    
+            nc_put_att_text( fpImage, 
+                             NCDFVarID, 
+                             "spatial_ref",
+                             strlen( pszWKT ),
+                             pszWKT );
+            nc_put_att_text( fpImage, 
+                             NCDFVarID, 
+                             "GeoTransform",
+                             strlen( szGeoTransform ),
+                             szGeoTransform );
+            nc_put_att_text( fpImage, 
+                             NCDFVarID, 
+                             GRD_MAPPING_NAME,
+                             30,
+                             "Geographics Coordinate System" );
+            nc_put_att_text( fpImage, 
+                             NCDFVarID, 
+                             LNG_NAME,
+                             14,
+                             "Grid_latitude" );
+	}
+	
+        /* write geo info to global attribute GDAL_Geographics */
+        *szGdalGeo = '\0';
+
+        sprintf( szTemp, "Northernmost_Northing = %.15g \n",
+                 dfNN );
+        strcat( szGdalGeo, szTemp );
+        sprintf( szTemp, "Southernmost_Northing = %.15g \n",
+                 dfSN );
+        strcat( szGdalGeo, szTemp );
+        sprintf( szTemp, "Easternmost_Easting = %.15g \n",
+                 dfEE );
+        strcat( szGdalGeo, szTemp );
+        sprintf( szTemp, "Westernmost_Easting = %.15g \n",
+                 dfWE );
+        strcat( szGdalGeo, szTemp );
+
+        sprintf( szTemp, "spatial_ref = %s \n",
+                 pszWKT );
+        strcat( szGdalGeo, szTemp );
+        sprintf( szTemp, "GeoTransform = %s \n",
+                 szGeoTransform );
+        strcat( szGdalGeo, szTemp );
+        sprintf( szTemp, "%s = Geographics Coordinate System \n",
+                 GRD_MAPPING_NAME );
+        strcat( szGdalGeo, szTemp );
+        sprintf( szTemp, "%s = Grid_latitude",
+                 LNG_NAME );
+        strcat( szGdalGeo, szTemp );
 
         nc_put_att_text( fpImage, 
-                         NCDFVarID, 
-                         GRD_MAPPING_NAME,
-                         30,
-                         "Geographics Coordinate System" );
-        nc_put_att_text( fpImage, 
-                         NCDFVarID, 
-                         LNG_NAME,
-                         14,
-                         "Grid_latitude" );
+                         NC_GLOBAL, 
+                         "GDAL_Geographics",
+                         strlen( szGdalGeo ),
+                         szGdalGeo );
+
+        /* -------------------------------------------------------------------- */
+        /*      Write latitude attributes                                       */
+        /* -------------------------------------------------------------------- */
+        anLatLonDims[0] = nYDimID;
+        status = nc_def_var( fpImage, "lat", NC_DOUBLE, 1, anLatLonDims, &NCDFVarID );
+        nc_put_att_text( fpImage,
+                         NCDFVarID,
+                         "standard_name",
+                         8,
+                         "latitude" );
+        nc_put_att_text( fpImage,
+                         NCDFVarID,
+                         "long_name",
+                         8,
+                         "latitude" );
+        nc_put_att_text( fpImage,
+                         NCDFVarID,
+                         "units",
+                         13,
+                         "degrees_north" );
+        nc_put_att_text( fpImage,
+                         NCDFVarID,
+                         "axis",
+                         1,
+                         "Y" );
+        
+        /* -------------------------------------------------------------------- */
+        /*      Write latitude values                                           */
+        /* -------------------------------------------------------------------- */
+        // dfY0 = adfGeoTransform[3];
+        /* invert latitude values */
+        dfY0 = ( adfGeoTransform[5] * nYSize ) + adfGeoTransform[3];
+        dfDY = adfGeoTransform[5];
+        
+        pafLonLat = (double *) CPLMalloc( nYSize * sizeof( double ) );
+        for( i=0; i<nYSize; i++ ) {
+                /* The data point is centered inside the pixel */
+                // pafLonLat[i] = dfY0 + (i+0.5)*dfDY ;
+                /* invert latitude values */
+                pafLonLat[i] = dfY0 - (i+0.5)*dfDY ;
+        }
+
+        startLonLat[0] = 0;
+        countLonLat[0] = nYSize;
+
+        /* Temporarily switch to data mode and write data */
+        status = nc_enddef( fpImage );
+        status = nc_put_vara_double( fpImage, NCDFVarID, startLonLat,
+                                            countLonLat, pafLonLat);
+        status = nc_redef( fpImage );
+
+
+        /* -------------------------------------------------------------------- */
+        /*      Write longitude attributes                                      */
+        /* -------------------------------------------------------------------- */
+        anLatLonDims[0] = nXDimID;
+        status = nc_def_var( fpImage, "lon", NC_DOUBLE,
+                                 1, anLatLonDims, &NCDFVarID );
+        nc_put_att_text( fpImage,
+                         NCDFVarID,
+                         "standard_name",
+                         9,
+                         "longitude" );
+        nc_put_att_text( fpImage,
+                         NCDFVarID,
+                         "long_name",
+                         9,
+                         "longitude" );
+        nc_put_att_text( fpImage,
+                         NCDFVarID,
+                         "units",
+                         12,
+                         "degrees_east" );
+        nc_put_att_text( fpImage,
+                         NCDFVarID,
+                         "axis",
+                         1,
+                         "X" );
+
+        /* -------------------------------------------------------------------- */
+        /*      Write longitude values                                          */      
+        /* -------------------------------------------------------------------- */
+        dfX0 = adfGeoTransform[0];
+        dfDX = adfGeoTransform[1];
+
+        pafLonLat = (double *) CPLRealloc( pafLonLat, nXSize * sizeof( double ) );
+        for( i=0; i<nXSize; i++ ) {
+                /* The data point is centered inside the pixel */
+                pafLonLat[i] = dfX0 + (i+0.5)*dfDX ;
+        }
+
+        startLonLat[0] = 0;
+        countLonLat[0] = nXSize;
+
+        /* Temporarily switch to data mode and write data */
+        status = nc_enddef( fpImage );
+        status = nc_put_vara_double( fpImage, NCDFVarID, startLonLat,
+                                        countLonLat, pafLonLat);
+        status = nc_redef( fpImage );
+
+        /* free lonlat values */
+        CPLFree( pafLonLat );
     }
 
-    if( oSRS.IsProjected() )
+    else if( oSRS.IsProjected() )
     {
         const char *pszParamStr, *pszParamVal;
         const OGR_SRSNode *poPROJCS = oSRS.GetAttrNode( "PROJCS" );
@@ -2734,7 +2889,12 @@
 /*      Write Data from Band i                                          */
 /* -------------------------------------------------------------------- */
                 start[0]=iLine;
-                start[1]=0;
+		/* invert latitude values */
+		/* must be done for projected also */
+		/* NOTE: only tested for BYTE and Int16 data types */
+        	if( oSRS.IsGeographic() ) 
+		    start[0]=nYSize - iLine - 1;
+		start[1]=0;
                 count[0]=1;
                 count[1]=nXSize;
 
@@ -2745,6 +2905,10 @@
                 status = nc_put_vara_uchar (fpImage, NCDFVarID, start,
                                             count, pabScanline);
 
+		if (status != NC_NOERR)
+		    CPLError( CE_Failure, CPLE_AppDefined, 
+			      "netCDF write failed: %s", 
+			      nc_strerror( status ) );
 
 /* -------------------------------------------------------------------- */
 /*      Put NetCDF file back in define mode.                            */
@@ -2781,6 +2945,11 @@
                                             0,0);
 
                 start[0]=iLine;
+		/* invert latitude values */
+		/* must be done for projected also */
+		/* NOTE: only tested for BYTE and Int16 data types */
+        	if( oSRS.IsGeographic() ) 
+		    start[0]=nYSize - iLine - 1;
                 start[1]=0;
                 count[0]=1;
                 count[1]=nXSize;
@@ -2823,6 +2992,11 @@
                                             0,0);
 
                 start[0]=iLine;
+		/* invert latitude values */
+		/* must be done for projected also */
+		/* NOTE: only tested for BYTE and Int16 data types */
+        	if( oSRS.IsGeographic() ) 
+		    start[0]=nYSize - iLine - 1;
                 start[1]=0;
                 count[0]=1;
                 count[1]=nXSize;
@@ -2864,6 +3038,11 @@
                                             0,0);
 
                 start[0]=iLine;
+		/* invert latitude values */
+		/* must be done for projected also */
+		/* NOTE: only tested for BYTE and Int16 data types */
+        	if( oSRS.IsGeographic() ) 
+		    start[0]=nYSize - iLine - 1;
                 start[1]=0;
                 count[0]=1;
                 count[1]=nXSize;
@@ -2905,6 +3084,11 @@
                                             0,0);
 
                 start[0]=iLine;
+		/* invert latitude values */
+		/* must be done for projected also */
+		/* NOTE: only tested for BYTE and Int16 data types */
+        	if( oSRS.IsGeographic() ) 
+		    start[0]=nYSize - iLine - 1;
                 start[1]=0;
                 count[0]=1;
                 count[1]=nXSize;