title: KML in QGis with additional Data

text: i´ve developed a android app that extracts the exif-information of a photo and creates a KML-file (Placemark) - some additional information like camera-type, focal length,... are stored in the tag "ExtendedData". When i try to import the KML-file in qgis the attribute table only shows the name and description field. is it possible to import a the KML-file in qgis without loosing any information stored in "ExtendedData"? The KML-File is like this:               <?xml version="1.0" encoding="UTF-8"?>     <kml xmlns="http://www.opengis.net/kml/2.2">     <Placemark>     <name>20130404_132050</name>     <ExtendedData>     <Data name="Orientierung">     <displayName>Orientierung</displayName>     <value>1</value>     </Data>     <Data name="Brennweite">     <displayName>Brennweite</displayName>     <value>279/100</value>     </Data>     </ExtendedData>     <Point>     <coordinates>15.448561,47.11721,0</coordinates>     </Point>     </Placemark>     </kml>      Attributes in Qgis: ![enter image description here](http://i.stack.imgur.com/GY5Hb.jpg) Hope you can help me! thanks, michael EDIT: ok works fine now - i did not had the newest version of QGIS (2.0.1), the KML-File gets importet now correctly with all the attributes in "ExtendedData" ![enter image description here](http://i.stack.imgur.com/bzoVl.png)


title: Working with address data extracted from *.osm

text: My goal is get all the addresses from a certain country, ideally - from all countries, but for now - only from one. I downloaded an *.osm file from here http://download.geofabrik.de/ for a certain territory which contained a few countries and extracted the data from that file using Osm2pgsql to a Postgres database. There were no errors, all went well. But now I can't get along with those tables in the db, I can't figure out how actually to **query all the addresses** from a country I want?


title: Using a different projection to calculate areas using the Field Calculator

text: I have a dataset that is in a WGS84_UTM projected coordinate system, can I force the Field Calculator to use another coordinate system when I use it to calculate areas using pythons `!shape.area@squareyards!` ? My question is similar to checking "Use coordinate system of the data frame" in calculate geometry when the data frames projection is different; I would use this method but it does not have the option to calculate Square Yards. I could reproject the data to a new feature class, then calculate it and perform a join, or calculate the area in square feet using the data frames projection and dividing by 9, but since I will be doing this often the field calculator route would be the best way to go.


title: How to delete label style saved in ArcView?

text: The answer to my question may be very simple, but I can't seem to find it myself with online searches I have symbology and text/label styles (etc) that I've saved for easy access when I'm working on a Map. This is mostly for reasons of consistency. If I wanted to delete one, how would I go about accomplishing this? EXAMPLE (see image below): Highlighted in yellow, are duplicate label styles I saved at some point. I want to do some house cleaning. How do I delete one of them (for dummies please)? ![enter image description here](http://i.stack.imgur.com/SqMNt.jpg)


title: Convert .asc format into shapefile

text: I've got some data from "http://srtm.csi.cgiar.org/SELECTION/inputCoord.asp" and it is in *.asc form. To process them I need to convert them to shapefile (need to merge them all into one shapefile and add attribute table, on and on ...), but I don't know how. Some web sites such as "http://converter.mygeodata.eu/ " do convert my data but I'm talking about 2.5 GB of data and it is impossible for me to upload this much data ... how can I convert them inside QGIS or ArcGIS ??? Is there another software that can convert my *.asc file into *.shp ??? A lot of thanks


title: How to extract information from a buffer in ArcMap

text: I have a buffer in ArcMap 10 extending from a polyline, with numerous 5 km intervals (multiple buffer tool). I also have many layers already in raster format in ArcMap, including slope angle, geology, precipitation data etc. Basically I want to be able to get a value for each of these layers for each buffer polygon. The area is a transect through Nepal, and thus there are approximately 30 separate buffers extending 155 km towards the Himalayas. Ultimately I am aiming to produce a graph with distance up Nepal (x-axis) and then elevation (y-axis) and another variable such as precipitation (opposite y-axis).


title: Prevent OpenLayers ClickEvent from firing on Popup closeButton

text: My Code to add popups as selectControl to a Layer. Works fine!                     map.addLayer(window.markerLayer);          // SelectControl für Marker --> Popup     window.selectControl = new OpenLayers.Control.SelectFeature(             [ window.markerLayer ]);     map.addControl(window.selectControl);          window.selectControl.activate();               // code of: http://docs.openlayers.org/library/overlays.html#displaying-popups      function onPopupClose(evt) {               // 'this' is the popup     selectControl.unselect(this.feature);      } function onFeatureSelect(evt) {               feature = evt.feature;          var html_popup = "some html code";          popup = new OpenLayers.Popup.FramedCloud("featurePopup", feature.geometry             .getBounds().getCenterLonLat(), new OpenLayers.Size(250, 220),             html_popup, null, true, onPopupClose);          feature.popup = popup;     popup.feature = feature;     map.addPopup(popup);      } function onFeatureUnselect(evt) { feature = evt.feature; if (feature.popup) { popup.feature = null; map.removePopup(feature.popup); feature.popup.destroy(); feature.popup = null; } } ..But by clicking the closeButton on the Popup i am creating a new Marker because of this line:                   map.events.register("click", map, newMarkerClickEvent);      the Method: newMarkerClickEvent creates a new marker on the xy-Position of the event... How can i prevent from firing the click-Event to the map on clicking the Popup-Button ?? Thanks!


title: Creating legends on fly in Mapserver

text: I've got a spatial table in Oracle and a WMS-T in Mapserver, defined in a layer that query that table. The table has several columns and I pass the column to select as a parameter in OpenLayers, so the layer definition in mapserver is the following:               LAYER         NAME "layer_example"         CONNECTION "oracle/user@oracleservice/oracle"         CONNECTIONTYPE ORACLESPATIAL         DATA "geom FROM  (SELECT E.%type% as value, E.YEAR as YEAR, E.SHAPE_60M as geom FROM EMISSIONS E                           WHERE E.level=%level% ) USING NONE"              TYPE POLYGON              STATUS OFF         TRANSPARENCY 100              VALIDATION             'level' '^[0-9]$'             'type'       '^[a-zA-Z0-9]+$' #column name         END              PROJECTION                 "init=epsg:4326"         END         UNITS METERS         SIZEUNITS PIXELS         TOLERANCE 0         TOLERANCEUNITS PIXELS         CLASSITEM value         ### LEGEND DEFINITION HERE ###              METADATA                 "ows_title" "emissions"             "ows_srs" "EPSG:3035 EPSG:4326 EPSG:900913"             "ows_timeformat" "YYYY"             "ows_timeextent" "2000/2011"             "ows_timeitem" "YEAR"                 "ows_timedefault" "2000"             "wms_enable_request" "*"             "wms_metadataurl_format" "text/xml"             "wms_metadataurl_type" "ISO19115:2003"              END     END      The problem is the legend! Since the columns in the table has got different classes of values for their legends, how can I define the legend on fly? I tried to define legends in different files and then INCLUDE "legend_%type%.map", but the variable substitution doesn't work for INCLUDE. Actually, it would be great if Mapserver had this feature. Anyway, some ideas and suggestions? Hope the question is clear!


title: Converting map image to GIS

text: I’m very new to GIS, and I’ve been reading a lot of manuals and information, but I thought I’d ask a couple of questions here. One of my major hobbies is drawing maps of a completely imaginary country; it’s part of a massive project of mine. I’m getting to the point where I’m reaching practical limitations in the way I’m drawing the maps using traditional graphics software. However, I’ve already done quite a lot of work on these SVG maps in Inkscape, and I really would prefer not to start completely from scratch: the current file represents several months, if not a couple of years worth of work. They aren’t georeferenced, but they are drawn to scale (1:10000). Is there a method to convert or import the SVG into a GIS application, such as QGIS? I’ve been playing around with it all week to no avail. Ideally I want to just import the file, layers and all, and pick up drawing where I left off, but that seems to be all-but-impossible. Even just an easy way to add georeferencing to the exported image would be a good start. EDIT: Here’s one of the maps in question. It’s an export from Inkscape at 60 dpi; the original SVG file is ~ 3×3 metres(!) in size and tends to crash browsers. I’ve been documenting my progress with an album too. Although it’s unfinished, I figure if I’m planning to switch to QGIS, the earlier I do so the less re-mapping I’m probably going to have to do too. The bottom left corner is the point 0°E 40°N, and the map covers an area of roughly 1100 km² (33×33 km).


title: Openlayers - How to put street view icon on google layer?

text: I have a basic OL map with google layer. http://jsfiddle.net/Alophind/bDzg7/ I've read in similiar questions that the following css disable the google street view icon:               img {max-width: 100%;}      So I added (just for case) this css into my code:               img {max-width: none;}      And I still have no google street view icon. What else can I do ?


title: Slow tiles generation for leaflet js

text: I am currently trying to use leaflet to display images (jpeg or png) instead of maps, I have followed the instructions on this blog post to create the tiles using gdal2tiles.py here is the code we used               gdal2tile.py -p raster -z 0-5 -w none "SOURCE_PATH.png" "DESTINATION_FOLDER"      The command produces the necessary tiles for leaflet, but the problem is that it takes around 2 minutes to create tiles for the attached image (size: 777KB, resolution: 2592x1728) and a disastrous number of tiles of around 1600 file. I am trying to optimize the time and number of files to generate tiles for leaflet and looking for suggestions, please help.


title: Joining spatial and not spatial table in oracle

text: I have spatial table with rows of data representing polygon(A). I have another table that contains values for polygon field(B). Now I would like to join these two table and create a spatial table. I selected spatial column GEOLOC and data from another table. I want to make view out of this which would be mappable in mapinfo. But ![enter image description here](http://i.stack.imgur.com/7BiPI.png) **Making table mappable:** ![enter image description here](http://i.stack.imgur.com/OUMqp.png) **While making table mappable:** ![enter image description here](http://i.stack.imgur.com/tqnX6.png) I tried this with MS-SQL geometry colum and it works but fails with oracle spatial.


title: Placing a cost raster in Arcmap Table of Contents in 10.2.2

text: We have a set of routines that have been used for writing rasters to the ArcMap table of contents from our tools. The code that always worked is this:                   /// <summary>         /// Used for testing only.  Do not execute this in published server code.         /// </summary>         /// <param name="rasterDataset"></param>         /// <param name="layerName"></param>         internal static void PutToMap(IGeoDataset rasterDataset, string layerName)         {             if (rasterDataset == null)             {                 throw new ArgumentNullException("rasterDataset");             }             RasterLayer layer = new RasterLayer();             IApplication baseApp = (IApplication)Activator.CreateInstance(Type.GetTypeFromProgID("esriFramework.AppRef"));             IMxDocument doc = (IMxDocument)(baseApp.Document);                  IMap activeMap = doc.ActiveView.FocusMap;             if (rasterDataset is Raster)             {                 layer.CreateFromDataset(((IRaster2)(Raster)rasterDataset).RasterDataset);             }             else             {                 layer.CreateFromDataset((RasterDataset) rasterDataset);             }             layer.Name = layerName;             layer.Visible = true;             activeMap.AddLayer(layer);         }              /// <summary>         /// Used for testing only.  Do not execute this in published server code.         /// </summary>         /// <param name="layer">The raster layer or subclass that will be displayed</param>         /// <param name="layerName">The name to be put in the TOC</param>         internal static void PutToMap(IRasterLayer layer, string layerName)         {             IApplication baseApp = (IApplication)Activator.CreateInstance(Type.GetTypeFromProgID("esriFramework.AppRef"));             IMxDocument doc = (IMxDocument)(baseApp.Document);                  IMap activeMap = doc.ActiveView.FocusMap;                  layer.Name = layerName;             layer.Visible = true;             activeMap.AddLayer(layer);         }              internal static void PutToMap(IRaster raster, string layername)         {             IRasterLayer layer = new RasterLayerClass();             layer.CreateFromRaster(raster);             PutToMap(layer, layername);         }      We are now trying to debug a problem but when I try to write the cost raster out, these methods return an error when I try getting the singleton for ArcRef. The Activator.CreateInstance throws an exception in trying to load the module. Has anyone done this in 10.2.2 and is there a way to get things running that I'm missing??? I have tried changing the CreateInstance from Type.GetTypeFromProgID to Type.GetTypeFromCLSID and both approaches throws an exception. It returns an error of 8000ffff or something like that.


title: Can I open a form in layer action?

text: I want to open a form (made in QTDesign) in the layer action. It is possible?


title: FME: Identify a string within a field and populate another field with that same string

text: I have a shp which needs to be split according to attributes in a field. The problem is that the field contains other Characters and dates. For Example the Field is called "Landtype" and contains the following "Woodland_20120502". The Date changes for each data record so I can´t just select all records with this value. I need to first select all records that have the word "Woodland" and then populate another field with the word "Woodland" in FME or ArcGIS modelbuilder or Python. Could anyone advise me on how I would best do this?


title: Multiple Leaflet maps on website in different divs

text: I'm pretty new to leaflet.js and I'm trying to figure out how to assign the same map with the same set of options and layers to a different HTML container instead of having to remove and add a new one every time? I used to deal with Open layers 2 and I had map.render(div); option every time I wanted to set map to another div. Is there solution to this? Many thanks!


title: Import shapefile in CAD

text: I'm trying to import shapefiles into a DWG using diferent layers names and during my web search I found this post. I know you did this 2 years ago but my problem is very similar and I wonder if you could help me. Basically I've got a shapefiles of pipelines and I want to use the attribute name to sort the features into different CAD levels in teh same layer. I tried the import and although soem people say that "you can set the import to use data field" I can't see any option to do it and I'm really struggling. Any help appreciated. I'm running on AutoCAD Map 3D 2010. Thank you very much Patricia


title: OpenLayers vs Leaflet for specific app

text: Could anybody give me some advice on which one should I use? OpenLayers or Leaflet I will implement OSM on a web-app for taxi management. The requirements for the map part are:   1. Tracking of the drivers cars using a tablet through my web app.   2. Conversion of a given address to lat-long.   3. Showing a route between start and destination lat-long.   4. Completely free of charge Thanks in advance for your help. Lori


title: Can CartoDB provide Web Mapping Service (WMS)?

text: Does cartodb implements a wms service? If not, is there a way to contribute to the project to implement it?


title: How to get the name of a subtype field in a feature class in Python?

text: I have a geodatabase with many feature classes. Some of the feature classes have subtype fields. Their names vary feature class by feature class. I now want to find the names of those subtype fields. What I have so far is:               for fc in arcpy.ListFeatureClasses("*","",dataset):         subtypefields = arcpy.da.ListSubtypes(fc)         for subtype in subtypefields:      Now, I don't know how to continue from here. Any help will be appreciated. Thanks!


title: SLD style descriptor doesn't validate

text: I am trying to style a geotiff in OpenGeo Boundless...it doesn't validate... http://schemas.opengis.net/sld/1.0.0/StyledLayerDescriptor.xsd" version="1.0.0"> raster_layer PM 2.5 PM 2.5 Layers Color Scheme for PM 2.5 Rasters Feature


title: How to implement Google Earth style zoom in/out and pan animation in ArcGIS Engine applications?

text: I am developing an application using ArcGIS Engine 10 and the C# language. I wish to add some more animation capabilities into my application. In other words, I want to add animation to zoom in, zoom out and pan like the animation in Google Earth's pan and zoom in/out. I use OpenStreetMap (OSM) maps for the base map and also have other local layers in my application. Is there any way to use animations like Google Earth animation for zoom in/out and pan to act with my map ?


title: ArcGIS Query Layer performance documentation or feedback

text: I first tried asking this question in the ArcGIS forums but as is common I got a total of 0 responses. I'm willing to accept this may be due to a lack of knowledge on the subject, but I'm optimistic. I will cross-pollinate responses or just delete the ArcGIS post if it remains untouched. I've just started using Query Layers for read-only maps served through ArcGIS Server (10 SP2) reading data from Oracle 10g. I love it! I've always felt oddly uneasy about SDE and this seems a good sidestepping opportunity. In my tests performance looks good but I would like to see some documentation (ideally) or other user feedback on Query Layer performance on live services or large datasets and how it compares to other DBMS options (SDE, SDE DC). The only docs I've seen make a big deal out of the data being queried each time a user action occurs, but I would expect this to happen in all types of connections. Does SDE do some magical caching that a Query Layer doesn't? Does it better handle concurrent users in some way? If anything I would expect Query Layers to improve performance at runtime (lower overhead) and reduce development / setup time (fewer SDE configuration headaches). Any experience or informed thoughts on this would be very much appreciated. Please assume I am not concerned with versioning / multi-user editing, or any of SDE's other non performance-related benefits.


title: How to get the layer list back?

text: Sorry for the probably stupid question, but I closed a side panel I had, where I had a list of layers and which I used to switch on and off layers, and I cannot find how to get it back. I have tried all the buttons I could think of, but I cannot see it again. Can someone help? Cecilia


title: Open source methods for kriging?

text: I have a point dataset which I'd like to Krige, ideally using an open-source software package. If possible, I'd also like to choose the semi-variogram model during the process to improve the estimation.


title: GIS Shapefiles map to Zipcodes to districts

text: I have found GIS shape files from the following link: http://elections.il.gov/votinginformation/2011congressionaldistricts.aspx What would be the easiest way if possible to process/tool that I need so that we can map/lookup the list of zip+4 codes to the shapefiles to find the senator and representatives name, address, zip, district and county.


title: How to create callout/balloon labels in QGIS 2.0?

text: I'd like to use callout/balloon labels for a vector layer in my GIS using QGIS 2.0. The only way I can see to do this is to position an ordinary label where I want it and then draw an empty text annotation box behind the label. I've tried using an HTML annotation but I can't work our how to link to a feature attribute using HTML - I understand that I could use [% "attribute name" %] in the label (this results in a blank label), but how do I include the layer name in the html? I'd really appreciate some help. My apologies if my question isn't clear - I'd be more than happy to provide further details if necessary. Thanks


title: How to convert Ext.data.TreeStore loading from file to loading from .net webservice?

text: I have a Ext.data.TreeStore which read data from text file : ![enter image description here](http://i.stack.imgur.com/AdMc3.png)               this.store = Ext.create('Ext.data.TreeStore', {             model: "GeoExt.data.LayerTreeModel",             proxy: new Ext.data.HttpProxy({                 url: 'data/test.json',                 reader: {                     type: 'json',                     root: 'children',                     idProperty: 'Id'                 }             }),             folderSort: true         });      This is the file content:               [         {             text: "Category1",             leaf: false,             expanded: true,             checked: false,             children: [                 {                     text: "A1",                     layer: "A2",                     name: "A3",                     leaf: true,                     checked: false,                     nodeType: "gx_layer"                 },                 {                     text: "B1",                     layer: "B2",                     name: "B3",                     leaf: true,                     checked: false,                     nodeType: "gx_layer"                 }                 ],             nodeType: "gx_layer"         },         {             text: "Category2",             leaf: false,             expanded: true,             checked: false,             children: [                 {                     text: "C1",                     layer: "C2",                     name: "C3",                     leaf: true,                     checked: false,                     nodeType: "gx_layer"                 }],             nodeType: "gx_layer"         }     ]      I've noticed in firebug that the get doesn't return JSON struct , but plain text. but it works. ![enter image description here](http://i.stack.imgur.com/LU5sa.png) I want to get the same result but getting the tree nodes from webservice (which brings it back from postgres). The first simplest thing I tried was (also tried text/json):               [WebMethod]                     public void GetJSON()             {                 HttpContext.Current.Response.Clear();                 HttpContext.Current.Response.ContentType = "text/plain";                 HttpContext.Current.Response.Charset = "utf-8";                      string res = File.ReadAllText(@"D:\test.json");                 HttpContext.Current.Response.Write(res);             }           And I'm getting error :               Request format is unrecognized for URL unexpectedly ending in '/getJSON'.      which make sense because the file isn't JSON (but this is how I managed to make it work). I would like help with 2 things :   1. What is the correct JSON format I need to use to get the same result as I got in this example ?   2. How do I serve the same data from webservice ? Thank you !


title: PHP Mapscript Installation Problem on Mac (Mavericks)

text: After installing PHP and Mapserver and all the frameworks needed from http://www.kyngchaos.com, I can access Mapserver via cgi-bin, but I got the following error for PHP: 500 Internal Server Error. I've tried changing the permissions for the php.ini, httpd.conf and the folder: "Library/WebServer/Documents/" to 777 and this error remains the same. I am completely newbie in Mac.


title: How to create a shapefile with quantity of points?

text: I need a shapefile able to show the amount of houses in every 100m2 of a region. I don't want any density calculation, I need just the amount of point in each cell in a new shapefile. How can I do it? Thank you very much for your help!


title: How to overlay an NOAA nautical map as a layer in QGIS

text: I am new to QGIS and am doing some work for a friend. He is hoping to overlay an NOAA nautical map as a layer in QGIS. My only real experience is with ArcGIS. He is looking to do something similar to this: ![enter image description here](http://i.stack.imgur.com/wP8gr.png) Does anyone have any starting points for this? Any suggestions or directions in terms of how to do this, tutorials for this type of thing, etc. would be greatly appreciated. Here is the original link for the image: http://earthnc.com/online-nautical- charts


title: Calculating shape perimeter/area per zone in ArcGis 10

text: How can I calculate geometrical features of polygon (perimeter and area) per zone (similar to zonal stats for rasters)? I have a square fishnet (approx 20x35, 1 square kilometer each) and I need to calculate forest area and boundary length per each zone.


title: GeoServer/OpenLayers: Overlay WMS on Bing Map

text: I am coming from the ESRI side and learning GeoServer/OpenLayers. It seems the only way to overlay one layer on Bing or the like is to use WFS. In my opinion, WFS transfers a lot of data to the client side and the client has to render it. It's not that efficient when all I need is to have one layer on the top of the base map. allOverlays=true is too extreme for my case. Is there a way to overlay it as WMS? Similar to that in ArcGIS Server. Thanks Here is my code to try to overlay wms on Bing. The wms layer is not displaying. What's wrong?                       var map;              function init(){             var restrictedExtent = new OpenLayers.Bounds(-180, -90,                  180, 90);                  map = new OpenLayers.Map("map");             //map.allOverlays = true;                   map.addControl(new OpenLayers.Control.LayerSwitcher());                  var shaded = new OpenLayers.Layer.VirtualEarth("Shaded", {                 type: VEMapStyle.Shaded             });             var hybrid = new OpenLayers.Layer.VirtualEarth("Hybrid", {                 type: VEMapStyle.Hybrid             });             var aerial = new OpenLayers.Layer.VirtualEarth("Aerial", {                 type: VEMapStyle.Aerial             });                  var states = new OpenLayers.Layer.WMS(                 "States", "http://localhost:8080/geoserver/wms",                 {                     layers: 'state',                     styles: '',                     format: 'image/png',                     transparent: true                 }, {                     opacity: 0.5,                     singleTile: false,                      ratio: 1,                      isBaseLayer: true                 }              );                  map.addLayers([aerial, hybrid, shaded, states]);                  map.setCenter(new OpenLayers.LonLat(-110, 45), 3);         }


title: Is it possible to mimic ArcPublisher with PDF or similar?

text: I want to provide mapping output in a format that mimics the delivery provided by ArcPublisher and wonder if anybody has any experience in this. For instance, the PDF format now supports layers which can turned on and off, but can QGIS or GDAL export in such a way to enable this to happen rather than a basic PDF? Alternatively, perhaps there is another FOSS product which can provide a similar export. I want to avoid the client having to download and install software if at all possible (I am assuming that they will have a PDF viewer available). All suggestions welcome!


title: QGIS Composer Prints Blank Page

text: I am trying to print a map that I have made in QGIS print Composer. When I try to print the map or export it to a PDF (or any other format) as a raster all I get is a blank page. When I try to print not as a raster (unchecking the "print as raster option") it will print but obviosly not the way I want it to. eg. Layers arent transparent nor do other bits look right. I don't understand why this would happen. Why wont do I have a blank page when trying to print as a raster? I am running QGIS 2.0 on a UBUNTU 13.10/64bit/16GiB ram


title: IMapServerObjects.get_Map says object is null reference.

text: So I am trying to access a map layer, which I can connect to.               try                     {                                     agsconn.Connect();                         Log.Write("Connected.", "myclass");                         if (agsconn.IsConnected)                         {         ESRI.ArcGIS.Server.IServerObjectManager som = agsconn.ServerObjectManager;                                 IServerContext serverContext = som.CreateServerContext(mapDocName, "MapServer");                                 IMapServer mapServer = serverContext.ServerObject as IMapServer;                                 Log.Write("Server context established MapName=" + mapServer.DefaultMapName, "myclass");                                 ESRI.ArcGIS.Carto.IMapServerObjects mapSO = mapServer as ESRI.ArcGIS.Carto.IMapServerObjects;                                 Log.Write("IMapServerObjects created", "myclass");                                 IMap map = mapSO.get_Map(mapServer.DefaultMapName);                                 Log.Write("Layers map established", "myclass");                                 //int layerIndex = -1;                                 IFeatureLayer fLayer = map.get_Layer(0) as IFeatureLayer;      The map was created with direct connect, I don't know if that makes a difference. So far, the arcgis connects, and everything goes fine, but it stops here: **IMap map = mapSO.get_Map(mapServer.DefaultMapName);** I checked mapServer.DefaultMapName, and it returns "Layers" so that's all correct. mapSO also works because mapServer and mapSO are the same object basically. So how is there a null here? My log prints out an exception and never reaches "layers map established". The exception:               Object reference not set to an instance of an object.      I've tried a few other maps too and it doesn't work. But this code is straight from sample code I've seen on ESRI website, and it works for ESRI. I'm using ArcGIS 9.3.1.


title: Use ArcPy to rename MXD layers

text: I have a MXD with some layers that have parentheses in them. I'd like to rename them via script. Is this a valid approach?               import re     map = arcpy.mapping.MapDocument("CURRENT")     layers = arcpy.mapping.ListLayers(map)     for layer in layers:         if re.search("[()]", layer.name):             newName = re.sub("[()]", "", layer.name)             arcpy.Rename_management(layer, newName)      `Rename_management` is giving me an error: > Runtime error : ERROR 000840: The value is not a Data Element.


title: Use field names in combobox

text: I am currently stuck on a problem surrounding ArcGIS Explorer and comboboxes What I am trying to do is use 3 comboboxes to do the following:   1. View all features contained within ArcExplorer's table of Contents   2. View the selected layers containing field names    3. View unique attributes for the fields All of these would go in 3 seperate comboboxes. Would anyone have any suggestions on working through this? Thanks


title: Generalisation Strategies for Building Outlines

text: I have a dataset of building footprints vectorised from a raster dataset. The dataset is currently just a vector representation of the raster data. I would like a polygon dataset that only has vertices at the corners of the building. I have attempted to use the Simplify command in PostGIS with limited success. Is there a better way?


title: esriFramework.appRef is null

text: Trying to use the following code with an ArcGlobe Control               Type t = Type.GetTypeFromProgID("esriFramework.AppRef");     System.Object obj = Activator.CreateInstance(t);      Howevewr, t is coming back as null. Any ideas as to why? I do have the ESRI.ArcGIS.Framework reference in my workspace. I was hoping to access this so I could then use some ArcGlobe Commands via ICommandBar which I wanted to access via               IApplication myApp = obj as IApplication;     ICommandBars myCommandBars = myApplciation.Document.CommandBars;     myCommandItem = commandBars.Find("ID" , false);     myCommandItem.Execute();


title: ArcGIS: Create a line layer from a point layer and csv data file

text: **UPDATED DECEMBER 20** I have a Regular point layer and a csv file that has the same ids as the point layer. The csv file is organized as               ID     From    To     Weight     --     ----    --     ------     01     1234    1235   2     02     1234    1236   4     03     1234    1237   1     04     1235    1234   4     05     1236    1235   3          // etc      the point layer is organized in the following way:               MY_ID     -----     1234     1235     1236     1237      the From / To fields correspond to a my_ID in the point layer. I'd like to create a new line layer, such that the From, To, Weight correspond to both the point and csv files


title: ArcGlobe annotation

text: Hoping someone here can point me in the right direction. All i want to do is put some labels on a globe in ArcGlobe. I've spent time looking here : http://resources.arcgis.com and here http://help.arcgis.com/en/sdk/10.0/arcobjects_net/conceptualhelp/index.html Howver, I can't seem to find any solid examples. I've looked in the ArcGIS Snippet Finder when using Visual Studio as well. Also tried looking in the examples that include source code. Is there somewhere else to look - am I overlooking something really simple? Any help would be great.


title: Creating GeoTiff from Tiff+WorldFile+Proj

text: I'm currently using the global mapper export function to convert my tiffs+worldfile+projection file to geotiff. Is there a commandline tool (for Linux) that lets me do the same thing? "gdaltranslate" doesn't seem to use the projection file.


title: How to get Name of IRaster object using ArcObjects?

text: I have an IRaster object as output from a function. I need to be able to get the name of that raster. Is there a way that I can do that? I know that I can't mask an IRaster as an IRasterDataset. If I could, I would be able to get the Name property from that (I think). Thanks!


title: Routing Application in Windows

text: I am making an routing web application similar to Google Maps 'Get Direction' feature. But because of the limited control and data on GMaps, I decided to use OpenStreetMap instead. This also made me install PostGIS-PostgreSQl database, added routing functionality using pgRouting and now, I'm planning to use a Geoserver. All of these are done in Ubuntu Linux 9.04. My problem is I have difficulty using Linux and I'm more comfortable with Windows. As much as possible, I want to use Windows, but some of the things that I need are tested only in Linux (like osm2pgrouting, which I need to import OSM data to my database). Is there a way I can make this project in Windows? I mean a way I can implement the project using Windows, to avoid using Linux? Thank you in advance!


title: Geoserver Routing with OSM data and PostGIS

text: Does Geoserver have a routing functionality for the OSM data in PostGIS?


title: Has anyone successfully used Generate Near Table to create a table of distances between all polygons in a layer?

text: I'm trying to find the shortest edge-to-edge distance between polygons in a layer. An example is a table of the distances from county to county in a state. I'm attempting to use the generate near table tool. I set the proximity to 0 and made sure to de-select only the nearest neighbor. It runs for a few hrs, but when its finished, the ouput only returns the distance of the 1 nearest neighbor, which is 0 because they are touching each other. Any help would be much appreciated! -Kevin


title: Calculate distance in kilometers between two points

text: I have a few questions regarding the calculation of the distance between 2 postal codes. Take a look at the following code.               ###     # load the different packages we need     ###          require("rgdal")     require("maptools")     require("ggplot2")     require("plyr")     require("rgeos")          ###     # load the shape file     ###          belgium = readOGR(dsn="XXX", layer="NAME")     belgium.points = fortify(belgium, region="POSTCODE")      head(belgium.points)     tail(belgium.points)      This gives ;                     long      lat     id     1 151584.6 175068.7    1000     2 150829.0 174739.9    1000     3 150614.4 174580.1    1000     4 150178.8 174208.4    1000     5 149715.7 173823.0    1000     6 149126.5 173269.7    1000                long      lat    id     51839 83871.46 215487.5 9992     51840 84519.46 215217.5 9992     51841 83925.51 214823.3 9992     51842 83578.20 214704.9 9992     51843 82796.53 215022.7 9992     51844 81811.86 214716.5 9992      What I now want to do is calculating the distance between the centroid of each postal code. The postal code is given by the variable id in the dataset belgium.points. **Is the code below valid? Can I convert it to kilometers?** (Please if you need additional info just ask me.               lat = sapply(split(belgium.points, belgium.points$id), function(x) mean(x$lat))/1000     long = sapply(split(belgium.points, belgium.points$id), function(x) mean(x$long))/1000     lat.dist = outer(lat,lat,function(lat1,lat2) (lat1-lat2)^2)     long.dist = outer(long,long,function(long1,long2) (long1-long2)^2)     spatial.dist = lat.dist + long.dist


title: Trimble Net R9 Base Station Receiver File Structure

text: We just installed our new GPS Base Station Net R9 receiver that replaced our old 5700 receiver. The old 5700 used Trimble's GPS Base software, which the new Net R9 does not. My first problem that I can't seem to figure out is the file structure. With the old GPS Base software we were able to set the data logging system to write the files in folders that went Year -> Month -> Day -> hour. It was a nice clean structure. The new Net R9 logs all files in one folder where we have it set for every hour. This is going to be a massive folder after a while. Does anyone have any experience with the Net R9 file structure and how to get it to write like the old GPS Base software? Below are some pictures. Thanks in advance!! Trimble 5700 GPS Base File Structure Net R9 File Structure


title: How to adjust spacing and alignment between legend symbol and text in Qgis

text: I'm getting some weird and seemingly random spaces between the symbol and text labels in the legend for the Qgis print composer. How can I set the symbols and text to align left/right or divide the symbols into one centered column and the text in another? Also, in the legend currently, the third entry, "10m interim..." is not displaying the symbol at the proper size. It should be the largest yellow circle. ![Screenshot. Note different sized symbols and erratic spacings between symbols and labels.](http://i.stack.imgur.com/M1Vgb.jpg)


title: Is there an OSM High Level Dataset (pbf) available?

text: I'm currently trying to setup an offline Nominatim geocoder VM and I only need to have access high level data (countries, states and cities). Is there an OSM high level dataset (.pdf package) available?   A lighter version of Planet.pbf would be perfect for my needs (ref: geofabrik). Thanks


title: shapefile damaged after editing in ArcMap 10.0

text: I drew a layer on Autocad and exported it to a shapefile. I edited this layer on Arcmap. After editing the layer it is corrupted. I use Shapechk to fix it however I will need to edit plenty of layers and shapechk wouldn't be practical for me.i need an automated method


title: What kind of coordinate system is this?

text: I have a csv with a list of the world sea ports and their coordinates. The coordinates are in this format: 5234N 00226W How can i convert this to decimal Latitude / Longitude ? Data source


title: Easy way to Create Vector Grid polygons?

text: We're currently using MapInfo & MapBasic to analyze daily tests results for planning purposes.For each test,exported result may contain empty values (say 999,we know Mapinfo has no NULL values),so we decided to maintain each result (by point values) by averaging and cleaning empty results for each parameters inside a square(piece of grid),that way results may not repeat at the same coordinates and separation of empty values.(Or you can say get results by grid where each grid contain how many points and joining results for each cell by averaging) In this scenario,firstly we thought about creating these grids via **gridmaker.mbx** which is free tool for **Mapinfo** ,but when the area is too _large_ or the cell size is too _small_ or even if grid maker has no chance to imply our algoritm that like in the following We approach this problem by creating rectangles by 10 to 5000 meter values which are square cell sizes,and insert these values onto these tests via MBR or Convex polygon objects by checking out whether cell intersects these generated object from test. The problem starts here even if we process these cells,generating results **may last weeks,months even**.So my question is likely in the title,so i ask a different and easy way to approach this problem. Solutions,scenarios,algorithms you answer may not be solved only in mapinfo,but other tools as well. Best Regards


title: Calculating vector statistics in GRASS

text: I have a points layer, each of points have an year attribute and I want to count how many points in each year. I use GRASS-GIS and its v.univar function:               v.univar -e map=mod14_hmao column=year where="year=2007"      It is works for one year only, but I want get report like this:               2001 - 1000 points     2002 - 2000 points     ...      What function should I use?


title: How to fix distortions in rendered OpenStreetMap maps?

text: I use OpenStreetMap in order to generate maps and add some cool stuff in it. Actually I have good results but I seem to have a problem with reporting a good image ratio when I create my generated file. I put OpenStreetMap data in a PostgreSQL and I use PostGis to work with it in the EPSG:4326 projection. Everything is fine except this ratio problem. Actually I've tried almost every things (except the good one indeed)... that's why I'm looking for your help. I think I forgot something. In my program I use the min lon & min lat (as origin of the selection) ans max lon & max lat (as selection max coordinates)... basically I though a correct ratio would have been given by something like that (using N as size factor):               N = 10000     x_size = maxlon - minlon     y_size = maxlat - minlat     file_width = x_size * N     file_height = y_size * N      The most intriguing part is that this computation works for some parts of the world. In my examples below, the three cities use this ratio computation, but even if Montreal and Lyon don't work, Stockholm seems to work pretty well... The two others seem stretched. Lyon ( http://www.openstreetmap.org/?lat=45.7731&lon=4.8544&zoom=14&layers=M ): ![enter image description here](http://i.stack.imgur.com/6JOEn.png) Stockholm ( http://www.openstreetmap.org/?lat=59.3278&lon=18.0616&zoom=14&layers=M ): ![enter image description here](http://i.stack.imgur.com/Nw9KR.png) Montreal ( http://www.openstreetmap.org/?lat=45.5088&lon=-73.5878&zoom=13&layers=M ): ![enter image description here](http://i.stack.imgur.com/v9JqR.jpg) I think it's a noob problem... something silly enough to be barly invisible after hours of work on it. So I think you will find my problem easily. Thank you.


title: Limits to processing files from a list of rasters

text: I have written a script which selects all files that are within or subset by an area of interest polygon. This works fine but when the file list it's working has more than 150 records it stops comparing the extents after the 150th file. Does anyone know why? the log files and script are here in googledocs as they are very long. FileList = is a log of all 450 files which the program uses to iterate through. This is generated by FileList = arcpy.ListRasters() NoSubset – is a log of files which aren’t subset by the AOI. Note it finishes at 150 Success-Subset – is a log of all files that subset –note nothing beyond file Success-Within – is a log of all the files within the AOI. 20 files are logged when there should be 54. The ESRI log clearly shows that it stops processing beyond a certain point. I have included the tbx in the zip. Any suggestions for improving the coding would be appreciated as well. Note that I have turned the copy file and subset file functions off and am just using the generated csv's (link to files) rather than creating gigs of new files.


title: How to fill a hole in a multi-polygon layer on QGIS 2?

text: I need to fill a hole in a multi-polygon layer on QGIS 2.0. I've set the snapping option for the layer as follows :   * mode: to vertex and segment   * tolerance 4 pixels   * avoid intersection "enabled"   * topological editing "enabled"  When I try to edit a rough shape around the hole and confirm it right clicking, the system told me "The feature cannot be added because its geometry collapsed due to intersection avoidance" ![enter image description here](http://i.stack.imgur.com/sI5DC.jpg) * * * Update: I have a multipolygon file with a hole. I'm not able to fill the hole because if I manually edit it the system tells me that it is impossible to add it because the geometry is collapsed due to intersection avoidance. If I try to use the specifc tool it tells me that it's impossible due to the snapping settings. These are the images: ![specific tool](http://i.stack.imgur.com/MK4qA.jpg) ![error message about snapping setting](http://i.stack.imgur.com/mnTeU.jpg) ![snapping setting](http://i.stack.imgur.com/mffl8.jpg) ![enter image description here](http://i.stack.imgur.com/2GqVa.jpg) ![enter image description here](http://i.stack.imgur.com/8sZD5.jpg)


How do I export labels as editable text objects to be used in Illustrator? QGIS 1.8 Thanks.


title: Is there a API which return the geojson data for a countries borders?

text: Is there a simple free REST API where I can send the name or code of a country and get in return geojson data which highlights it's borders ? I'm writing an application where I will have a changing list of countries and I have to outline their surface/borders on a openlayers map. I've perfected/learned how to do it using sample geojson from this github but I haven't been able to find a RESTful api which serves me similar data? Does anyone know of such a service, bonus points if I can send multiple country names at once so I don't need to make a lot of calls.


title: why can I only define 6 custom CRS?

text: Using Q GIS 1.8.0 - I am trying to create several custom CRS definitions for the purpose of rectifying several split areas of a census shapefile. In the "Custom Coordinate Reference System Definition" window I seem to be limited to creating no more than 6 custom CRS. If I advance the arrows to 6 of 6 and enter a new definition it overwrites USER 100006: CRS. This is repeated on any of the other CRS definitions (1 of 6, USER 100000, etc...). Can I only create 6???


title: QGIS - Importing USGS DRGs

text: I have asked about this before. If I import a USGS DRG raster for Anchorage, Alaska using EPSG:26706 - NAD27 / UTM zone 6N all the UTM coordinates are correct. Now if I want to look at the coordinates in latitude/longitude and use ESPG:4267 NAD27 the image is rectified and stretched so the East/West looks twice as long as the North/South. The latitude/longitude of the coordinate points are correct but the map just looks funny. What am I doing wrong?


title: How can I merge a LOT of polylines automatically?

text: I have lots of MultiPolylines, three-dimensional arrays where the bottom level is an array of lat/lngs, the second level is an array of those points (a single polyline), and the top level is the group of lines. Like this:               [         [             [50.12312,0.2541],             [50.12322,0.2547],             ...         ],                  [52.412,1.142141],             [51.9812,0.21237],             ...         ],         ...     ]      The polylines in each group are generally contiguous, and I would like to merge them so that connecting lines within a group become a single line. How can I do this programmatically without some sort of desktop software? Is there Python,PHP,JavaScript,etc. code out there for this task or do I have to write it from scratch?


title: Run an ERDAS Imagine pmdl (or bcf) via Python?

text: How (or is it) possible to run an Erdas Imagine .pmdl model from a Python script? I have successfully run .eml models but cannot get the syntax right for .pmdl Here is what the line I want to run looks like in a BCF file:               modeler -nq $IMAGINE_HOME/etc/models/hpf_resolutionmerge2.pmdl -meter -state "$(ms)" "$(pan)" "$(pansharp)" 1 Unsigned_16_bit Integer Continuous useall 11 120 "None" 0.6 1 0.4 1:7 28 '$(mask.Ulx)' '$(mask.Uly)' '$(mask.Lrx)' '$(mask.Lry)' Map      ms, pan, pansharp, and mask are custom variables for raster files. Alternatively, can the BCF file be launched from outside Erdas?


title: What is the correct way to handle assembly resolution in .NET ArcGIS add-ins?

text: This has been a thorn in my side for some time now -- in .NET ArcGIS add-ins, assembly dependencies packaged with the add-in are often not resolved correctly. I am still not sure why and under which circumstances this happens but there are a couple of examples (1, 2) on this site, and I have run into several others involving reflection and WPF/XAML resources. Presumably this is happening because add-ins run in the ArcGIS application's `AppDomain`, whose assembly probing path is that of the folder containing the application executable, not the add-in. If you use `fuslogvw` to log assembly binding failures you'll see that it is probing e.g. `C:\Program Files\ArcGIS\Desktop10.0\Bin` for your assemblies, which will obviously never work with assemblies packaged in your add-in, which are extracted on app startup to e.g. `%USERPROFILE%\Local Settings\Application Data\ESRI\Desktop10.0\AssemblyCache\<{Add-In Guid}>`. Some workarounds I have seen include handling the `AppDomain.AssemblyResolve` event and explicitly forcing a referenced assembly to load by referencing a type defined within that assembly, e.g. `var dummy = typeof(SomeTypeInReferencedAssembly);`. The latter is what I have settled on doing since it is simple enough, but it seems tacky. The former seems fragile, especially when you have the potential for multiple add-ins needing to do the same thing. Is there not potential for them to mishandle the others' AssemblyResolve events, since they are all in the same AppDomain? Why did ESRI not provide an AppDomain for each add-in so they automatically get the correct probing path and avoid all this trouble? What practice does ESRI recommend in this situation?


title: How can I easily detect gaps in a polyline using arcmap 10.0 (arcview license)?

text: I'm making routes by selecting links from a streets layer which I export as a new temporarely layer in arcmap 10. This layer which I call route_temp unfortunately can be very long and contain gaps, for instance when i forget to select very small links which are difficult to see because they are overlaid by other layers. But sometimes the route contains gaps which have to be there because a cluster is in between it. What I would like to be able of is to visiualize all the gaps in the route, so that I can judge mysel which gaps have to be filled and which gaps have to stay in the route. Can anybody help me with this question, knows a tool for this, maybe a fieldcalculator-script or something else? Thank you in advance


title: When to use ModelBuilder over Python scripting and vice versa?

text: I am new to Python scripting in ArcGIS but not new to ModelBuilder. I would like to know the benefits of Python vs ModelBuilder. When should we write Python script for GIS automation rather than use ModelBuilder? What are the automation capabilities of Python that we cannot find in ModelBuilder?


title: Why can't I save my custom CRS?

text: I am unable to create custom CRSes in QGIS. The 'save' button doesn't seem to have any effect. My best guess is that I don't have permission to save in the folder or database where custom CRSes go. I have looked at the screenshots of other users who are having trouble with this, and I usually see "1 of 1" in the custom CRS browser. On my system, I see "* of 0" instead ![enter image description here](http://i.stack.imgur.com/VRfTb.png) Any assistance appreciated. QGIS version: 1.8.0 Ubuntu 12.04 LTS


title: run sql on Qgis

text: I want to do simple demo to our manager with quantum gis. i have buildings from my city with wgs84 projection on postgre sql database and also have non gis table which have who leaves this buildings (people name,surmane,birth date,etc...) my demo is ; select some building with selection tool from qgis and run some sql from non gis table and show how many people leaves in selected buildings? My sql is ready for work but how can i do this scenario with quantum gis? thanks


title: Alternatives for ERDAS Imagine

text: I am having trouble with the ArcGIS 10.1 and ERDAS Imagine 9.2 license conflict. From what I found out it isn't that rare, but interestingly enough why is there such a conflict in the first place? Therefore I am looking for some alternatives, heard about PCI's Geomatica but never tried it, any open source solutions?


title: Unexpected EOF while looking for matching

text: I was wondering if someone could help me with this error. I am new to writing shell script so here is my problem. When I run the following script it gives the error, line 31: unexpected EOF while looking for matching line 33: syntax error: unexpected end of file               g.remove rast=tmp     i=1     while [ $i !=366 ]     do          j=1951          r.mapcalc "tmptot = 0.0"          r.mapcalc "tmpsum = 0.0"          r.mapcalc "tmptot2 = 0.0"               while [$j !=2008 ]          do              g.copy rast=prec$j.$i,tmp              r.mapcalc "tmp2 = if(isnull(tmp),0.0,tmp) + tmptot"              r.mapcalc "tmptot = tm2"              r.mapcalc "tm3 = if(isnull(tmp),0.0,tmp*tmp) + tmptot2"              r.mapcalc "tmptot2 = tm3"              r.mapcalc "tmp4 = if(isnull(tmp),0.0,1) + tmpsum"              r.mapcalc "tmpsum = tmp4"              g.remove rast=tmp,tmp2,tmp3,tmp4              j=`expr $j + 1~          done              r.mapcalc "tmp5 = tmptot / tmpsum"         r.mapcalc "tmp6 = (tmptot * tmptot) / tmpsum"         r.mapcalc "tmpvar = sqrt((tmptot2 - tmp6) / (tmpsum - 1))"         g.rename rast=tmp5,av$i         g.rename rast=tmpvar,sd$i         g.rename rast=tmpsum,sum$i         g.remove rast=tmp6             i=`expr $i + 1`     done


title: Installing Demo Datasets

text: I installed QGIS without selecting the demo datasets. Now I realise that was a mistake. Is there a way to install these after QGIS installation or do I have to uninstall completely and re-install? Any advice appreciated.


title: how does arcgis server manage tile caches with an unlimited scale choice?

text: I am familiar with tile caching with geowebcache, but not with ArcGIS server. I imagine the basic concept is the same. First determine how many scales with be cached, set the image quality and format, and off you go. I would then usually design my own web client and set the scale choice to those set in the tile cache. Now I am using ArcGIS server and do not have the ability to set the client scales. So the client users can will basically hardly ever use the caches set to 1:1000, 1:2500, 1:5000, 1:10000,1:2500 etc...they will always be somewhere in between these set scales. Is this how it will work, or am I missing an important ESRI concept?


title: Successfully extracting values from a GRIB2 file...is it possible?

text: Let me say first that I'm starting to work with GRIB2 files, I've read the docs, and more or less know how they are structured. My goal is to be able to parse some GRIB2 files from NOAA, more precisely, those from the WaveWatch III model. So please forgive my ignorance if some of what I say sounds too naive :D I'm not attached to any technology, although I feel comfortable with Java and Python. To be more precise, you can find the GRIB2 file I'm using as a sample at this link. I'm also using the ToolsUI from Unidata for a more "user friendly" inspection of the GRIB2 file. As you can see, there are some variables in the provided GRIB2, but I'm only interested in 2 of them: Wind_direction_from_which_blowing_surface Direction_of_swell_waves_ordered_sequence_of_data So far, I've tried this code, which is based on this question also, and I'm able to open the GRIB2 file:               Grib2Input grib2Input;     String pathfile = "glw.grl.WDIR.grb2";          //Create RandomAccessFile     ucar.unidata.io.RandomAccessFile gribfile = null;         try {         gribfile = new ucar.unidata.io.RandomAccessFile(pathfile, "rw");     }      catch (IOException e) {         e.printStackTrace();     }          //Create grib2Input file     grib2Input = new Grib2Input(gribfile);     Grib2Data gd = new Grib2Data(gribfile);     try {         grib2Input.scan(false,  false);         List<Grib2Record> records = grib2Input.getRecords();         float[] data = null;         for (int i = 0; i < records.size(); i++){             //First get record information            System.out.println("record " + i +  "-> Param. Categ.: " +   records.get(i).getPDS().getPdsVars().getParameterCategory());         System.out.println("record " + i +  "-> Param. Number: " + records.get(i).getPDS().getPdsVars().getParameterNumber());         try {         data = gd.getData(records.get(i).getGdsOffset(),records.get(i).getPdsOffset(), records.get(i).getId().getRefTime());         }      catch (Exception e) {         e.printStackTrace();     }     //Second print some data     for (int j = 0; j < 50000; j+=5000){         System.out.println("data[" + j + "] -> " + data[j]);     }          System.out.println("--");     }      } catch (IOException e) { e.printStackTrace(); } Now the output, for the first 2 records, is this:               record 0-> Param. Categ.: 2     record 0-> Param. Number: 0     data[0] -> NaN     data[5000] -> 64.34     data[10000] -> NaN     data[15000] -> NaN     data[20000] -> 38.13     data[25000] -> NaN     data[30000] -> NaN     data[35000] -> NaN     data[40000] -> 66.55     data[45000] -> 107.06     --     record 1-> Param. Categ.: 0     record 1-> Param. Number: 7     data[0] -> NaN     data[5000] -> NaN     data[10000] -> NaN     data[15000] -> NaN     data[20000] -> NaN     data[25000] -> NaN     data[30000] -> NaN     data[35000] -> NaN     data[40000] -> NaN     data[45000] -> 328.5      Well, here are my question: 1.- Is there a better way to do this with NetCDF-Java library? It looks to me that for bigger size GRIBS this method is going to take ages. Besides, I should remove all NaN values from the data array for each record 2.- Should I go for more GRIB2 specific libraries like pygrib or pynio? If so...any working examples on parsing GRIB2 files with those? 3.- I also tried using wgrib2 command line tool...like...generating a CSV file with that tool and then parse the CSV file. Thing is that the generated CSV file is not very comprehensive either Thanks!!! Alejandro **EDIT: You can find an answer to this (using Python, tough)here. I must also say that I was able to do this too using Netcdf4Java, although with GRIB-API from ECMWF is way simpler.**


title: Is there an alternative for ERDAS IMAGINE AutoSync (eventuelly for QuantumGis)?

text: I had the possibility to work with a Demo-Version of the ERDAS IMAGINE module AutoSync. I have different geo-rasters (i.e. geoTIFFs) produced by UAV photography and which have to be matched the best way possible in order to carry out some analysis. If possible, even inner distortions should be corrected in oder to match the reference raster. So, the AutoSync module above does the trick, BUT as a small startup company, we are not able to pay the prices the company asks for the license. So I'm searching an alternative to this software. Does someone know a good alternative? Ideal would be open source, but even a proprietary software would be OK, if the price is accessible. Well, thanks for your answer in advance, Cheers Patribu


title: Geoserver WFS gives error but return features

text: I've noticed that after couple of requests geoserver becomes very slow. This is how I do the WFS request :               filter_buildings = new OpenLayers.Filter.FeatureId({ fids: featuresArr });         pins = new OpenLayers.Layer.Vector("Selected", {                     strategies: [new OpenLayers.Strategy.Fixed()],                     styleMap: highlightStyleMap,                     renderers: renderer,                     protocol :new OpenLayers.Protocol.WFS({                         version: "1.1.0",                         defaultFilter: filter_buildings,                         url: geoserver + "/wfs",                         featurePrefix: "ZoharHadar", //geoserver worspace name                         featureType: "Buildings", //geoserver Layer Name                                                 geometryName: "geom", // field in Feature Type details with type "Geometry"                         srsName: "EPSG:900913"                     })                 });      After inspecting the log I've noticed this :                   2872388 [btpool0-10] INFO  org.geoserver.wms  -      Request: getServiceInfo     2872397 [btpool0-10] INFO  org.geoserver.wms  -      Request: getFeatureInfo         PropertyNames = null         FeatureCount = 1         QueryLayers = [org.geoserver.wms.MapLayerInfo@6cc5ad04]         GetMapRequest =      GetMap Request      version: 1.1.1      output format: image/png      width height: 886,412      bbox: ReferencedEnvelope[3868350.809246 : 3868879.896801, 3767361.003943 : 3767607.035628]      layers: ZoharHadar:Buildings      styles: BuildingsPost         XPixel = 331         YPixel = 326         InfoFormat = application/vnd.ogc.gml         Exceptions = application/vnd.ogc.se_xml         Version = 1.1.1         Request = GetFeatureInfo         RawKvp = {INFO_FORMAT=application/vnd.ogc.gml, BBOX=3868350.809246,3767361.003943,3868879.896801,3767607.035628, QUERY_LAYERS=ZoharHadar:Buildings, SERVICE=WMS, HEIGHT=412, REQUEST=GetFeatureInfo, STYLES=, WIDTH=886, FEATURE_COUNT=1, VERSION=1.1.1, FORMAT=image/png, LAYERS=ZoharHadar:Buildings, Y=326, X=331, SRS=EPSG:900913}         Get = false         BaseUrl = http://localhost:8080/geoserver/         RequestCharset = null     2872651 [btpool0-0] INFO  org.geoserver.wfs  -      Request: getServiceInfo     2872670 [btpool0-0] WARN  org.geoserver.monitor  - Error handling request object     java.lang.ClassCastException: com.vividsolutions.jts.geom.Envelope cannot be cast to org.geotools.geometry.jts.ReferencedEnvelope         at org.geoserver.monitor.ows.wfs.WFSRequestObjectHandler.getBBox(WFSRequestObjectHandler.java:179)         at org.geoserver.monitor.ows.RequestObjectHandler.handle(RequestObjectHandler.java:50)         at org.geoserver.monitor.ows.MonitorCallback.operationDispatched(MonitorCallback.java:93)         at org.geoserver.ows.Dispatcher.fireOperationDispatchedCallback(Dispatcher.java:727)         at org.geoserver.ows.Dispatcher.dispatch(Dispatcher.java:722)         at org.geoserver.ows.Dispatcher.handleRequestInternal(Dispatcher.java:263)         at org.springframework.web.servlet.mvc.AbstractController.handleRequest(AbstractController.java:153)         at org.springframework.web.servlet.mvc.SimpleControllerHandlerAdapter.handle(SimpleControllerHandlerAdapter.java:48)         at org.springframework.web.servlet.DispatcherServlet.doDispatch(DispatcherServlet.java:923)         at org.springframework.web.servlet.DispatcherServlet.doService(DispatcherServlet.java:852)         at org.springframework.web.servlet.FrameworkServlet.processRequest(FrameworkServlet.java:882)         at org.springframework.web.servlet.FrameworkServlet.doPost(FrameworkServlet.java:789)         at javax.servlet.http.HttpServlet.service(HttpServlet.java:727)         at javax.servlet.http.HttpServlet.service(HttpServlet.java:820)         at org.mortbay.jetty.servlet.ServletHolder.handle(ServletHolder.java:487)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1093)         at org.geoserver.filters.ThreadLocalsCleanupFilter.doFilter(ThreadLocalsCleanupFilter.java:23)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1084)         at org.geoserver.filters.SpringDelegatingFilter$Chain.doFilter(SpringDelegatingFilter.java:74)         at org.geoserver.monitor.MonitorFilter.doFilter(MonitorFilter.java:137)         at org.geoserver.filters.SpringDelegatingFilter$Chain.doFilter(SpringDelegatingFilter.java:70)         at org.geoserver.filters.SpringDelegatingFilter.doFilter(SpringDelegatingFilter.java:45)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1084)         at org.geoserver.platform.AdvancedDispatchFilter.doFilter(AdvancedDispatchFilter.java:49)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1084)         at org.vfny.geoserver.filters.SetCharacterEncodingFilter.doFilter(SetCharacterEncodingFilter.java:109)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1084)         at org.springframework.security.web.FilterChainProxy$VirtualFilterChain.doFilter(FilterChainProxy.java:311)         at org.geoserver.security.filter.GeoServerCompositeFilter$NestedFilterChain.doFilter(GeoServerCompositeFilter.java:68)         at org.springframework.security.web.access.intercept.FilterSecurityInterceptor.invoke(FilterSecurityInterceptor.java:116)         at org.springframework.security.web.access.intercept.FilterSecurityInterceptor.doFilter(FilterSecurityInterceptor.java:83)         at org.geoserver.security.filter.GeoServerCompositeFilter$NestedFilterChain.doFilter(GeoServerCompositeFilter.java:72)         at org.geoserver.security.filter.GeoServerCompositeFilter.doFilter(GeoServerCompositeFilter.java:91)         at org.springframework.security.web.FilterChainProxy$VirtualFilterChain.doFilter(FilterChainProxy.java:323)         at org.geoserver.security.filter.GeoServerCompositeFilter$NestedFilterChain.doFilter(GeoServerCompositeFilter.java:68)         at org.springframework.security.web.access.ExceptionTranslationFilter.doFilter(ExceptionTranslationFilter.java:113)         at org.geoserver.security.filter.GeoServerCompositeFilter$NestedFilterChain.doFilter(GeoServerCompositeFilter.java:72)         at org.geoserver.security.filter.GeoServerCompositeFilter.doFilter(GeoServerCompositeFilter.java:91)         at org.springframework.security.web.FilterChainProxy$VirtualFilterChain.doFilter(FilterChainProxy.java:323)         at org.geoserver.security.filter.GeoServerAnonymousAuthenticationFilter.doFilter(GeoServerAnonymousAuthenticationFilter.java:53)         at org.springframework.security.web.FilterChainProxy$VirtualFilterChain.doFilter(FilterChainProxy.java:323)         at org.geoserver.security.filter.GeoServerCompositeFilter$NestedFilterChain.doFilter(GeoServerCompositeFilter.java:68)         at org.springframework.security.web.authentication.www.BasicAuthenticationFilter.doFilter(BasicAuthenticationFilter.java:150)         at org.geoserver.security.filter.GeoServerCompositeFilter$NestedFilterChain.doFilter(GeoServerCompositeFilter.java:72)         at org.geoserver.security.filter.GeoServerCompositeFilter.doFilter(GeoServerCompositeFilter.java:91)         at org.geoserver.security.filter.GeoServerBasicAuthenticationFilter.doFilter(GeoServerBasicAuthenticationFilter.java:82)         at org.springframework.security.web.FilterChainProxy$VirtualFilterChain.doFilter(FilterChainProxy.java:323)         at org.geoserver.security.filter.GeoServerCompositeFilter$NestedFilterChain.doFilter(GeoServerCompositeFilter.java:68)         at org.springframework.security.web.context.SecurityContextPersistenceFilter.doFilter(SecurityContextPersistenceFilter.java:87)         at org.geoserver.security.filter.GeoServerCompositeFilter$NestedFilterChain.doFilter(GeoServerCompositeFilter.java:72)         at org.geoserver.security.filter.GeoServerCompositeFilter.doFilter(GeoServerCompositeFilter.java:91)         at org.springframework.security.web.FilterChainProxy$VirtualFilterChain.doFilter(FilterChainProxy.java:323)         at org.springframework.security.web.FilterChainProxy.doFilter(FilterChainProxy.java:173)         at org.geoserver.security.GeoServerSecurityFilterChainProxy.doFilter(GeoServerSecurityFilterChainProxy.java:97)         at org.springframework.web.filter.DelegatingFilterProxy.invokeDelegate(DelegatingFilterProxy.java:346)         at org.springframework.web.filter.DelegatingFilterProxy.doFilter(DelegatingFilterProxy.java:259)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1084)         at org.geoserver.filters.LoggingFilter.doFilter(LoggingFilter.java:71)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1084)         at org.geoserver.filters.GZIPFilter.doFilter(GZIPFilter.java:41)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1084)         at org.geoserver.filters.SessionDebugFilter.doFilter(SessionDebugFilter.java:46)         at org.mortbay.jetty.servlet.ServletHandler$CachedChain.doFilter(ServletHandler.java:1084)         at org.mortbay.jetty.servlet.ServletHandler.handle(ServletHandler.java:360)         at org.mortbay.jetty.security.SecurityHandler.handle(SecurityHandler.java:216)         at org.mortbay.jetty.servlet.SessionHandler.handle(SessionHandler.java:181)         at org.mortbay.jetty.handler.ContextHandler.handle(ContextHandler.java:726)         at org.mortbay.jetty.webapp.WebAppContext.handle(WebAppContext.java:405)         at org.mortbay.jetty.handler.ContextHandlerCollection.handle(ContextHandlerCollection.java:206)         at org.mortbay.jetty.handler.HandlerCollection.handle(HandlerCollection.java:114)         at org.mortbay.jetty.handler.HandlerWrapper.handle(HandlerWrapper.java:152)         at org.mortbay.jetty.Server.handle(Server.java:324)         at org.mortbay.jetty.HttpConnection.handleRequest(HttpConnection.java:505)         at org.mortbay.jetty.HttpConnection$RequestHandler.content(HttpConnection.java:842)         at org.mortbay.jetty.HttpParser.parseNext(HttpParser.java:648)         at org.mortbay.jetty.HttpParser.parseAvailable(HttpParser.java:205)         at org.mortbay.jetty.HttpConnection.handle(HttpConnection.java:380)         at org.mortbay.io.nio.SelectChannelEndPoint.run(SelectChannelEndPoint.java:395)         at org.mortbay.thread.BoundedThreadPool$PoolThread.run(BoundedThreadPool.java:450)     2872679 [btpool0-0] INFO  org.geoserver.wfs  -      Request: getFeature         service = WFS         version = 1.1.0         baseUrl = http://localhost:8080/geoserver/         query[0]:             filter = [ Buildings.880254 ]             srsName = EPSG:900913             typeName[0] = {http://zoharhadar.opengeo.org}Buildings         outputFormat = text/xml; subtype=gml/3.1.1         resultType = results      What does the log means ?Am I doing something wrong?


title: Graduated colouring Voronoi polygons? (QGis)

text: I've got a huge database of depth points for the whole Irish sea which I'd like to be a nice smooth surface, or at least have a depth value associated with the entire area (i.e. include the inter-point spaces). I've tried interpolation but can't get it to work well, so I thought to make Voronoi polygons from the points, and then graduated colour them. Problem: In 2.1 master, a blank column is created in the Voronoi polygons shapefile, and switching from single symbol: ![single symbol](http://i.stack.imgur.com/q9zlz.png) towards graduated results in a complete blank: ![complete blank](http://i.stack.imgur.com/3xeEu.png). This even happens with the csv file, from dots: ![dots](http://i.stack.imgur.com/EUHe3.png) towards blank: ![blank](http://i.stack.imgur.com/livdG.png). Andre Joost got around a similar problem by using semicolon delimited files, but it feels to me like this is a genuine bug in master, especially given the same files don't blank the Voronoi shapefile in DuFour: ![don't blank the Voronoi shapefile in DuFour](http://i.stack.imgur.com/8Plyh.png). Has anyone else experienced this? I'm about to raise it on QGis hub. Question: can I colour the Voronoi polygons using the depth value? The Voronoi shapefile has depth as a field but it doesn't do anything - possibly it doesn't know how to deal with the fact that each Voronoi polygon has a number of depths equal to its number of edges (usually 4 in my case here). Can I colour them using the depth value from the points file somehow?


title: GetCapabilities for NED 1/9 sec and 1 foot orthoimagery?

text: I'm having trouble locating the WMS servers for these, specifically for San Francisco, which makes me think this level of data is only available tiled via email links from the USGS. Is that true? I've tried this WMS for NED without luck (unless I'm not recognizing the 1/9 sec DEM), and this WMS for the 1 foot ortho shows the 1 foot data set but not getting data back from it.


title: How to define Z dimension on table with PostGIS for import in TileMill?

text: I'm trying to reproduce the idea behind this blogpost from Mapbox, but with other data, for example the `ne_10m_populated_places_simple` dataset from Natural Earth. By adding a z-dimension to the data (eg. the `scalerank` column, size of city), the idea is to be able to set a height in `CartoCSS` in `Tilemil`l and get a 3D-representation of the data, with x,y as geo-point and z as value (like in Mapbox's example). The result could perhaps be something like this (except this map is oblique, different story): ![enter image description here](http://i.stack.imgur.com/ftFi6.png) What do I need to do with `PostGIS` (or other) so that the z-dimension is set for value and correctly converted so I can set the `height` in `CartoCSS`?


title: Moving GPS Points with a VB script - Bearing works well but distance is wrong

text: I am using a script I found online created in 2003 by one of the commentators here. It is a great script and exactly what I have been looking for and I'm very glad to have found it. It simply takes a points file and using two fields in the table, bearing and distance it moves the points based on these two fields. When I run it, the points always seem to move to the correct bearing but the distances are never correct. For example, I may need to move a point 30 meters away from the original position at a bearing of 85 degrees. The points move much more than that (hundreds of kilometers or more). Is this a projection issue? When I enter 30 in my table as my distance, do some projections interpret this as 30 degrees? How do I get ArcMap to recognize that this is 30 meters? If I can solve this I will be extremely happy. My points file covers very large sections of geography covering many states. Thank you The script is below               Option Explicit     Sub MovePoints()         ' move points in layer(0) based on bearing and distance         Dim pMxDoc As IMxDocument         Set pMxDoc = ThisDocument              Dim pFLayer As IFeatureLayer         Set pFLayer = pMxDoc.FocusMap.Layer(0)              Dim pFCur As IFeatureCursor         Set pFCur = pFLayer.FeatureClass.Search(Nothing, False)              Dim lBFld As Long, lDFld As Long         lBFld = pFCur.FindField("Bearing")         If lBFld = -1 Then             MsgBox "bearing field not found"             Exit Sub         End If              lDFld = pFCur.FindField("Distance")         If lDFld = -1 Then             MsgBox "distance field not found"             Exit Sub         End If              Dim pFeat As IFeature         Set pFeat = pFCur.NextFeature         Do Until pFeat Is Nothing             Dim pConstPoint As IConstructPoint             Set pConstPoint = New Point             Dim dAngle As Double, dDist As Double, PI As Double             PI = Atn(1#) * 4             dAngle = Azi2Cart(pFeat.Value(lBFld)) * PI / 180             dDist = pFeat.Value(lDFld)             pConstPoint.ConstructAngleDistance pFeat.ShapeCopy, _                                                dAngle, _                                                dDist             ' maybe set distance to 0 here             ' to avoid moving it again             'pFeat.Value(lDFld) = 0             Set pFeat.Shape = pConstPoint             pFeat.Store             Set pFeat = pFCur.NextFeature         Loop         Dim pAV As IActiveView         Set pAV = pMxDoc.FocusMap         pAV.Refresh     End Sub          Function Azi2Cart(ByVal dCompassAngle As Double) As Double         Azi2Cart = -dCompassAngle + 90     End Function


title: How do I interface with bundles using ESRI's Javascript API?

text: I've found this example... http://resources.esri.com/help/9.3/arcgisserver/apis/javascript/arcgis/help/jssamples_start.htm#jssamples/layers_custom_tiled.html Which uses an exploded cache. It notes at the beginning of the article... "You would most likely use ArcGISTileMapServiceLayer when accessing a service from 9.3 and beyond" Digging into the ArcGISTiledMapServiceLayer object it seems like it only connects to services via the rest url. I need to accomplish what is achieved in the Portland example (the link above) but with cache bundles instead of the exploded cache. **Potentially Useful Background** : I need to make an entirely self contained web map app with the javascript api that would use bundles that are stored with the app (as opposed to out on a server with arcgis installed) to render the map.


title: VB code for automatic refresh in arcmap with timer?

text: We have a shapefile of electrial network.The SCADA data of the electrical network is stored in oracle 10g express edition.We have to link oracle table with the arcmap and automatically refresh the map with the values which will be updated in oracle table.Do suggest a method and a code to execute this in Arcmap


title: How to Manually add a Z variable in Qgis onto a vector layer

text: I am trying to manually add one constant Z-variable (e.g. 200meters) across an entire Vector Layer in Qgis so that I may eventually use the Interpolation plugin to generate a Tin. I have not been able to find any information on how to do this on the website. Is this possible. Any information would be really helpful. Thanks! Anina


title: Is it possible to create Runtime content progamatically with ArcObjects

text: I'm exploring different options for creating an offline mobile application with Esri's technologies. The option that fits my needs the best is to use the ArcGIS Runtime SDK with offline maps created with ArcGIS for Desktop(In ArcMap : File->ShareAs->Runtime Content). The next step is now to automate the Runtime Content creation. Is there a way, using ArcObjects, to create Runtime content progamatically?


title: Querying Arcgis file geodatabase from sql

text: I need to compare data in a sql db with data stored in an arcgis file geodatabase. I'm not interested in the spatial attributes, only the added text attributes. What's the easiest way to do this? I can get the data exported to csv. But I can see this exercise will be repeated. And it would really make my life easier if I could also update data to clean it up. Any suggestions are greatly appreciated. Edit: I believe we run arcgis 10.1 or 10.2 I don't have arcgis desktop myself. Not sure if there are license limitations but installing it wouldn't be a problem. I didn't want to influence the answer but I've been reading about gdal/ogr. And wonder if that's an option. Also I have no idea how hard it would be to install arcsde. And whether that could read/convert a fgdb on the fly.


title: How to symbolize a layer on multiple fields?

text: I'm working in ArcGIS 10. I have a layer that contains **8 attributes** and each one of these fields can be categorized into **3 sub categories**. I need to display them on a single map. * * * _**Update_** I tried to add the same layer multiple times (8) and then symbolized each category. when I display the results are good but symbols now **superimpose on one another** , and hence the categories cannot be differentiated. The following is an image that contains a sample of what I am working on. "Area" corresponds to the point referring the region and 1A,1B and 1C are the fields of each feature. Each one of them (1A,1B,1C...,1H) is further categorized into 3 categories as shown in the table. I have added same file multiple times to get all the variables displayed on a single map, which are being displayed as well but they tend to super impose over one another. I want all the symbols to be shown separately. ![enter image description here](http://i.stack.imgur.com/C6N1O.png) Can anyone provide any ideas on how to do this?


title: Point features disappear when adding new data

text: Recently while adding features to a point layer i noticed that once i had saved my edits the already existing points disappeared. Although not visable on the map all of them and their data were still in the attribute table. Any guesses on how to make all of the points viewable again? Thanks, L


title: Loading a raster into a PostGIS 2.0 database on Windows

text: I'm trying to figure out how to load a raster into a PostGIS2.0 database (I have asked previous questions on this topic here and here). I am trying to use the `raster2pgsql.exe` program provided with PostGIS2.0. After figuring out that the command prompt in windows needs to be run as an administrator (In Windows 7 to run the command line as an administrator, type `cmd` into the search bar and hit `ctrl` \+ `shift` \+ `enter`) to enable `raster2pgsql.exe` to function I have sort-of managed to load a raster into my database. I have a raster file called `ras_test.tif` that I temporarily placed in the `bin` folder of my `postgresql` installation. Using the following code I convert and load this raster:               C:\Program Files (x86)\PostgreSQL\9.1\bin>raster2pgsql -s 102003 ras_test.tif -t> elev.sql     Processing 1/1: ras_test.tif          C:\Program Files (x86)\PostgreSQL\9.1\bin>psql.exe -p 5434 -U postgres -d test2 -f elev.sql     BEGIN     psql:elev.sql:2: NOTICE:  CREATE TABLE will create implicit sequence "-t_rid_seq" for serial column "-t.rid"     psql:elev.sql:2: NOTICE:  CREATE TABLE / PRIMARY KEY will create implicit index "-t_pkey" for table "-t"     CREATE TABLE     INSERT 0 1     COMMIT      When I view this table in PostGIS it looks like this: ![enter image description here](http://i.stack.imgur.com/sCYLT.png) However, I haven't been able to view this in QGIS, and I am not sure if I loaded it correctly as there appears to be no data in this file. It seems like I have loaded the filename in as a raster, rather than the data content. _Have I made any obvious errors that are preventing me from loading a raster into my database?_ * * * The PostGIS documentation provides an example of how to load a raster, but I don't understand which arguments are optional, and am still unclear on what I should use if I want to use the default schema. For example, in the following example from the documentation:               raster2pgsql -s 4236 -I -C -M *.tif -F -t myschema.demelevation > elev.sql      psql -d gisdb -f elev.sql        1. Do I have to provide an SRID? `-s 4236`   2. Are the arguments `-I -C -M` all optional?    3. `-t` appears to be tile size; Do I need to specify this if I do not have a custom schema?   4. Can I just leave out `myschema.demelevation` ? * * * EDIT: I've included the outcome of the suggestions below:               C:\Program Files (x86)\PostgreSQL\9.1\bin>raster2pgsql -s 102003 -t 300x300 ras_test.tif ras_test | psql.exe -U postgres      -d raster_analysis -h localhost -p 5434     Processing 1/1: ras_test.tif     BEGIN     NOTICE:  CREATE TABLE will create implicit sequence "ras_test_rid_seq" for serial column "ras_test.rid"     NOTICE:  CREATE TABLE / PRIMARY KEY will create implicit index "ras_test_pkey" for table "ras_test"     CREATE TABLE     INSERT 0 1     INSERT 0 1     INSERT 0 1     INSERT 0 1     COMMIT      This results in a table with two columns `rid` and `rast`. `rid` has four values, nad `rast` has none. When I try using some more arguments:               C:\Program Files (x86)\PostgreSQL\9.1\bin>raster2pgsql -I -C -e -Y -F -s 102003 -t 300x300 ras_test.tif ras_test1 | psql     .exe -U postgres -d raster_analysis -h localhost -p 5434     Processing 1/1: ras_test.tif     NOTICE:  CREATE TABLE will create implicit sequence "ras_test1_rid_seq" for serial column "ras_test1.rid"     NOTICE:  CREATE TABLE / PRIMARY KEY will create implicit index "ras_test1_pkey" for table "ras_test1"     CREATE TABLE     CREATE INDEX     ANALYZE     NOTICE:  Adding SRID constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding scale-X constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding scale-Y constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding blocksize-X constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding blocksize-Y constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding alignment constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding number of bands constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding pixel type constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding nodata value constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Unable to add constraint "enforce_nodata_values_rast"     CONTEXT:  PL/pgSQL function "_add_raster_constraint_nodata_values" line 40 at RETURN     PL/pgSQL function "addrasterconstraints" line 94 at assignment     PL/pgSQL function "addrasterconstraints" line 49 at RETURN     WARNING:  Unable to add constraint: 'nodata_values'.  Skipping     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN     NOTICE:  Adding maximum extent constraint     CONTEXT:  PL/pgSQL function "addrasterconstraints" line 49 at RETURN      addrasterconstraints     ----------------------      t     (1 row)      I get the following output. This results in a new table with the following structure: ![enter image description here](http://i.stack.imgur.com/0XaGK.jpg) I assume that this is not a properly loaded raster as I cannot view the data. Are there any other options that I can try? EDIT: This last attempt did work, I just was not accessing the raster correctly.


title: Why ArcGIS ONLINE OVER ArcGIS Server?

text: I do not understand why ESRI has two products serving similar purpose. Why one should have ArcGIS Server and why one should choose ArcGIS Online? Thanks.


title: ValueError when reading a raster band as a numpy array

text: Here's the problem shown in Python shell:               In [32]: from osgeo import gdal          In [33]: geotiff = '/Users/calvin/work/luc/proj_public/media/GLO     /Users/calvin/work/luc/proj_public/media/GLOBCOVER_L4_200901_200912_V2.3.color.tif     /Users/calvin/work/luc/proj_public/media/GLOBCOVER_L4_200901_200912_V2.3.color_orig.tif     /Users/calvin/work/luc/proj_public/media/GLOBCOVER_L4_200901_200912_V2.3.tif     /Users/calvin/work/luc/proj_public/media/GLOBCOVER_L4_200901_200912_V2.3_CLA_QL.tif          In [33]: geotiff = '/Users/calvin/work/luc/proj_public/media/GLOBCOVER_L4_200901_200912_V2.3.color.tif'          In [34]: src = gdal.Open(geotiff)          In [35]: band = src.GetRasterBand(1)          In [36]: band     Out[36]: <osgeo.gdal.Band; proxy of <Swig Object of type 'GDALRasterBandShadow *' at 0x10632b8f0> >          In [37]: ar = band.ReadAsArray()     ---------------------------------------------------------------------------     ValueError                                Traceback (most recent call last)     /Users/calvin/work/luc/<ipython-input-37-d7dde170464f> in <module>()     ----> 1 ar = band.ReadAsArray()          /opt/local/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/osgeo/gdal.pyc in ReadAsArray(self, xoff, yoff, win_xsize, win_ysize, buf_xsize, buf_ysize, buf_obj)         893         return gdalnumeric.BandReadAsArray( self, xoff, yoff,         894                                             win_xsize, win_ysize,     --> 895                                             buf_xsize, buf_ysize, buf_obj )         896          897     def WriteArray(self, array, xoff=0, yoff=0):          /opt/local/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/osgeo/gdal_array.pyc in BandReadAsArray(band, xoff, yoff, win_xsize, win_ysize, buf_xsize, buf_ysize, buf_obj)         225                                     buf_xsize, buf_ysize, datatype )         226         ar = numpy.fromstring(band_str,dtype=typecode)     --> 227         ar = numpy.reshape(ar, [buf_ysize,buf_xsize])         228          229         return ar          /opt/local/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/fromnumeric.pyc in reshape(a, newshape, order)         169     except AttributeError:         170         return _wrapit(a, 'reshape', newshape, order=order)     --> 171     return reshape(newshape, order=order)         172          173           ValueError: total size of new array must be unchanged


title: GDAL local files VS ftp files

text: I can successfully get raster information when the raster file is located on my local file system:               import gdal     from gdalconst import *     filename= "MOD11A1.A2012193.h11v11.005.2012196013534.hdf"     dataset = gdal.Open( filename, GA_ReadOnly )     print ( dataset.GetDriver().ShortName,dataset.GetDriver().LongName)      Output:               ('HDF4', 'Hierarchical Data Format Release 4')      But, when I try this with the same raster file, which this time is located on a ftp server, I get an error:               import gdal        from gdalconst import *     filename="ftp://ladsweb.nascom.nasa.gov/allData/5/MOD11A1/2012/193/MOD11A1.A2012193.h11v11.005.2012196013534.hdf"     dataset = gdal.Open( filename, GA_ReadOnly )     print ( dataset.GetDriver().ShortName,dataset.GetDriver().LongName)      Output with error:               ERROR 4: `/vsimem/http_1/MOD11A1.A2012193.h11v11.005.2012196013534.hdf' not recognised as a supported file format.          ('HDF4', 'Hierarchical Data Format Release 4')      Why is this?


title: postGIS (?) - ST_Difference preserving attributes

text: ## EDITED QUESTION: PROPER SYNTAX FOR ST_DIFFERENCE BETWEEN TABLES? I have two huge polygon datasets from shapefiles, call them 'baselevel' coverage and 'deluxe' coverage. They each contain over 2.4 million polygons, with many polygons describing complex coastal features. The baselevel dataset has complete coverage over the area of interest but less detailed attributes, while the deluxe dataset has gaps and less extensive spatial coverage but more detailed attributes. I want to create a third dataset which uses all the deluxe polygons, and uses the baselevel coverage to fill the gaps and extend the coverage. I've been approaching this as two steps: produce a 'fill' dataset which is just the base-level coverage minus all the deluxe polygons, then union the 'baselevel- fill' dataset with the 'deluxe' dataset. I can do this fairly efficiently in grass via v.overlay (create the fill layer with operator=NOT, then union with deluxe using operator=AND), after dicing the shapefiles into much smaller tiles. ### Final Update I'm pasting my final syntax below, though I'd still welcome any optimization suggestions ;-) based on Paul Ramseys answer, as well as some other stackexchange answers and postgis tutorials GRASS is pretty excellent, and could handle my polygons where a lot of other python, ogr, and qgis approaches. However the fastest I could complete a similar operation on my (very small) test dataset was 154 seconds. Code below completed in under 28 seconds :-)               CREATE TABLE deluxe_plus_baselvl AS     WITH deluxe_cutter AS (      SELECT ST_Union(d.geom) AS geom, b.gid      FROM baselevel b JOIN deluxe d ON ST_Intersects(b.geom, d.geom)      GROUP BY b.gid      ),     baselevel_cut AS (      SELECT ST_Difference(b.geom, d.geom) AS geom, b.gid      FROM baselevel b JOIN deluxe_cutter d ON b.gid = d.gid      )     SELECT 'baselevel' AS type, b.geom, b.gid     FROM baselevel_cut b     UNION ALL     SELECT 'deluxe' AS type, d.geom, d.gid     FROM deluxe d     UNION ALL     SELECT 'baselevel' AS type, b.geom, b.gid     FROM baselevel b LEFT JOIN deluxe d ON     ST_Intersects(b.geom,d.geom)     where d.gid is null;      Thank you for your help!


title: WFS connect SSL handshake fails

text: I'm running geoserver on my own server at home. Everything works great usually but when I try to connect to WFS from QGIS over https I get a warning saying SSL Handshake Failed. This is probably because I have a self-signed certificate. Do I have to pay for an authorized certificate or is there a way around this. P.S. https works fine in the web browser I just have to click to accept the certificate. I really don't want to pay $70 a year for a certificate just to use my GIS server for my own personal experimental use. Thank you.


title: How to find closest point without crossing obstacles?

text: I have a set of points and I'd like to be able to find the nearest point from another layer. The problem is that I need to find the closest point without crossing an obstacle (such as a river, highway, etc). Any idea how to do this?


title: How to restore gdalinfo metadata after ImageMagick processing?

text: I first print out the file info :               $gdalinfo crop.origin.tif      result:               + gdalinfo crop.origin.tif     Driver: GTiff/GeoTIFF     Files: crop.origin.tif     Size is 1920, 1950     Coordinate System is `'     Origin = (66.991666666666674,37.508333333333354)     Pixel Size = (0.016666666666667,-0.016666666666667)     Metadata:       NC_GLOBAL#Conventions=COARDS/CF-1.0       NC_GLOBAL#GMT_version=4.4.0       NC_GLOBAL#history=grdreformat ETOPO1_Ice_g_gdal.grd ETOPO1_Ice_g_gmt4.grd=ni       NC_GLOBAL#node_offset=0       NC_GLOBAL#title=ETOPO1_Ice_g_gmt4.grd       x#actual_range=-180, 180       x#long_name=Longitude       x#units=degrees       y#actual_range=-90, 90       y#long_name=Latitude       y#units=degrees       z#_FillValue=-2147483648       z#actual_range=-10898, 8271       z#long_name=z     Image Structure Metadata:       INTERLEAVE=BAND     Corner Coordinates:     Upper Left  (  66.9916667,  37.5083333)      Lower Left  (  66.9916667,   5.0083333)      Upper Right (  98.9916667,  37.5083333)      Lower Right (  98.9916667,   5.0083333)      Center      (  82.9916667,  21.2583333)      Band 1 Block=1920x2 Type=Int16, ColorInterp=Gray       NoData Value=-2147483648       Metadata:         NETCDF_VARNAME=z      **Given I then process my raster tif via ImageMagick** , here, for a resizing :               convert crop.origin.tif     -resize $(WIDTH) crop.small.tif      the file become a non-georeferenced *.tif image, losing its geotif metadata. So, **I reinject georeferencing & altitudes** :               gdal_translate -a_ullr 44 44 33 33 crop.small.tif crop.tif     gdalinfo -mm crop.tif      but the result is now different :               $gdalinfo crop.tif      result:               Driver: GTiff/GeoTIFF     Files: crop.tif     Size is 1280, 1300     Coordinate System is `'     Origin = (67.000000000000000,37.500000000000000)     Pixel Size = (0.025000000000000,-0.025000000000000)     Metadata:       TIFFTAG_DOCUMENTNAME=crop.small.tif       TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)       TIFFTAG_XRESOLUTION=72       TIFFTAG_YRESOLUTION=72     Image Structure Metadata:       INTERLEAVE=BAND     Corner Coordinates:     Upper Left  (  67.0000000,  37.5000000)      Lower Left  (  67.0000000,   5.0000000)      Upper Right (  99.0000000,  37.5000000)      Lower Right (  99.0000000,   5.0000000)      Center      (  83.0000000,  21.2500000)      Band 1 Block=1280x3 Type=UInt16, ColorInterp=Gray      There is a huge lost of data, especially into the `Metadata` section. And most of them could be recalculated. ![enter image description here](http://i.stack.imgur.com/ZPsse.png) So, **how to restore full metadata upon an ImageMagick processed file ?** * * * Please help me to create the `ImageMagick` tag.


title: Popups not showing on layer/feature/marker groups

text: I had two problems. One was that, once I got lots of markers in one layer, and lots of polylines in another, someof my markers, starting from the north, started to vanish.   I noted in the doco that the zindex is derived initially from the latitude,so I bumped the offset up to +1000 for my markers, and -1000 for my polylines, and things look ok. But I need the popups. This suggests they actually dont work https://github.com/Leaflet/Leaflet.markercluster/issues/180 So I've had to dump leaflets groups, and do my own layer management. Further more, the popups seem to vanish when I use optional baselayers. is this right ?


title: Avoid adding new generated intermediate layers into ArcMap during Python script running?

text: I am running my script to do geoprocessing in the python window of ArcMap. Since the geoprocessing includes lots of intermediate files, and these files were added in ArcMap at the same time of the processing. I thought it is a waste of time to display all these intermediate files in ArcMap. Can anyone help to suggest a way to only operate on these intermediate files in background and without adding them to ArcMap?