Q&A - Question

Here you find a c# procedure that takes a TGIS_LayerGRD and an array of class and returns a TGIS_LayerVector, hopeing It may help someone.

Usage inside the main program having four classes: 0-0.5;0.5-1;1-2;>2.

Please notice it may work on Grid files having header in the form of

ncols       
nrows       
xllcorner   
yllcorner   
cellsize     
nodata_value

and does not distinguish between no-data-value and 0 (due GetGrid function).
 

			// Vettore con le classi valMin,valMax,classe
			float[][] classi = new float[4][];
			for (int i = 0; i < 4; i++)
			{
				classi[i] = new float[3];
			}
			// Classe 1
			classi[0][0] = 0;
			classi[0][1] = 0.5f;
			classi[0][2] = 1;

			// Classe 2
			classi[1][0] = 0.5f;
			classi[1][1] = 1;
			classi[1][2] = 2;

			// Classe 3
			classi[2][0] = 1;
			classi[2][1] = 2;
			classi[2][2] = 3;

			// Classe 4
			classi[3][0] = 2;
			classi[3][1] = 999999;
			classi[3][2] = 4;

			TGIS_LayerVector vectorizedRaster = Vettorizza(raster, classi);

		private TGIS_LayerVector Vettorizza(TGIS_LayerGRD locRaster, float[][] arrClassi)
		{
			TGIS_LayerVector vettoriale;

			// Carico il locRaster nella matrice
			double xllcorner = 100;
			double yllcorner = 1000;
			double cellsize = 10;
			int ncols = 16;
			int nrows = 16;
			float[][] grid = RiempiMatrice();

			vettoriale = new TGIS_LayerVector();
			vettoriale.AddField("valore", TGIS_FieldType.Float, 12, 6);
			vettoriale.SetCSByEPSG(locRaster.CS.EPSG);
			locRaster.Viewer.Ref.Add(vettoriale);
			vettoriale.Lock();
			locRaster.Viewer.Ref.Lock();

			TGIS_Topology tpl = new TGIS_Topology();

			string s = "";

			for (int c = 0; c < ncols; c++)
			{
				for (int r = 0; r < nrows; r++)
				{

					for(int k=0;k< arrClassi.Length; k++)
					{
						if (grid[r][c] > arrClassi[k][0] && grid[r][c] <= arrClassi[k][1])
						{
							CreaQuadratiContigui(r, c, arrClassi[k][0], arrClassi[k][1], arrClassi[k][2]);
							UnisciQuadrati();
						}
					}
				}
			}
			vettoriale.Unlock();
			locRaster.Viewer.Ref.Unlock();

			void UnisciQuadrati()
			{

				TGIS_Shape shpUnione;
				shpUnione = tpl.UnionOnList(lShp, true);
				vettoriale.AddShape(shpUnione);
				for (int k = lShp.Count - 1; k >= 0; k--)
				{
					TGIS_Shape shpToDelete = (TGIS_Shape)lShp[k];
					shpToDelete.Delete();
				}
				lShp.Clear();
			}

			float[][] RiempiMatrice()
			{
				xllcorner = 100;
				yllcorner = 1000;
				cellsize = 10;
				ncols = 16;
				nrows = 16;

				using (StreamReader sr = File.OpenText(locRaster.Path))
				{
					string x = "";
					int righe = 0;
					while ((x = sr.ReadLine()) != null)
					{
						if (x.ToLower().IndexOf("xllcorner") != -1)
						{
							xllcorner = float.Parse(x.Substring(11, x.Length - 11));
						}
						if (x.ToLower().IndexOf("yllcorner") != -1)
						{
							yllcorner = float.Parse(x.Substring(11, x.Length - 11));
						}
						if (x.ToLower().IndexOf("cellsize") != -1)
						{
							cellsize = float.Parse(x.Substring(11, x.Length - 11));
						}
						if (x.ToLower().IndexOf("ncols") != -1)
						{
							ncols = int.Parse(x.Substring(11, x.Length - 11));
						}
						if (x.ToLower().IndexOf("nrows") != -1)
						{
							nrows = int.Parse(x.Substring(11, x.Length - 11));
						}

						// esco dopo aver letto lo header
						righe++;
						if (righe > 6)
						{
							break;
						}

					}
				}

				float[][] locGrid = new float[nrows][];
				for (int i = 0; i < nrows; i++)
				{
					locGrid[i] = new float[ncols];
				}

				yllcorner = yllcorner + (cellsize * nrows);

				locRaster.GetGrid(locRaster.Extent, locGrid);

				return locGrid;

			}

			void CreaQuadratiContigui(int r, int c, float valMin, float valMax, float classe)
			{
				// Creo il quadrato GIS
				TGIS_Shape shpQuadrato = vettoriale.CreateShape(TGIS_ShapeType.Polygon);
				TGIS_Point ptgQ;

				shpQuadrato = vettoriale.CreateShape(TGIS_ShapeType.Polygon);
				shpQuadrato.AddPart();

				ptgQ.X = xllcorner + (cellsize * c);
				ptgQ.Y = yllcorner - (cellsize * r);
				shpQuadrato.AddPoint(ptgQ);

				ptgQ.X = xllcorner + (cellsize * c);
				ptgQ.Y = yllcorner - (cellsize * (r + 1));
				shpQuadrato.AddPoint(ptgQ);

				ptgQ.X = xllcorner + (cellsize * (c + 1));
				ptgQ.Y = yllcorner - (cellsize * (r + 1));
				shpQuadrato.AddPoint(ptgQ);

				ptgQ.X = xllcorner + (cellsize * (c + 1));
				ptgQ.Y = yllcorner - (cellsize * r);
				shpQuadrato.AddPoint(ptgQ);

				shpQuadrato.SetField("valore", classe);
				lShp.Add(shpQuadrato);

				// invalido la cella corrente
				grid[r][c] = -1;

				// verifico se nelle quattro direzioni ho celle buone
				// Top
				if (r > 0 && grid[r - 1][c] > valMin && grid[r - 1][c] <= valMax)
				{
					CreaQuadratiContigui(r - 1, c, valMin, valMax, classe);
				}
				// Right
				if (c < ncols - 1 && grid[r][c + 1] > valMin && grid[r][c + 1] <= valMax)
				{
					CreaQuadratiContigui(r, c + 1, valMin, valMax, classe);
				}
				// Bottom
				if (r < nrows - 1 && grid[r + 1][c] > valMin && grid[r + 1][c] <= valMax)
				{
					CreaQuadratiContigui(r + 1, c, valMin, valMax, classe);
				}
				// Left
				if (c > 1 && grid[r][c - 1] > valMin && grid[r][c - 1] <= valMax)
				{
					CreaQuadratiContigui(r, c - 1, valMin, valMax, classe);
				}
			};

			return vettoriale;
		}

private TatukGIS.RTL.TObjectList<object> lShp = new TatukGIS.RTL.TObjectList<object>();