title: Why does the QGIS Server WMS layer not load in QGIS Desktop?

text: I have created a .QGS file with one PostGIS layer which is the alaska.shp file that comes from the samples. It was loaded via the shp2pgsql utility. When looking at my qgs project referencing the postgis layer, I see the alaska shape, everything appears good. Then I copy my .QGS file to my cgi-bin directory (along with the qgis_mapserve.fcgi file) to set up my WMS. Then I create a new project in Quantum GIS (1.7.2) and add WMS layer. I connect to my WMS and select the layer I want to import...so far so good. But after I add the WMS layer to my new QGIS project, I see nothing. The layer is checked "on" to render, and "render" is checked in the bottom right of Quantum GIS. Is there some way to tell if I have something set up wrong?


title: How to assign a grid value to points?

text: I'm sure there is an easy way to do what I am trying to do but I cant seem figure it out. All that i would like to do it create a grid (say 10mx10m) with unique values (i.e. grid 1,2,3,4) and then overlay a point file and assign a grid value to each point. All points within the first grid have field grid ID = 1 all points within the second grid have field grid ID = 2, etc. I have created grids with Xtools Pro> Create Fishnet which have unique FID but haven't figured out how to get that associated with the point data. Thanks, Corey


title: Spatial indexing of a point pattern in R

text: I have a point pattern and a tessellation (grid) in R. Is there a way to spatially index each of my points? In other words, my tessellation grid has _n_ quadrats, I would like to know for each of my points which of these _(1:n)_ quadrats they are associated with. If not, is there another way to do so in R? Many thanks for your suggestions!


title: raster standard deviation stretch python script

text: I am looking for ways of applying a permanent standard deviation stretch to rasters via python script. Are there any arcpy raster stretch functions? I was looking for a reference but could not find any. Are there any other python libraries that can perform the stretch? Any help much appreciated! Many thanks, Meryloo


title: PostGIS does not export shortest path to Shapefile

text: I am trying to draw a bunch of lines overlaying the globe map in a WGS84 projection. To this effect I have created a PostGIS DB with a table that has one column containing linestring data. My problem is that when I export it to WGS84 shapefile (either using pgsql2shp or adding a PostGIS layer directly in QGIS) the lines around the edges are drawn the wrong way. _For example:_ Suppose you have a line that connects Tokyo to Los Angeles. I would expect it to go over the pacific, but it does not. Instead I see a line go over Asia, Europe, Atlantic, US and finally end at LA. What I think is happening is that it's trying to keep the line contiguous in the WGS84 projection (which cuts the map in the middle of the Pacific Ocean), instead of picking the shortest route. How can I configure either PostGIS or one of the tools to do the right thing? _The right thing being:_ A line segment goes from Tokyo to the **right edge** of the map, and another line segment picks up at the **left edge** of the map (same coordinate) and goes to LA.


title: Is it best to use multipolygon or polygon features for large intersect operations?

text: I have two large layers that I want to do intersect operations. One is municipalities and the other is landuse. Right now, the landuse layer is of MultiPolygon type that is, for the entire state, I have only one entry in the table with the class "Forest", but composed of several small polygons. If I wanted to find all the forest areas in my city, what would be more efficient in PostGIS:   1. Do an intersect between my municipality and the FOREST multipolygon that takes up the entire state or,   2. Convert the Forest multipolygon into a new layer with several polygons and then do the Intersect My gut feeling is that if I break the multipolygon feature into single polygons, the intersect might be faster because it would only occur in the polygons that have a bounding box that overlaps my municipality. If I have multipolygon features, the intersect would have to go through the entire state to identify only the parts that fall inside the municipality. This is more or less what I'm trying yo get. I want the light green polygons. Does it matter it the forest polygons (F) are stored as multipolygons or singlepolygons? I'm going to be processing a very large dataset (Amazon) ![enter image description here](http://i.stack.imgur.com/RNctu.png)


title: QGIS developer help: custom icon

text: I'm a new QGIS developer and I have some difficulties to implement the follow feature: I want to track a live gps receiver and draw a custom icon (like png) in a map layer. I'm using a rubber band to move the gps point, and I found some plugins, but I need codify this feature for my application. perhaps have another way? Thank you.


title: Question for MOD35 cloud mask

text: I have a question about processing MOD35 cloud mask image. I worked reprojection on MOD35 cloud mask image using MODIS Swath Tool. The processed results of MOD35 cloud mask are consist of 6 band images (band0, band1, band2, band3, band4, band5). And each band images is 255 bits. I know that MOD35 cloud mask images is 48bit having a lot of information about clouds(1-2 : unobstructed FOV donfidence flag, 5bit : snow/ Ice background Flag, 9bit : thin cirrus dected ). what are the meaning of the band 0 to band 6 of processed result of MOD35 ? Is it right that create 48 bit image using all bands (band 0 to band 6)? I want to extract information from MOD35 cloud mask about (1-2 : unobstructed FOV donfidence flag, 5bit : snow/ Ice background Flag, 9bit : thin cirrus dected). How should handle this image? Ask an answer. Thank you :D


title: Install PostGIS and TIGER Data in Ubuntu 12.04

text: Could someone write a brief albeit dumbed down idiot's guide to installing postgis and loading national Tiger data on ubuntu? I've tried a few guides, namely, http://wiki.bitnami.com/@api/deki/pages/302/pdf, but I'm not having much luck. I apologize for the open ended nature of this question.


title: SpatialLite BulkLoad and Merge

text: Via the U.S. Census I have downloaded all of the census group shapefiles. These are approximately 19-20GB across 56 shapefiles. I want to merge these to do some testing. Clearly, shape file format is not going to work, so I am using spatialite. I have read that spatialite journalling causes performance issues when trying to load lots of data. Two questions:   1. Can I turn off journaling and just bulk load the data. I'll build stats and indices later.   2. Can I merge multiple shapefiles into one table within spatialite? I know I could use OGR and write a script to pump the data into a table one shape file at a time, but is this possible direct via spatial lite. Disclaimer:   1. I know that this will be an ~8,200,000 row table. I am not running this in production or planning on running many queries on this. This is a research project.   2. I could also go with a full blown PostgreSQL install and have this running, but I want to dig a bit deeper into Spatialite.


title: INSERT INTO Error

text: Hey guys I am trying to get the results of a query into a column in my table. The error I received is this ERROR: Quota exceeded by 16KB. I am using cartodb API,               INSERT INTO address (distance)     Select ST_Distance_Sphere(c.the_geom,(SELECT the_geom FROM lod_poly ORDER BY the_geom <->      c.the_geom LIMIT 1)) From address c      Is this a correct query for what I am trying to do? One other question I have is how can I convert meters into feet within the query ?


title: Calculate distance of downstream trace (geometric network) using ModelBuilder iterator

text: I'm using ArcGIS 10.1 and need to calculate the distance of down stream traces for each of several points in ModelBuilder. Using the Trace Geometric Network too, I have a feature class of points for the "flags" and a feature class of points for the "barriers". I want the trace to stop at the first barrier.   Since each point needs its own distance calculation, I'm imaging I need to use the iterator tool to accomplish this task. I need some guidance on how to go about this. Do I use the row iterator? Where does the iterator go? My end product should be a record for each flag showing the total downstream trace distance.   Another challenge I have is that I need to link the distance calculations for each point back to the flag, since the results of the trace calculations are independent of the flags (attribute wise), how can this link be made?   Any input would be helpful. Thank you! Edit: The attached graphic depicts what I've started out with. This is my first foray into Modelbuilder and the graphic is just a very rough start. ![ModelBuilder](http://i.stack.imgur.com/udVqX.png)


title: real-time weather forecast raster for Europe

text: I am trying to find a numeric weather forecast map (air pressure, temperature, snowfall, precipitation..) that I could display in my web map client. For the client I'm currently using OpenLayers with the 'web mercator' projection. I searched for "numeric weather forecast models" and found the GFS (global forecast system) and ECMWF. There are also some websites like: http://www.weatheronline.co.uk/cgi- bin/expertcharts?LANG=en&MENU=0000000000&CONT=euro&MODELL=gfs&MODELLTYP=1&BASE=-&VAR=prec&HH=3&ZOOM=0&ARCHIV=0&RES=0&WMO=&PERIOD= Here the weather forecast raster map is available as a .gif image and no info about the projection is provided, so I can't use it directly in the OpenLayers map client. I also found some weather forecast WMS services for USA, but nothing for Europe. Do you know about any free data source where the weather forecast map (ideally GFS or ECMWF model) available in a more GIS-compatible format? I would prefer the WMS but I also could use other formats (NetCDF, GeoTiff, ASCII grid..) and re-publish them with MapServer or GeoServer as a WMS. **Edit:** When I ask the question "where to find a WMS for my topic" It seems a good idea to do a search at http://spatineo.com. I typed in some keywords (GFS, ECMWF, weather forecast) and I found some interesting WMS services with weather forecast data. Spatineo.com also shows if the service is up and if the data access is free or restricted.


title: ArcGIS Server 10.0 Java - Upgrading JAI

text: I am developing a SOE for ArcGIS Server 10 SP 5 in Java. The SOE makes use of some GeoTools functionality - specifically the HeatmapProcess. I have written a small api (with all GeoTools, no ArcObjects) and packaged it as a standalone jar. The standalone jar works perfect in projects outside of ArcGIS Server. However, when used inside the SOE, I am eventually get hit with a "Provider com.sun.media.jai.imageioimpl.ImageReadWriteSpi not found" error message. I am traditionally a .NET guy, so Java's class loading mechanisms are a little alien to me. I am fairly confident I have the class-path of my SOE set properly - as all other GeoTools classes instantiate just fine. Looking at the ESRI install, I can see that some version of JAI is included in the java/lib/shared directory - but it is an older version than what I have included in the class-path of my SOE. Again admitting my lack of knowledge about all things Java to include the JVM, I am assuming that the version of JAI shipped with ArcGIS Server is taking precedence over the version I have included in the class-path of my SOE - through how ArcGIS Server starts the JVM. So I have a couple of questions for you Java guru's out there:   1. Is there some way I can update the version of JAI loaded by ArcGIS Server?   2. Is there some way I can make ArcGIS Server start the JVM giving precedence to the class-path I define?


Can ogr2ogr be used to append multiple DBF files? It works with SHP files, but I can't see how to do it with DBF files on their own. It ought to be possible, considering DBF is a subset of SHP.


title: select column in csv file in Python

text: I would like to study a range of data using the same technique in ArcMAP (it will take too much time to do it one by one ;) !). I've been just recently playing with ArcMAP and Python (which can help me to generate an automatic sequence). I also read that it is easier to import .csv file into ArcMAP. I'm importing the csv file using csv.reader (just coordinate of points without header, two colums x and y) like:               file = csv.read(...)      Then I didn't managed to print the file, because i obtained this message:               <_csv.writer object at 0x1426F7B0>      And I want to put the first column in the x coordinate and the second column in the y coordinate in a point shapefile? I read a lot of information how to create file. But I'm stuck in taking the first and second column of my file. Can you give me some advices? TA


title: CRS Conversion problem

text: My question is quite simple, this is probably why I couldn't find a conclusive answer anywhere. I have two datasets representing the same area in Romania, but in different CRS, A is local, B is WGS 84. **dataset A with:** Stereo_70 Authority: Custom Projection: Double_Stereographic False_Easting: 500000,0 False_Northing: 500000,0 Longitude_Of_Origin: 25,0 Scale_Factor: 0,99975 Latitude_Of_Origin: 46,0 Linear Unit: Meter (1,0) Geographic Coordinate System: GCS_Dealul_Piscului_1970 Angular Unit: Degree (0,0174532925199433) Prime Meridian: Greenwich (0,0) Datum: D_Dealul_Piscului_1970 Spheroid: Krasovsky_1940 Semimajor Axis: 6378245,0 Semiminor Axis: 6356863,018773047 Inverse Flattening: 298,3 **and dataset B with:** WGS 84 WKID: 4326 Authority: EPSG Angular Unit: degree (0,0174532925199433) Prime Meridian: Greenwich (0,0) Datum: WGS_1984 Spheroid: WGS 84 Semimajor Axis: 6378137,0 Semiminor Axis: 6356752,314245179 Inverse Flattening: 298,257223563 In order to convert B to A, I imported A's CRS to B=> Datasets did not overlap, they are miles away from each other. Could anyone explain what exactly is the problem. Thanks in advance


title: Zip Code-Level Choropleth Using Leaflet.JS

text: Brand new to Leaflet.js, but I'm enjoying learning it. However, I have a question I've been working on that even the Google hasn't been able to answer just yet. Is it possible to create a choropleth using zip-code level data? I've obviously seen plenty of choropleth examples using state-level data and even a few using county-level data. If so, what would be a good way to go about it using Leaflet? Do I need to obtain coordinates for the shapes of each zip code I'm interested in visualizing? Any help would be greatly appreciated.


title: Grass 7 - where is r.regression.multi hiding?

text: I compiled Grass 7 on Ubuntu and all works fine - except, I cannot find r.regression.multi. I installed it, yet it is nowhere to be found in the drop- down menus; r.regression.line is there and works, however.


title: identify subgraphs within planet_osm_line table

text: I want to identify subgraphs within the planet_osm_line table of the osm2pgsql-schema. those subgraphs are classified by certain tags i want to check. My first attempt was to do this recursively with a sql-function like this:               CREATE FUNCTION private.find_subgraphs() RETURNS void AS $$     DECLARE         good_way RECORD;         id_rec RECORD;         next_subgraph_rec RECORD;         id BIGINT;         next_subgraph_id BIGINT;     BEGIN         SELECT osm_id FROM private.good_ways WHERE subgraph_id IS NULL LIMIT 1 INTO id_rec;         id := id_rec.osm_id;         WHILE (id IS NOT NULL) LOOP             SELECT nextval('private.subgraph_seq') AS subgraph_id INTO next_subgraph_rec;             next_subgraph_id := next_subgraph_rec.subgraph_id;                  WITH RECURSIVE waysForSubgraph(osm_id, previous_way, visited_ways, done) AS (                  SELECT osm_id, public.planet_osm_line.way, ARRAY[ROW(osm_id, public.planet_osm_line.way)], false                      FROM public.planet_osm_line                     WHERE osm_id = id                 UNION ALL                 SELECT                         public.planet_osm_line.osm_id,                         public.planet_osm_line.way,                         visited_ways || ROW(public.planet_osm_line.osm_id, way),                         ROW(public.planet_osm_line.osm_id, public.planet_osm_line.way) = ANY (visited_ways)                 FROM waysForSubgraph, public.planet_osm_line JOIN private.good_ways                     ON public.planet_osm_line.osm_id = private.good_ways.osm_id                 WHERE ST_Intersects(previous_way, public.planet_osm_line.way) AND NOT done             )             UPDATE private.good_ways                 SET subgraph_id = next_subgraph_id                  WHERE osm_id IN (SELECT osm_id FROM waysForSubgraph);                  SELECT osm_id FROM private.good_ways WHERE subgraph_id IS NULL LIMIT 1 INTO id_rec;             id := id_rec.osm_id;         END LOOP;         RAISE NOTICE 'subgraph building done';       END;     $$ LANGUAGE plpgsql;      I know, that the assignment of `id` and `next_subgraph_id` is not perfect, but works. **Now to the problem:** For a small number of ways (like 10-20) this seems to work just fine. But when i run this on a sample of about 300 ways it looks like an infinite loop. Does anybody have suggestions on where the actual problem is? Or does anybody have a better approach on finding subgraphs like the ones i am looking for? Thanks in advance! **edit (in response to Mapper):** For now i am working on the lines of Germany (a little bit filtered). But here i am looking for ways that match certain criteria (like "is wihtin a forest (polygon)") and now i want to build subgraphs out of them. Of course, i am well aware of this being ressource hungry. But in this scenario it would be ok to run for one or two hours. Can you maybe a little bit more precise about that "smal set of subgraphs (Contraction hierarchies)" please? Thank you for your answer! **edit alogirthm:** For each `good_way` that is not in another subgraph already i look recursively for intersecting `private.good_ways` and than update their `subgraph_id` in `private.good_ways` with a common id i retrieve from the sequence `private.subgraph_seq`. The table good_ways has the `osm_id` of the ways that intersect forrests and the like (that match criteria for good ways). The `way` of a `good_way` is retrieved via the join with `planet_oms_line`. I hope you see what i try here. The problem seems to be, that the recursion is very very slow when the subgraph contains a lot of lines. But at least it should finish at some point, but it does not seem that way.


title: How to programatically create point and use it in Select by Location query?

text: I have a section grid for a US state. I also have a quarter/quarter section grid which splits the section grid into 16 equal squares. I have written a Python script that finds the XY coordinates of the centroid of each quarter/quarter section, so 16 of them. I need to create a temporary point feature at that centroid and then run a Select by Location on that point to see which quarter/quarter section it falls within and then I will be adding an attribute to that quarter/quarter section the point falls in. How can I create this temporary point feature to be used in a Select by Location? Is there a way to create a single point in memory?


title: Arcpy SelectLayerByAttribute for loop with "ADD_TO_SELECTION"

text: I am a little bit confused with the "ADD_TO_SELECTION" possibility. I process a point layer. First, the idea was to process the points with the different GBAUPS attribute separately. And it worked like this:               periode_field = "GBAUPS"     valueList = []     rows = arcpy.SearchCursor(input_gws_layer)          for row in rows:       valueList.append(row.getValue(periode_field))     uniqueSet = set(valueList)     uniqueList = list(uniqueSet)     uniqueList.sort()     print uniqueList          for period in uniqueList:            period = int(period)       bauperiod = arcpy.Select_analysis(input_gws_layer, "input_gws_layer_%i" % period, "GBAUPS = %i" % period)            ...... and the processing ......      But I would like now to add continuously the points from the previous layers to the next one (first just one value of the attribute GBAUPS, then this one and the second one etc.) I've tried like this, but the result is exactly the same as in the previous case.               for period in uniqueList:            period = int(period)            # selection       arcpy.MakeFeatureLayer_management(input_gws_layer, "input_gws_layer_lyr")       arcpy.SelectLayerByAttribute_management("input_gws_layer_lyr", "ADD_TO_SELECTION", "GBAUPS = %i" % period)       arcpy.CopyFeatures_management("input_gws_layer_lyr", "input_gws_layer_lyr_%i" % period)       ...... and the processing ......      I am not sure if I do something wrong in the SelectLayerByAttribute_management or in the CopyFeatures_management. Any help is very much appreciated! thank you!


title: Trouble intersecting a polyline and contours (Hack index / stream length gradient index)

text: Working on creating polylines for measurement using SL index (hack index). So far I have the main trunk of a stream, and contours at 10m. What I am trying to do is a) get the elevation at each intersection of the stream (polyline) and a contour for the full length of the stream b)get the distance along the stream between each contour intersection So far I've been able use spatial join to get the contact between the stream and the contour but the distance is off (checked by using measure tool on the polyline of the full stream and summing them in the attribute table) Any suggestions on how to make this work? Thanks!


title: Dissolve/aggregate polygons with ogr2ogr or gpc?

text: I have a large number of polygons representing the boundaries of US counties which I need to merge into as few polygons as possible. Is it possible to do this with non-commercial software, or software of a reasonable price? PostGIS is not available, nor is Java, but C is. Thanks in anticipation Lee


title: Execute OpenLayers WPS request using WPS Client

text: I am having trouble executing WPS processes with the OpenLayers WPS Client. I can execute the processes in Geoserver, however when I try to execute the process from my OpenLayers test app I keep getting errors in my browsers javascript console relating to JSON and xmlhttprequest. My JS code for creating the client, geometry and executing the process is below - In the Javascript console I can see the WPS XML and it all looks well formed and good to go. Is the issue with my code trying to decipher the results? JS console image below. ![enter image description here](http://i.stack.imgur.com/QgLQ1.jpg)                   runoff = OpenLayers.Geometry.fromWKT('Point(148 -36)');          client = new OpenLayers.WPSClient({         servers: {             local: '/geoserver/wps'         }     });          client.execute({             server: 'local',             process: 'gs:PointBuffers',             inputs: {                 center: runoff,                 crs: 'EPSG:4326',                 distance: 1000             },             success: function(outputs) {                 for (var i=0, ii=outputs.result.length; i<ii; ++i) {                 alert(outputs.result[i].geometry.toString());             }             }         });      Thanks Andrew


title: How to edit PostgreSQL geodatabase table using ArcPy?

text: I'm trying to update a table using a python script and it's not working. The table resides in a postgres geodatabase and is versioned. I have also verified my user permissions (I'm the admin). The version is public (i.e. not protected). The table properties show that the current user has select, insert, update, and delete privileges. Here's my command               with arcpy.da.Editor(arcpy.env.workspace) as edit:          with arcpy.da.InsertCursor("service", ["service_id"]) as cur:             cur.insertRow([0])      The error:               Runtime error      Traceback (most recent call last):       File "<string>", line 3, in <module>     RuntimeError: The requested operation is invalid on a closed state [ugdb.sde.service]


title: How to batch rename feature datasets and then rename feature classes based on the new dataset name?

text: I need to create a python script (or model) that would first rename all the feature datasets in a gdb. The feature dataset is in a set format and I need to change the prefix of the dataset. For example: from "H2611111:xxx" to "gen26_xxx" \- where xxx is the string that I need to keep and prefix of "H2611111" is a constant that I need to replace with "gen26_". Once this is completed, I would then need to rename each feature class inside each renamed dataset. The feature classes needs to be prefixed with the newly renamed dataset separated by an underscore. For example: in dataset gen26_xxx: from "test1" to "gen26_xxx_test1" from "test2" to "gen26_xxx_test2" in dataset gen26_zzz: from "new1" to "gen26_zzz_new1" I am new to python (and model building), if anyone has any helpful insights, that would be greatly appreciated. I am using ArcGIS 10 and 9.3.1 Thank you


title: DrawPolygon problem in DrawDynamicLayer

text: I have a class that inherits from BaseDynamicLayer and I override the DrawDynamicLayer() method. Is there something special that needs to be done for a successful inDynamicDisplay.DrawPolygon(poly) call? The method crashes every time I call it with a "Exception from HRESULT: 0x80040218." I can use the same point collection to create a Polyline and call DrawPolyline(polyLine) without a problem. I can also successfully use all the other draw commands, like DrawMarker() & DrawRectangle().


title: How to use Select by Attributes query with semicolons?

text: I am trying to select multiple (50 or more) unique numeric records using ArcGIS 10's **"Select by Attributes"** using semicolons ';' after each record which it verifies as good to go. Please see **screen shot** of where it verifies my expression "successful". Yet when OK is hit it **only selects the first record.** Do we really have to enter **_"xyz" = 1 OR "xyz" = 2 OR "xyz" = 3_** etc. etc. to have multiple records selected?


title: Idaho's southern border looks curved in ArcMap. How do I make Idaho look flat?

text: I have a polygon of Idaho's borders in ArcMap. The southern border of the state is slightly curved. I realize that the Earth is curved in reality, but I want Idaho to look flat. Is there a way to do this in ArcMap? I feel like a change in projection or coordinate system might help, but I don't know what those actually do, nor how to do it. Any help is appreciated.


title: PHP tile server

text: I am newbie in PHP. I am using tileserver.php to serve my *.mbtiles as web map. I have configure it and when I am trying to access the localhost it replying like following ![enter image description here](http://i.stack.imgur.com/Lp1oc.jpg) Now please help me to visualize the available map in browser. What should I do to visualize the map in browser? I have tried the following codes with OpenLayers                               <title>OpenLayers Example</title>


title: Map Server 2013 GetFeature request failing without providing enough of a clue

text: This is a Map Server support question more then a GIS question... I am trying to buffer features using Map Server 2013 but receiving the following error:               **An exception occurred in FDO component. String does not represent a valid filter.**      But I cannot see anything wrong with the request below. The capabilities service boasts the following spatial capabilities (note _Intersects_ is listed):               <ogc:Spatial_Capabilities>       <ogc:GeometryOperands>         <ogc:GeometryOperand>gml:Envelope</ogc:GeometryOperand>         <ogc:GeometryOperand>gml:LineString</ogc:GeometryOperand>         <ogc:GeometryOperand>gml:Point</ogc:GeometryOperand>         <ogc:GeometryOperand>gml:Polygon</ogc:GeometryOperand>       </ogc:GeometryOperands>       <ogc:SpatialOperators>         <ogc:SpatialOperator name="BBOX" />         <ogc:SpatialOperator name="Contains" />         <ogc:SpatialOperator name="Crosses" />         <ogc:SpatialOperator name="Disjoint" />         <ogc:SpatialOperator name="DWithin" />         <ogc:SpatialOperator name="Equals" />         <ogc:SpatialOperator name="Intersects" />         <ogc:SpatialOperator name="Overlaps" />         <ogc:SpatialOperator name="Touches" />         <ogc:SpatialOperator name="Within" />       </ogc:SpatialOperators>     </ogc:Spatial_Capabilities>      This is the actual request which is failing:               <GetFeature xmlns:wfs="http://www.opengis.net/wfs" service="WFS"     version="1.1.0" maxFeatures="500" xmlns="http://www.opengis.net/wfs">       <wfs:Query typeName="ns63133509:ST_WASTE_CONTAINER" SRSNAME="EPSG:26910">         <Filter xmlns="http://www.opengis.net/ogc">           <Intersects>             <PropertyName>Geometry</PropertyName>             <Polygon srsName="EPSG:26910" xmlns="http://www.opengis.net/gml">               <exterior>                 <LinearRing>                   <posList>489272.945 5457256.531 489271.98426402 5457246.7764839                   489269.138976626 5457237.39682838 489264.518480615                   5457228.75248835 489258.300339059 5457221.17566094                   489250.723511651 5457214.95751939 489242.079171618                   5457210.33702337 489232.699516101 5457207.49173598 489222.945                   5457206.531 489213.190483899 5457207.49173598 489203.810828382                   5457210.33702337 489195.166488349 5457214.95751939                   489187.589660941 5457221.17566094 489181.371519385                   5457228.75248835 489176.751023374 5457237.39682838                   489173.90573598 5457246.7764839 489172.945 5457256.531                   489173.90573598 5457266.2855161 489176.751023374 5457275.66517162                   489181.371519385 5457284.30951165 489187.589660941                   5457291.88633906 489195.166488349 5457298.10448062                   489203.810828382 5457302.72497663 489213.190483899                   5457305.57026402 489222.945 5457306.531 489232.699516101                   5457305.57026402 489242.079171618 5457302.72497663                   489250.723511651 5457298.10448062 489258.300339059                   5457291.88633906 489264.518480615 5457284.30951165                   489269.138976626 5457275.66517162 489271.98426402 5457266.2855161                   489272.945 5457256.531</posList>                 </LinearRing>               </exterior>             </Polygon>           </Intersects>         </Filter>       </wfs:Query>     </GetFeature>      The DescribeFeatureType looks like this:               <?xml version="1.0" encoding="utf-8"?>     <xs:schema xmlns:xs="http://www.w3.org/2001/XMLSchema"     targetNamespace="http://fdo.osgeo.org/schemas/feature/ns63133509"     xmlns:fdo="http://fdo.osgeo.org/schemas" xmlns:gml="http://www.opengis.net/gml"     xmlns:ns63133509="http://fdo.osgeo.org/schemas/feature/ns63133509"     elementFormDefault="qualified" attributeFormDefault="unqualified">       <xs:element name="ST_WASTE_CONTAINER"       type="ns63133509:ST_WASTE_CONTAINERType" abstract="false"       substitutionGroup="gml:_Feature">         <xs:key name="ST_WASTE_CONTAINERKey">           <xs:selector xpath=".//ST_WASTE_CONTAINER" />           <xs:field xpath="PRIMARYINDEX" />         </xs:key>       </xs:element>       <xs:complexType name="ST_WASTE_CONTAINERType" abstract="false"       fdo:geometryName="Geometry">         <xs:complexContent>           <xs:extension base="gml:AbstractFeatureType">             <xs:sequence>               <xs:element name="PRIMARYINDEX">                 <xs:simpleType>                   <xs:restriction base="xs:int" />                 </xs:simpleType>               </xs:element>               <xs:element name="COMPKEY" minOccurs="0">                 <xs:simpleType>                   <xs:restriction base="xs:string">                     <xs:maxLength value="30" />                   </xs:restriction>                 </xs:simpleType>               </xs:element>               <xs:element name="FID" minOccurs="0">                 <xs:simpleType>                   <xs:restriction base="xs:double" />                 </xs:simpleType>               </xs:element>               <xs:element name="UNITID" minOccurs="0">                 <xs:simpleType>                   <xs:restriction base="xs:string">                     <xs:maxLength value="30" />                   </xs:restriction>                 </xs:simpleType>               </xs:element>               <xs:element name="Geometry" type="gml:AbstractGeometryType"               fdo:hasMeasure="false" fdo:hasElevation="false"               fdo:srsName="UTM83-10" fdo:geometricTypes="point"               fdo:geometryTypes="point multipoint" />             </xs:sequence>           </xs:extension>         </xs:complexContent>       </xs:complexType>     </xs:schema>      UPDATE: Error log isn't very helpful but does provide a stack trace:               <2012-11-26T14:50:00>   9400        10.42.100.93    WfsUser      Error: An exception occurred in FDO component.             String does not represent a valid filter.       StackTrace:       - MgFeatureServiceHandler.ProcessOperation() line 83 file d:\build\ims2013\build_52.8\ent\os\server\src\services\feature\FeatureServiceHandler.cpp       - MgOpGetWfsFeature.Execute() line 208 file d:\build\ims2013\build_52.8\ent\os\server\src\services\feature\OpGetWfsFeature.cpp       - MgServerFeatureService.GetWfsFeature() line 2017 file d:\build\ims2013\build_52.8\ent\os\server\src\services\feature\ServerFeatureService.cpp       - MgServerFeatureService.SelectFeatures() line 462 file d:\build\ims2013\build_52.8\ent\os\server\src\services\feature\ServerFeatureService.cpp       - MgServerSelectFeatures.SelectFeatures() line 306 file d:\build\ims2013\build_52.8\ent\os\server\src\services\feature\ServerSelectFeatures.cpp


title: How do I select all polygons in a layer that border 2+ polygons in the same layer?

text: My data consists of equal-sized polygons. I need to identify and select all the polygons that border two or more polygons.


title: How to create a polyline-based "heatmap" from GPS tracks?

text: This winter I am planning to track my downhill skiing/snowboarding using a GPS. Most of my riding will occur at the same resort. I would like to be able to create a sort of "heatmap" that shows the amount of runs that I have made in a given area. As I add more and more GPS traces to my database, my goal would be to see a sort of linear heatmap of the most traveled areas. Given the nature of downhill skiing, you would expect the uphill chairlift lines will be the "hottest" because it will be the only places visited over-and-over again. Given that 1) my track will not be the same every time and 2) the area covered by following one "run" may be a few hundred feet wide, what might be the best way to analyze this "linear" data to create a sort of heatmap? My thoughts were to buffer the lines, then intersect the polys to get a sort of Venn Diagram thing going. My preference is to use open-source technologies. I've got QGIS and PostGIS loaded and available. **UPDATE** : In regards to @blah238's response, I was thinking of something that might be able to "collect" the number of passes ("runs") through an area, and then symbolize by the count. Conceptually, this would be similar to ArcGIS "Collect Events" (but for lines, not points) or Collapse Dual Lines To Centerline (but for multiple line in roughly the same area). A more visual example of a similar concept might be a traffic-flow map, where highly-congested areas would equate to "highly-traveled" ski runs/areas: ![Google Traffic Map](http://i.stack.imgur.com/J27W8.jpg) I've read the following questions that may give some ideas, but they don't really address what I am trying to accomplish: Clustering Trajectories (GPS data of (x,y) points) and Mining the data Managing error with GPS routes (theoretical framework?)


title: How successful are current metadata implementations in achieving its objectives

text: Metadata infrastructure is well established as an essential practice within the GIS community. Quoting from the FGDC: > Metadata helps people who use geospatial data find the data they need and > determine how best to use it. Metadata supports producers in locating and > using their own data resources and data consumers in locating and using data > resources produced by others. In your experience how effectively do current standards actually help to achieve these objectives? I am particularly interested in analysis of ISO 19115:2003 and how it compares to over methods of capturing this kind of information.


title: Is there a GIS to locate exact opposite locations on earth?

text: I was wandering if there was a GIS that offers to locate the exact opposite location of key locations (countries, even cities ?) on earth ?


title: Is the ArcGIS 10 SP2 an improvement worth installing? How buggy is it?

text: I've downloaded the service pack but have not yet installed it. Today I stumbled over a post on the ESRI forums where problems with reference scale and labeling are discussed. Apparently, this is now a known bug (NIM067942). Normally, I would install a Service Pack to correct existing issues and to possibly improve performance but it seems that new bug(s) were introduced as well. This particular bug could potentially be a nuisance. How buggy is it? Will I benefit from installing it? What are your experiences with this particular service pack?


title: Raster overlay in openlayers

text: I am using openlayers and recently started to learn it from the 'Open Layers CookBook'. I have been able to get Google Maps and WMS layers as base layers and even a vector layer overlayed correctly on top of the base layers at the exact location. Lastly, I wanted to display my Interpolation Raster image which I generated in QGIS to be displayed over the vector layer at the precise location. **But the 'Interpolation Raster' seems to load as an image and does not come up at the correct location, the way my vector layer does without any trouble. Can anyone point me in the right direction ?**. I have been searching two books on layers and on google for the solution to this problem but it seems to be not talked about much. The code for my webpage using openlayers is as below:               var map;          function init() {       // Initiating 'map' Variable and Adding Controls       map = new OpenLayers.Map('map_element', {         controls: [           new OpenLayers.Control.Navigation(),           new OpenLayers.Control.PanZoomBar(),           new OpenLayers.Control.ZoomBox(),           new OpenLayers.Control.ZoomToMaxExtent(),           new OpenLayers.Control.LayerSwitcher({}),           new OpenLayers.Control.MousePosition({             prefix: 'Long: ',             separator: ' | Lat: '           }),           new OpenLayers.Control.ScaleLine({}),           new OpenLayers.Control.OverviewMap()         ]       });            // BaseMap Layers       var wms = new OpenLayers.Layer.WMS(         'OpenLayers WMS', 'http://vmap0.tiles.osgeo.org/wms/vmap0', {           layers: 'basic'         }, {           projection: 'EPSG:4269',           maxExtent: new OpenLayers.Bounds(-180, -90,             180, 90),           maxResolution: 1.40625,           units: 'm'         }       );       var google_hybrid = new OpenLayers.Layer.Google(         "Google Hybrid", {           type: google.maps.MapTypeId.HYBRID         }       );       var google_physical = new OpenLayers.Layer.Google(         "Google Physical", {           type: google.maps.MapTypeId.TERRAIN         }       );       var google_satellite = new OpenLayers.Layer.Google(         "Google Satellite", {           type: google.maps.MapTypeId.SATELLITE         }       );       var google_streets = new OpenLayers.Layer.Google(         'Google Streets', {           type: google.maps.MapTypeId.ROADMAP,           projection: new OpenLayers.Projection("EPSG:3857"),           units: "dd",           maxResolution: 1.40625,           displayProjection: new OpenLayers.Projection("EPSG:4326")         }, {           "sphericalMercator": true         }       );       // Vector Layers       var vector_layer = new OpenLayers.Layer.Vector('Basic Vector Layer')       // Load Vector Data       var vector_KML = new OpenLayers.Layer.Vector("KML", {         //projection: 'EPSG: 4269',         protocol: new OpenLayers.Protocol.HTTP({           url: 'New_Sector_Polygon_KML.kml',           format: new OpenLayers.Format.KML({             extractStyles: true,             extractAttributes: true,             maxDepth: 2           })         }),         strategies: [new OpenLayers.Strategy.Fixed()]       });            var vector_KML_new = new OpenLayers.Layer.Vector(         "New KML", {           //projection: 'EPSG: 4269',           protocol: new OpenLayers.Protocol.HTTP({             url: 'New_Sector_Polygon_KML_new.kml',             format: new OpenLayers.Format.KML({               extractStyles: true,               extractAttributes: true,               maxDepth: 2             })           }),           strategies: [new OpenLayers.Strategy.Fixed()]         }       )            // Loading Image Layer---Interpolation       var image_layer = new OpenLayers.Layer.Image(         "IDW_Interpolation_RASTER",         'Export_CLIP_RENDERED_Image.png',         new OpenLayers.Bounds(-5, -5, 5, 5), {           numZoomLevels: 7,           maxResolution: .625         }, {           isBaseLayer: false         }       );            // Adding map Controls       map.addControl(new OpenLayers.Control.EditingToolbar(vector_layer));            // Adding BaseMap and Vector Layers       map.addLayers([         wms,         vector_layer, vector_KML, vector_KML_new, image_layer       ]);       var map_control_panel = new OpenLayers.Control.Panel({         div: document.getElementById('panel_div'),         defaultControl: new OpenLayers.Control.Navigation()       });            map.setCenter(new OpenLayers.LonLat(0, 0), 2);                           // Create some empty vector layers       var pointLayer = new       OpenLayers.Layer.Vector("Points");       var lineLayer = new OpenLayers.Layer.Vector("Lines");       var polygonLayer = new       OpenLayers.Layer.Vector("Polygon");       // Add layers to the map       map.addLayers([polygonLayer, lineLayer, pointLayer]);                 // Fill layers       initializePointLayer();       initializeLineLayer();       initializePolygonLayer();       // Create some random points.       function initializePointLayer() {         var pointFeatures = [];         for (var i = 0; i < 50; i++) {           var px = Math.random() * 360 - 180;           var py = Math.random() * 180 - 90;           var pointGeometry = new           OpenLayers.Geometry.Point(px, py);           var pointFeature = new           OpenLayers.Feature.Vector(pointGeometry);           pointFeatures.push(pointFeature);         }         pointLayer.addFeatures(pointFeatures);       }                 // Create some random lines       function initializeLineLayer() {         for (var j = 0; j < 1; j++) {           var pointGeometries = [];           for (var i = 0; i < 10; i++) {             var px = Math.random() * 240 - 120;             var py = Math.random() * 100 - 50;             var pointGeometry = new             OpenLayers.Geometry.Point(px, py);             pointGeometries.push(pointGeometry);           }           var lineGeometry = new OpenLayers.Geometry.           LineString(pointGeometries);           var lineFeature = new           OpenLayers.Feature.Vector(lineGeometry);           lineLayer.addFeatures(lineFeature);         }       }                      // Create some random polygons       function initializePolygonLayer() {         for (var j = 0; j < 1; j++) {           var pointGeometries = [];           for (var i = 0; i < 5; i++) {             var px = Math.random() * 240 - 180;             var py = Math.random() * 100 - 90;             var pointGeometry = new             OpenLayers.Geometry.Point(px, py);             pointGeometries.push(pointGeometry);           }           var linearGeometry = new OpenLayers.Geometry.           LinearRing(pointGeometries);           var polygonGeometry = new OpenLayers.           Geometry.Polygon([linearGeometry]);           var polygonFeature = new OpenLayers.           Feature.Vector(polygonGeometry);           polygonLayer.addFeatures(polygonFeature);         }       }                      //Trying to Create a circle >>>>>---------------------------------------->>>>>            // CIRCLE       var circle = OpenLayers.Geometry.Polygon.createRegularPolygon(         new OpenLayers.Geometry.Point(30, 53),         50,         30       );            var attributes_c = {         name: "my circle",         bar: "foo"       };            var feature_c = new OpenLayers.Feature.Vector(circle);       var layer_c = new OpenLayers.Layer.Vector("Circle");       layer_c.addFeatures(feature_c);       map.addLayer(layer_c);            // ----xxxxxxx------------------------------------------------------------------------xxxxxxxxxxxxx                      if (!map.getCenter()) {         map.zoomToMaxExtent();       }     }


title: Arcobjects - Properties Grid

text: Does anybody know a way, in VBA, to change the grid interval? For example: if my layout´s scale is 1:2000, my grid will have 200 meters interval and my layout´s scale is 1:5000, my grid will have 500 meters interval. Thanks, Alex


title: Are there any online WKT editors?

text: I've been using QGIS and it's great. However, we can't use that because it has to be installed (please, don't ask). Instead, what we need to do are very _simple_ administrative boundary drawings over a (insert fav map here) via some website. Our requirements are:   1. click to create points which then form a polygon   2. drag existing points to edit a polygon Yes, I know Google Maps and others can do all of this. We're hoping to leverage an existing site that can already offer this, though. We would also love to 'load' some WKT data, so we can edit that (if required), where load means paste into a textbox. Anyone know of any sites out there with this functionality?


title: Using Google Map Spatial Reference to plot HeatMap?

text: Basically I am trying to plot a Google heatmap based on the JSON objects returned. Here is my JavaScript:               function getHeatMap() {     var gradient = [                             'rgba(185, 185, 203, 0)',                             'rgba(145, 145, 192, 0)',                             'rgba(65, 65, 207, 0)',                             'rgba(30, 30, 229, 1)',                             'rgba(0, 185, 255, 1)',                             'rgba(0, 255, 215, 1)',                             'rgba(0, 255, 15, 1)',                             'rgba(0, 255, 0, 1)',                             'rgba(255, 255, 0, 1)',                             'rgba(255, 235, 0, 1)',                             'rgba(255, 0, 0, 1)'     ]     $.ajax({         url: "index.aspx/getBusCommuter",         type: "POST",         data: "",         contentType: "application/json; charset=utf-8",         dataType: "json",         success: function (data) {             var parsed = JSON.parse(data.d);             var data = [];             var max = 0;             $.each(parsed, function (i, jsondata) {                 var coordXicon = jsondata.BusStopX;                 var coordYicon = jsondata.BusStopY;                 var commuterAmt = jsondata.CommuterAmt;                 var point = new esri.geometry.Point({ "x": coordXicon, "y": coordYicon, "spatialReference": { "wkid": 3414 } });                      var count = commuterAmt;                 data.push({                     location: point,                     weight: count                 });                                   if (count > max) {                     max = count;                 }                      var heatmap = new google.maps.visualization.HeatmapLayer({                     data: data,                     radius: 30,                     opacity: 0.8,                     maxIntensity: max,                     setMap: map                 });                      heatmap.set('gradient', gradient);             });         },         error: function (request, state, errors) {         }     });     }       Here is the part where I initialise the baseMap:               function setMap() {     function init() {         require(             [                 "esri/map",                 "dojo/dom-construct",                 "dojo/domReady!"             ],             function              (                 Map,                 domConstruct             ) {                 map = Map("map-canvas",                 {                     //infoWindow: popup                 });                 map.setZoom(1);                 coreFunctions();             });          }     dojo.ready(init);     }          function coreFunctions() {     try {         addLayersToMap();     }     catch (err) {         console.log(err);     }     }     function addLayersToData() {         var layer = new esri.layers.ArcGISTiledMapServiceLayer("https://www.onemap.sg/ArcGIS/rest/services/BASEMAP/MapServer");         mapLayers.push(layer);         var layer2 = new esri.layers.ArcGISTiledMapServiceLayer("http://www.onemap.sg/ArcGIS/rest/services/LOT_VIEW/MapServer");         mapLayers.push(layer2);         layer2.hide();     }          function addLayersToMap() {         addLayersToData();         for (var a = 0; a < mapLayers.length; a++) {             map.addLayer(mapLayers[a]);         }         map.setZoom(1);     }      The JSON objects did returned some data. But when I try to plot it as a heatmap overlayer, the heatmap does not show. I am using spatial reference 3414 for the data to plot on OneMap. I wonder is it because of the difference of spatial reference that caused my heatmap does not show. **EDIT**               function getHeatMap() {     heatLayer = new HeatmapLayer({         config: {             "useLocalMaximum": true,             "radius": 40,             "gradient": {                 0.45: "rgb(000,000,255)",                 0.55: "rgb(000,255,255)",                 0.65: "rgb(000,255,000)",                 0.95: "rgb(255,255,000)",                 1.00: "rgb(255,000,000)"             }         },         "map": map,         "domNodeId": "heatLayer",         "opacity": 0.85     });          $.ajax({         url: "index.aspx/getBusCommuter",         type: "POST",         data: "",         contentType: "application/json; charset=utf-8",         dataType: "json",         success: function (data) {             var parsed = JSON.parse(data.d);             $.each(parsed, function (i, jsondata) {                 var coordXicon = jsondata.BusStopX;                 var coordYicon = jsondata.BusStopY;                 var data = [            {                attributes: {},                geometry: {                    spatialReference: { wkid: 102100 },                    type: "point",                    x: coordXicon,                    y: coordYicon                }            }                 ];                      heatLayer.setData(data);                  });             map.addLayer(heatLayer);         },         error: function (request, state, errors) {         }     });      }


title: How to convert a shapefile into a WMS? Wie kann man eine Shape-Datei in ein WMS umwandeln?

text: can anyone tell me by chance, how to transform a shape file into a WMS and then load it into mapbender? kann mir zufällig jemand sagen, wie man eine Shape-Datei in ein WMS bekommt, um diesen dann ins Mapbender zu laden?


title: ArcMap 10.1 How to make mass (in bulk) update by attribute value

text: I have feature class with 70,000+ records (file geodatabase). I want to make in bulk update on an attribute where the value is NULL. Select by attribute shows that I have 5,000 records like that. Doing it one by one is going to be a pain. In sql pseudo code it should be something like that:                UPDATE <MyAttributeName> = 'NewValue'      WHERE <MyAttributeName> IS NULL      How do I do that with ArcMap 10.1?


title: Feature Layer Web Mercator Service is too slow in IE and zooming does not work in FF

text: I am facing issue of Feature Layer being too slow and hangs browser , we changed the projection system to WebMercator performance increased in FF but in IE it sill hangs.when i zoom in FF it hangs what is the issue we have about 220 polygons which is not too much .Below is services and html code Rest Service               Layer: World_Countires_v3_WMS (ID: 0)          Name: World_Countires_v3_WMS           Display Field: COUNTRY           Type: Feature Layer           Geometry Type: esriGeometryPolygon           Description: World Countires Map (v3) with only few fields. No changes to the data from V1.           Definition Expression:           Copyright Text:           Default Visibility: true           MaxRecordCount: 1000          Supported Query Formats: JSON, AMF          Min Scale: 0          Max Scale: 0          Supports Advanced Queries: true          Supports Statistics: true          Has Labels: false          Can Modify Layer: true          Can Scale Symbols: false          Extent:     XMin: -2.00375070672E7     YMin: -3.02409719584E7     XMax: 2.0037507067199998E7     YMax: 1.8460513218900003E7     Spatial Reference: 102100  (3857)                Drawing Info:     Renderer:     Simple Renderer:     Symbol:     Style: esriSFSSolid      Color: [201, 242, 208, 255]      Outline:     Style: esriSLSSolid      Color: [110, 110, 110, 255]      Width: 0     Label:      Description:     Transparency: 0      Labeling Info:          HasZ: false           HasM: false           Has Attachments: false           HTML Popup Type: esriServerHTMLPopupTypeAsHTMLText           Type ID Field: null          Fields:     OBJECTID_1 ( type: esriFieldTypeOID , alias: OBJECTID_1 )     Shape ( type: esriFieldTypeGeometry , alias: Shape )     OBJECTID ( type: esriFieldTypeInteger , alias: OBJECTID )     COUNTRY ( type: esriFieldTypeString , alias: COUNTRY , length: 50 )      HTML code               <!DOCTYPE html>                     <meta http-equiv="Content-Type" content="text/html; charset=utf-8">         <meta http-equiv="X-UA-Compatible" content="IE=7, IE=9, IE=10">         <meta name="viewport" content="initial-scale=1, maximum-scale=1,user-scalable=no">         <title>Feature Layer Only Map</title>         <link rel="stylesheet" href="http://js.arcgis.com/3.6/js/esri/css/esri.css">


title: CartoDB Visualization refresh in legends does not work until you switch template

text: When changing the settings through the wizard, for example of the Bubble visuatlization, the legend does not update automatically. I managed to update the legend by changing the template to 'custom' then back to 'bubble' Is there another way to make it update?


title: Public Geoserver preview via third-party website?

text: A few years ago I remember seeing a website that would connect to a public Geoserver offering and allow you to preview all available layers overlaid on a Google Maps or similar base layer. Does anyone else remember that? Is it still around somewhere?


title: How to approximate point to segment distance on sphere?

text: I would like to approximate the distance of a line segment from a reference point on the map (coordinates are gives as lat/lon). In the Cartesian coordinate system I would simply find the closest point of the line segment to the reference point using the dot product and calculate the distance. I am new to working with geo-coordinates, and from what I have learnt, segments are called great circles on the Earth's surface. This site gives a clear formula for calculating what I would like (which is I think called the cross-track distance). However it seems an overkill for what I want to do. Firstly, I am interested in shorter distances only, maximum a few hundred kilometers (say 700 km). Second, I can do with an approximate value (say maximum 2-3% error). Third, I want the calculation as quick as possible, and the exact formula contains around 20 trig calculations. I would need something similar to _equirectangular approximation_ :               x = (lon2-lon1) * cos((lat1+lat2)/2)     y = lat2-lat1      d = sqrt(x*x + y*y) * R      Is there a way to generalize the above formula for point-segment distances? Or, alternatively, is there a handy approximation to find out the distance of a point and a great circle for small-to-mid distances with reasonable accuracy?


title: Architecture for java ee application

text: I´m new to GIS development. I have read about Openlayers, maps servers (i choose Geoserver), geodatabase (in my case i choose PostGis) and theory of gis in general. I´m developing a Java EE application that needs GIS functionality, but i´m having some problems defining the high level arquitecture for my requirements. Before i explain what alternatives i´m considering, i think it´s best if i give the main requirements that the present application must fulfill:   * It has to show a map on the browser, with several layers. Among the features that it has to show, are points that represents cities, or users, etc. In some of this features the user is allowed to click and a pop up should appear with related information, and effects like that.   * For some of the layers, the users may see only some information (i.e, some of the features) but not all, and this is based on some custom domain logic that its implemented server side with java EE.   * The users can modify some map information. They may create certain new features or modify existing ones, but all of these is based on custom logic implemented server side. They can also move through the map, and other users should be capable of seeing this updates. So, i have a lot of questions regarding how i should design something like this. My first question is:   * Should i use WMS or WFS services? As the users will be modifying things on the map, i first think of using WFS services. But, remember that all the operations that modify the map require specific domain logic in order to validate the update, before hitting in the data store. So perhaps i can access directly into the database from java and let the user retrieve the map from a WMS service. At the same time, i´m not so sure about this, because then i´ll be hitting the data base from two sides (Geoserver and jdbc or jpa or whatever) and i don´t know what kinds of problems could arise.   * Lets suppose i pick one of the previous services. How could be client-side logic be implemented? I think that with Ajax i can request the application to see if i can validate the operation on the map, and if the validation is succesful then make a call with Openlayers to the service in Geoserver (i don´t know if this is good, since the client ends up making two remote calls, one to the java EE app and other to Geoserver in order to retrieve the ¨results¨ (map)). Or how can i accomplish this? Thinking about this i´m not even sure about why i should be using a map server at all, since i have to do some geographical validations in the map before serving the geographical information.. perhaps i can just query Postgis and serve the contents directly to openlayers through the ajax call if the previous validations pass? or is it possible to pass the validation logic to geoserver somehow? As you can see, these are some fundamental questions and decisions i have to made before doing anything else, and i was hoping for advice and other alternatives as well. Thanks in advance and excuse my poor english!! :)


title: Python syntax error...cannot see the error!

text: I'm very new to python and have created a script for enabling a user to list file types in a workspace making use of the parameters function in ArcGIS so it's more user friendly. However I cannot get the script to work...I am struggling to finish the script off can anyone lend a hand? Python reports the syntax error in line 36 "f.close()" there is no white space around it so am struggling to see whats wrong...can you? Or what am I missing? Here's the code:               import arcpy     from arcpy import env     import os     env.workspace = ("C:\ArcGIS\Default.gdb")     ws = arcpy.GetParameter(0)     option = arcpy.GetParameterAsText(1)     outputFilePath = arcpy.GetParameterAsText(2)          dsList = []          if option == "List Feature Classes":         featureclasses = arcpy.ListFeatureClasses(ws)         for table in tables:             dsList.append(table)     elif option == "List Tables":         tables = arcpy.ListTables(ws)         for table in tables:             dsList.append(table)     elif option == "List Rasters":         rasters = arcpy.ListRasters(ws)         for table in tables:             dsList.append(table)     elif option == "List Workspaces":         workspaces = arcpy.ListWorkspaces(ws)         for table in tables:             dsList.append(table)     elif option == "List Files":         files = arcpy.ListFiles(ws)         for table in tables:             dsList.append(table)          f = open(outputFilePath, 'w')           for item in dsList:          f.close()


title: Animated lines and line styles in Javascript

text: I'm using Leaflet library and Mapbox. I need to create animated line from one city to another, slow motion "travelling" from one point to another point. Is that possible in javascript and how? I can add any tools if it's necessary. Also, I would need some zig zag line, not just straight line, do you have some advices regards to that? Is there some library or something similar where we can choose line styles maybe? Thank you.


title: Creating a buffer based on area (i.e. sqm) instead of distance

text: Is there a way in ArcGIS 10.1 or QGIS (or any add-on) to create buffer based on an area (i.e. sqm) rather than using buffer distance? I'd like to create a buffer of let's say 5% of the area within the buffer. The problem is that the feature is not regularly shaped. Thanks, Alex


title: Poor Performance when using Spatial Indexes in MySQL

text: Re-post of a question asked on Stack Overflow when it was suggested this would be a better forum. I'm trying a little experiment at pushing a data set which is not geo-spatial but fits it quite well and am finding the results somewhat unsettling. The data set is genomic data e.g. the Human Genome where we have a region of DNA where elements like genes occupy specific start and stop coordinates (our X axis). We have multiple regions of DNA (chromosomes) which occupy the Y axis. The goal is to bring back all the items which intersect two X coordinates along a single Y coordinate e.g. LineString(START 1, END 2). The theory seemed sound so I pushed it into an existing MySQL based genome project and came up with a table structure like:               CREATE TABLE `spatial_feature` (       `spatial_feature_id` int(10) unsigned NOT NULL AUTO_INCREMENT,       `external_id` int(10) unsigned NOT NULL,       `external_type` int(3) unsigned NOT NULL,       `location` geometry NOT NULL,       PRIMARY KEY (`spatial_feature_id`),       SPATIAL KEY `sf_location_idx` (`location`)     ) ENGINE=MyISAM;      **`external_id`** represents the identifier of the entity we have encoded into this table & **`external_type`** encodes the source of this. Everything looked good and I pushed in some preliminary data (30,000 rows) which seemed to work well. When this increased past the 3 million row mark MySQL refused to use the spatial index and was slower when it was forced to use it (40 seconds vs. 5 seconds using a full table scan). When more data was added the index started to be used but the performance penalty persisted. Forcing the index off brought the query down to 8 seconds. The query I'm using looks like:               select count(*)     from spatial_feature     where MBRIntersects(GeomFromText('LineString(7420023 1, 7420023 1)'), location);      The data going into this is be very dense along the Y dimensions (think of it like you've recorded the position of every building, telephone box, post box and pigeon on a very long road). I've done tests of how R-Indexes behave with this data in Java as well as others in the field have applied them to flat- file formats with success. However no one has applied them to databases AFAIK which is the goal of this test. Has anyone out there seen a similar behaviour when adding large quantities of data to a spatial model which is not very disparate along a particular axis? The problem persists if I reverse the coordinate usage. I'm running the following setup if that's a cause   * MacOS 10.6.6   * MySQL 5.1.46 Help!


title: Merge overlapping lines in ArcGIS

text: I have a polyline shapefile containing a number of lines, some of which overlap. By overlap, I mean that they follow exactly the same path for some of their route, but also have bits that are not shared. As it's difficult to explain in words, the following images should help: ![FirstExample](http://i.stack.imgur.com/KMPOs.png) ![SecondExample](http://i.stack.imgur.com/6S9z6.png) In the example above, you can see how the two selected lines overlap, and I want to merge them into one line feature. Is there a way to do this with a built-in ArcGIS tool? If not, is there an easy way to check for overlapping in Python? (and then an easy way to merge two lines in python?)


title: Examples of thematic mapping web services built with FOSS GIS

text: What are some current examples of web services serving thematic mapping based on free and open source software? The thematic mapping could be in the form of choropleth or point layers. It would be very interesting to see some good publicly accessible geo-information systems, with the emphasis as much on simple but well-implemented projects as on those that are technically accomplished.


title: How to add new projection using ArcObjects of ArcGIS for Desktop?

text: I want to use "Mars MOLA" { 86 3396000.0 3396000.0 "Mars MOLA" 0 0 0 0 0 0 0 } with my shapefiles.


title: Model Builder or a python script for pulling SDE data in a single GDB file

text: Would like to know if Model Builder or a Python Script would be better for pulling dozens of fc from ArcSDE into a single .gdb file ? The tool will need to be run on a regular basis.


title: name 'AddressLocator' is not defined, Python created Address Locator

text: When I run this tool, I get the error 'name 'locator_AddressLocator_n_' is not defined'. My guess is that this is an issue with the path where the locator is trying to save to? So, I am assuming I have not set the workspace or the scratch workspace? These are just guesses though. I will need to save each of the geocoding outputs, so should I just have the code create a new GDB and push the outputs to that location? How could this be accomplished?               import os     import arcpy     from arcpy import env          env.overwriteOutput = 1      Here is the geocoder code:               arcpy.CreateAddressLocator_geocoding("General - Single Field", inputFC+"'Primary Table'", "'Feature ID' '' VISIBLE NONE;*KeyField GeoLabel VISIBLE NONE;'Additional Field' '' VISIBLE NONE;'Altname JoinID' '' VISIBLE NONE", locator_AddressLocator_n_, "")      UPDATE: Created a GDB from paramter and set the workspace to that location, but I am still getting the same error when running the script. name 'AddressLocator' is not defined               inputGDBname = arcpy.GetParameterAsText(3)     inputGDBloc = arcpy.GetParameterAsText(4)          #Create GDB for address locators     arcpy.CreateFileGDB_management(inputGDBloc, inputGDBname)          env.workspace = inputGDBloc+"/"+inputGDBname      **UPDATE** New Error on the last line... RuntimeError: Object: Error in executing tool I have included all of my code this time, just so everything is visible to see what could be throwing the error.               import os     import arcpy     from arcpy import env          env.overwriteOutput = 1          #Create GDB for address locators     inputGDBname = arcpy.GetParameterAsText(0)     inputGDBloc = arcpy.GetParameterAsText(1)     arcpy.CreateFileGDB_management(inputGDBloc, inputGDBname)     env.workspace = os.path.join(inputGDBloc, inputGDBname)          #Get FC and Fields     inputFC = arcpy.GetParameterAsText(2)     inputField = arcpy.GetParameterAsText(3)     inputFieldLabel = "GeoLabel"     inputFieldName = arcpy.GetParameterAsText(4)          #Add the GeoLabel field to the Selected Feature Class     arcpy.AddField_management(inputFC, inputFieldLabel, "TEXT")          #Calculate the field values for the selected field and populate it with the inputFieldLabel + InputField     rows = arcpy.UpdateCursor(inputFC)     for row in rows:         idt = row.getValue(inputField)         exp = inputFieldName + ': ' + idt         row.setValue('GeoLabel', exp)         rows.updateRow(row)     del row, rows          Reference_Data = inputFC + ' Primary Table'     RefData = Reference_Data          arcpy.CreateAddressLocator_geocoding("General - Single Field", RefData, "'Feature ID' <None> VISIBLE NONE;*KeyField GeoLabel VISIBLE NONE;'Additional Field' <None> VISIBLE NONE;'Altname JoinID' <None> VISIBLE NONE", "Address_Loc_n_", "")


title: In ArcGIS how can I append the number of points that are located inside each individual block groups of the entire US?

text: I am using ArcGIS 10 and I have 2 shapefiles; (a) a polygon shapefile of block groups of the entire US, (b) a point shape file of customer locations throughout the US. What I want to know if how many points (b) are within each block group (a). What tool should I used to add a new field to the block group shapefile that will count and show the total number of points located inside each individual block group? I also would like to show that number as a percentage of the total number of points. Thank you for your help!


title: PHP script for publishing a layer

text: I want to publish a layer automatically in geoserver. Each layer is a database table in postgis db. I can see the table in geoserver store(already contains some layer tables of that particular db) as soon as I create it in postgis but publishing it should be done manually. Can I find any php script for just publishing the layer.... Does "Adding a postgis table" section in http://docs.geoserver.org/2.0.x/en/user/extensions/rest/rest-config-examples- curl.html also publishes the layer????


title: Drawing in ArcGIS Engine

text: I am trying to port over a large amount of code for drawing lines, bitmaps (with transparency), and text onto our map display. Originally it was done using Dynamic Display and `IDynamicMap` but that is no longer possible. Could anyone tell me the _right_ way to go about this without leveraging dynamic layers? I am completely confused... Sample code, links to tutorials, or advice would be greatly appreciated! Sorry for asking such a vague question, but the methods I have tried (`GraphicTracker`) have not been successful and I am fairly green when it comes to the ArcGIS Engine API. Thanks


title: Geoserver Rest folder not found

text: I am trying to learn geoserver's rest-api following http://docs.geoserver.org/2.0.x/en/user/extensions/rest/index.html. I have placed rest-2.0.3 jar file in geoserver/webmaps/geoserver/WEB-INF/lib folder and restarted the server. I don't know whether this completes the installation of rest. All the examples in curl rest shows the path as *http://localhost/geoserver/rest/workspaces..... but I am not able to find any rest folder and the link *http://localhost:8080/geoserver/rest giving me the page with "Geoserver Configuration API" . I have followed Using REST service with GeoServer but even this doesn't solve my problem. Can you please tell me how to find whether rest is installed or not.... *http://localhost:8080/geoserver/rest/layers/giant_polygon.html also giving me "Geoserver Configuration API".


title: Which one of these azimuth values is the right one? And why?

text: First of all, let's make clear that I want the Azimuth on the surface of the Earth, i.e. the angle between two locations, for example New York and Moscow. I am testing some azimuth calculations with my JS functions (shown below). For the points A(-170, -89) to B(10, 89), I get ~90º. JS function for Azimuth on sphere (from Wikipedia)               var dLon = lon2 - lon1;     var y = Math.sin(dLon) * Math.cos(lat2);     var x = Math.cos(lat1) * Math.sin(lat2) - Math.sin(lat1) * Math.cos(lat2) * Math.cos(dLon);     var angle = Math.atan2(y, x) * 180 / Math.PI;      JS function for Azimuth on oblate spheroid (from Wikipedia)               var dLon = lon2 - lon1;     var f = 1 / 298.257223563;    /* Flattening for WGS 84 */     var b = (1 - f) * (1 - f);     var tanLat2 = Math.tan(lat2);     var y = Math.sin(dLon);     var x;     if (lat1 === 0) {         var x = b * tanLat2;     } else {         var a = f * (2 - f);         var tanLat1 = Math.tan(lat1);         var c = 1 + b * tanLat2 * tanLat2;         var d = 1 + b * tanLat1 * tanLat1;         var t = b * Math.tan(lat2) / Math.tan(lat1) + a * Math.sqrt(c / d);         var x = (t - Math.cos(dLon)) * Math.sin(lat1);     }     var angle = Math.atan2(y, x) * 180 / Math.PI;      In Calculator 2, I get 90º. In PostGIS, I get 270º In Calculator 1, I get 180º. I know the Azimuth gets more and more distorted near the Poles, but that's exactly why I am testing at these spots. This variety of different solutions are confusing me. Could you please help me getting the right answer for this?


title: Creating Local Cache From ArcGIS Online Basemaps

text: I am contemplating creating a local caches from the standard ArcGIS Online Basemaps; streetmap, imagery and topo. I need to do this because the environment the web application will be hosted is DOES NOT have live connection to the internet. I have a couple of questions before I proceed:   1. Would creating a local tile cache from ArcGIS Online Basemaps violate any usage terms. Without doing any research (yet), I am guessing the local tile cache would be subject to the same usage restrictions of the ArcGIS Online basemaps?   2. Assuming this is legal, is there anything in place that would prevent me from doing so? Does anyone know if the ArcGIS Online servers start blocking an IP after _n_ tile requests in _n_ time period? Any other thoughts,comments,suggestions welcome.


title: geoserver-cURL-REST error

text: when I try to create a workspace/layer/store following http://docs.geoserver.org/stable/en/user/restconfig/rest-config-examples/rest- config-examples-curl.html I am getting "System cannot find the file specified".... My curl is working fine with http. What might be the solution?


title: For ArcGIS tool can Illustration of Item Description be used as Thumbnail in Contents of Catalog window?

text: I've made a Python script tool using ArcGIS Desktop 10.0 SP2 and have accessed its Item Description to Update its Illustration by browsing to a *.png file. The expected Illustration then appears in Item Description. However, when I use the Catalog window, to access Thumbnail view of the Toolset in which that tool sits, the Illustration does NOT appear as its Thumbnail. Am I perhaps missing a step?


title: Counting records that match attribute from another feature class?

text: I have geometry feature class of states that I would like to build a choropleth map from by using the quantity of people in that state. However, the data of what people in what state is in seperate table with other attributes like first and last name and so on. Is there a way to count all of the, say, Washington State people from the seperate table and append a "population" attribute at the end of the geometry feature class' attribute table so that I may then use this newly derived "population" attribute to use for the choropleth map?


title: Looking for Melbourne and Canberra metropolitan suburb shapefiles / ABS coded lists

text: I am looking for shapefiles or a list of suburb names for Melbourne and Canberra. I have ABS shapefiles, but I need to filter that to end up with only metropolitan suburbs. My shapefile has codes in it for Census Geographical Area Classification Listing \- I just need to know which suburbs belong to the metros.


title: How to convert SQL Server Rasters to ASCII Grid?

text: I have 25km by 25km squares of the world in a relational database. Each square has coordinates for the top left and lower right corner as latitudes and longitudes. Each square also has an index ilat and ilong value starting at the lower left corner with (0,0) to (maxlong, maxlat) for the top right corner. Each square also has measurements (X, Y, Z). I would like to transform this data into the ASCII Grid Format: link for each measurement. Here is some example data taken from the Internet:               ncols         419     nrows         407     xllcorner     823678.99884033     yllcorner     1151760.4752197     cellsize      14.314150820378     NODATA_value  -9999     -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999       I understand what ncol and nrows is but I am not too sure about xllcorner, yllcorner, cellsize. Could someone please advise me what these values would be in my situation? Do you reckon I can transform my data from my tile/square format into the ASCII Grid Format. * * * Actually my grid is as follows: Latitude is in the range [-90 ... 90] and longitude [-180 ... 180]. A square’s cell’s size is 2.5 degrees I have been told that corresponds to 25 square kilometres.


title: What do I need to know, to develop a plugin for QGIS?

text: I wish to develop a simple plugin for QGIS, which will open some data, and some specific processing. I have quite a it of development experience with ESRI's Arcojects, but haven't developed with QGIS before. I tried searching on the QGIS wiki, but only found this link : http://www.qgis.org/wiki/Writing_C++_Plugins which seems to be obsolete. There were also a few blog posts from 2006. I also saw something about Python plugins, but I am not sure whether I need to work with the C++ plugin or a Python plugin. Can Anyone guide me to the documentation, or tutorials that would help me get started


title: How do you set the position of a plugin window when you start the plugin?

text: I'd like to set the position of a plugin's window or dialog **relative to the map canvas** frame. For example, if I start a plugin, I want it to appear somewhere in the left upper corner of the map canvas. I can use absolute screen coordinates (choosing appropriate values for the first two parameters of the `setGeometry` method of a dialog), but that's not a flexible solution as users might adapt the size of the map canvas frame according to their wishes. FYI: I use the Plugin Builder plugin to take off... * * * Based on @Matthias Kuhn's concise answer, I found a solution which I share here below. Note that the code corresponds to the script that is produced by the Plugin Builder plugin. Also note that the argument of `mapToGlobal()` must be a `QPoint`.               def run():         self.mapcanvas_upperleft = self.iface.mapCanvas().mapToGlobal(QPoint(0, 0))         offset_x = 24         offset_y = 43         pos_x = self.mapcanvas_upperleft.x() + offset_x         pos_y = self.mapcanvas_upperleft.y() + offset_y         self.dlg.setGeometry(pos_x, pos_y, 0, 0)              # show the dialog         self.dlg.show()         # Do this a first step... then delete         print self.dlg.pos().x()         print self.dlg.pos().x()       In order to define the offsets, you can print the upper left corner of your dialog **before** you change its position with `setGeometry`. The `pos()` method returns a `QPoint`. (hint: in QGIS, click Plugins > Python Console to see the print output) As such, you know the difference between the global coordinates of 1) the upperleft corner of the map canvas and 2) the upperleft corner of the plugin dialog, which then helps you to define the offset you prefer.


title: How can I add menus to Flex Viewer 2?

text: I have a lot of widgets and my dashboard in Flex Viewer 2 has gotten cluttered. How can I add menus to the dashboard rather than buttons to organize the widgets? Any help will be appreciated.


title: How can I calculate the vertical distance to the nearest water?

text: Summary: I would like to create a raster where each cell contains the vertical distance to the nearest water (a decent approximation would be acceptable). The problem is that I'm not sure how to go about it. Current approach: I've looked through the documentation for my GIS (I'm using Grass 6.4 if that matters). And I can't find any direct way to perform this calculation. I've also searched the web without luck. So I dreamed up the following:   * convert my water vector to a water raster (1=water, null=everything else)   * replace the 1s in the water raster with the altitude above sea level   * perform "neighborhood" (moving window) analysis. Calculate the average altitude of water for a 100 meter radius for every cell in my region (ignoring nulls)   * repeat the above with 100m, 200m, 300m, ..., 2000m. One raster for each distance.   * using r.grow.distance calculate the distance between each cell and it's nearest water   * then, using map algebra, I look at the distance to nearest water for each cell and its own altitude and subtract the average altitude of water for the appropriately distanced raster. The problems with my method:   * It is very computationally demanding. I am working at a 10m resolution and I need to perform this modelling on around 60 million cells. I've been creating rasters on a fast desktop for days and days and I'm not even close to finished.   * At small radii, the calculation seems reasonable. If the distance to nearest water is 75 meters and I use a 100 m average altitude to water, the approximation is pretty good. But if the distance to nearest water is 1875 meters and I use the 1900 meter raster I will be averaging in the altitude to water of many cells far removed from the nearest water to my cell. The results are dissatisfying. I considered using some kind of smoothing process for the final raster but I wonder if there is a completely different approach that I am not seeing. I've been reading about predictive modeling projects where vertical distance to nearest water is a one of the inputs to the model so I know it can be done but I'm just not sure how. I'd appreciate any help I can get with this problem. Thanks.


title: How to can I run the OpenGeo Suite Community Edition under Apache Tomcat on Windows server 2003?

text: I have installed opengeo suite **community edition** on ubunutu 10.10 and it is working with apache tomcat and it all works great. So I decide to do the same on a Windows 2003 server at work. _Note: On the windows 2003 server machine I already have apache tomcat installed and working correctly on port 8080_ So what I did was download the opengeo suite **community edition** for windows from opengeo.org and run the installer. I change the default port to 8081 and everything works fine! However it installs using a different web server (JBOSS?? or J something?), and this server only seems to operate when I am logged onto the the machine which is undesirable. For example I can't access geoexplorer from a remote machine unless I am also logged onto the server. Has anyone got the openGeo suite **community edition** running with Tomcat? and how did you do it? Thanks Ando


title: Overlay polygons on raster colour map to extract area values for each colour within each polygon

text: Does anyone know an answer to the following: I have a set of polygons created as the MCP of animal movements. I also have a habitat colour map raster, with 13 different habitat types as different colours. I need to find a way of overlaying one on the other and calculating the area of each habitat type within each polygon. Possible? I think it has something to do with overlaying, but have not used it before, and usually asking those who know is far faster than self-discovery! Thanks in advance.


title: Producing a map of higher level regions using postal codes

text: Say I have a map of Australian postal codes (technically a postal code shapefile downloaded from the Australian Bureau of Statistics (ABS) Web Site). I want to use this map to construct a higher level map using proc gmap in SAS. For example, a number of postal regions will form a larger region and I only want to show this larger region. Is this possible? The reason I need to do this is that the ABS does not provide a map for this higher level region and so I need to do this myself. Of course, this is possible in MapInfo, but I want to do it in SAS.


title: Creating one DEM from multiple ASC files

text: I am struggling with this several days and now its time to ask how to do it in QGIS. I have very detailed 27 ASC grid files with overall size of 9 GB created from 1m contours of large area (600 m2). The goal is to merge all this data into one DEM for easier handling (output will be smaller i think ?) and later to create a slope, hillshade color shaded reliefs ,etc for whole area easier. Available Data:   1. I have simplified contours for whole area in SQL Sever which are much less in size than original ones.   2. ASC files are created from original contour data and I would like to create one DEM from that. I have tried next QGIS/GRASS solutions with simplified contours of a whole area:   * How to generate a DEM from a contour Shapefile?   * Is it possible to make a DEM from contours in QGIS?   * http://grasswiki.osgeo.org/wiki/Contour_lines_to_DEM but I had only success in rasterizing contours (v.to.rast.attr). That creates a 218 MB TIF with dimensions of 37945 x 30143. But, in the next step when I start 'r.surf.contour' it just hangs forever (last try was 24 hour long and then I quit, although process was very active and consumed almost all memory (6GB RAM). Is this too large dataset and this is normal execution time on modern PC configuration ? Other way would be to use ASC files and to merge it to one DEM. If someone could give me information how to do it efficently or some detailed workflow would be great. I would say that the purpose of DEM will be to create a shaded relief as a background for printing 1:5000 maps. **UPDATE:** Well I have converted all 27 ASC to GeoTiff and also realized there are 3 more grids (.grd) which has been made by GS Surfer 12. I converted them also in GeoTiff. Now when I try to build VRT next thing happens:   * The first of 27 GeoTifs are joined into VRT but I've lost all elevations. VRT have only one value (167, which is the lowest value from all dem's included).   * Those 3 files built from Surfer GRD have somehow different characteristics and it wont join with others. **               > Warning 6: gdalbuildvrt does not support heterogenous band     > characteristics. Skipping D:/DEM (GeoTIFF)/OF_IZO_28_dem.tif      ** gdalinfo for one tif in original group:               Driver: GTiff/GeoTIFF     Files: D:/DEM (GeoTIFF)/OF_IZO_27_dem.tif            D:/DEM (GeoTIFF)/OF_IZO_27_dem.tif.aux.xml     Size is 5001, 5001     Coordinate System is `'     Origin = (7584999.500000000000000,4795000.500000000000000)     Pixel Size = (1.000000000000000,-1.000000000000000)     Image Structure Metadata:       COMPRESSION=PACKBITS       INTERLEAVE=BAND     Corner Coordinates:     Upper Left  ( 7584999.500, 4795000.500)      Lower Left  ( 7584999.500, 4789999.500)      Upper Right ( 7590000.500, 4795000.500)      Lower Right ( 7590000.500, 4789999.500)      Center      ( 7587500.000, 4792500.000)      Band 1 Block=5001x1 Type=Float32, ColorInterp=Gray       Min=243.000 Max=834.072        Minimum=243.000, Maximum=834.072, Mean=503.092, StdDev=130.621       NoData Value=-1       Metadata:         STATISTICS_MAXIMUM=834.07208251953         STATISTICS_MEAN=503.0916431759         STATISTICS_MINIMUM=243         STATISTICS_STDDEV=130.62100147903      gdalinfo for one tif made from Surfer .grd>               Driver: GTiff/GeoTIFF     Files: D:\DEM (GeoTIFF)\OF_IZO_31_dem.tif            D:\DEM (GeoTIFF)\OF_IZO_31_dem.tif.aux.xml     Size is 2125, 3512     Coordinate System is `'     Origin = (7587875.569015772100000,4790000.499983053700000)     Pixel Size = (0.999968455743971,-0.999966106522373)     Image Structure Metadata:       COMPRESSION=PACKBITS       INTERLEAVE=BAND     Corner Coordinates:     Upper Left  ( 7587875.569, 4790000.500)      Lower Left  ( 7587875.569, 4786488.619)      Upper Right ( 7590000.502, 4790000.500)      Lower Right ( 7590000.502, 4786488.619)      Center      ( 7588938.036, 4788244.560)      Band 1 Block=2125x1 Type=Float64, ColorInterp=Gray       Min=409.000 Max=1239.595        Minimum=409.000, Maximum=1239.595, Mean=679.431, StdDev=154.121       NoData Value=1.7014100000000001e+038       Metadata:         STATISTICS_MAXIMUM=1239.5948358501         STATISTICS_MEAN=679.43099985765         STATISTICS_MINIMUM=408.9999998738         STATISTICS_STDDEV=154.12086745843      How I would correct this and get all grids in VRT with preserved elevations ?


title: Install pgRouting on Postgres 9.1 and Win 32

text: I referred some documents on GIS-stackexchange and over internet regarding installation of pgRouting on Postgres v9.1 and windows 32. I am not getting the correct version of pgRouting to install. Can anyone guide me to get the right pgRouting version for Postgres 9.1 and windows 32. I tried winnie builds and sanak builds but I am getting errors after running SQL query CREATE EXTENSION pgrouting;. The error are:   1. Using winnie build for Postgres 9.1 and win 64 ERROR: incompatible library "C:/Program Files/PostgreSQL/9.1/lib/librouting.dll": version mismatch   2. Using pgRouting for Postgres 9.2 and win 32 ERROR: incompatible library "C:/Program Files/PostgreSQL/9.1/lib/librouting.dll": version mismatch   3. Using Sanak build: ERROR: could not load library "C:/Program Files/PostgreSQL/9.1/lib/librouting_bd.dll": The specified module could not be found.


title: Can I speed up very long connection times to a PostGIS DB caused by complex views?

text: When connecting to a PostGIS database to add a postgis layer, all table and views are scanned by Qgis in order to define the geometry type (if I understand well the query performed). My problem is that I have complex views, quite long to be queried, that are not meant to be queried without a restricting WHERE clause. The automatic Qgis queries thus take more than half an hour to perform, and the application is frozen in the mean time. Is there a way to tell Qgis to drop the query if no answer was obtained before a timeout, and to assume the geometry type is unknown for these layers ? Thank you for your answers ! Diane


title: r raster package rasterize failing:no slots

text: I'm having trouble rasterizing a shapefile using the raster package in R.               shp<-shapefile(ZoneShape);     lcRas<-raster(lcRaster);          r<-raster(ncol=ncol(lcRas), nrow=nrow(lcRas), crs=CRS);     res(r)<-res(lcRas);     extent(r)<-extent(lcRas);          >r          class       : RasterLayer      dimensions  : 22610, 27959, 632152990  (nrow, ncol, ncell)     resolution  : 1, 1  (x, y)     extent      : 554739, 582698, 3837197, 3859807  (xmin, xmax, ymin, ymax)     coord. ref. : +proj=utm +zone=15 +ellps=GRS80 +datum=NAD83 +units=m +no_defs +towgs84=0,0,0          shp$GID<-1:nrow(shp);          > shp          class       : SpatialPolygonsDataFrame      nfeatures   : 1      extent      : 554838, 582597.6, 3837297, 3859707  (xmin, xmax, ymin, ymax)     coord. ref. : +proj=utm +zone=15 +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0      nvariables  : 14     names       : SP_ID,              NAME, Shape_Leng, Shape_Area, GID          >zoneRas<-rasterize(shp, r, "GID")      The following error is returned:               trying to get slot "coords" from an object of a basic class ("NULL") with no slots.      Can anyone see what I'm missing/screwing up here? Thanks


title: How to calculate distances in a point sequence?

text: I'm looking for help calculating the distance between a sequence of points that are in a single Shapefile in QGIS. Below is what my data look like and a blank distance column I've added to show how I'd like the distance to look. I want to know the distance between point 1 and 2, 2 and 3, etc. I'd like distance to be in meters or Km, but currently my Shapefile is in a projection whose units are decimal degrees.               ID  LAT         LON         TIME        DISTANCE     1   10.08527    124.59833   21:24:37    0     2   10.08523    124.59830   21:25:07    ?     3   10.08526    124.59832   21:25:37    ?     4   10.08526    124.59831   21:26:07    ?      A number of people have asked similar questions, but none quite get at what I want to do. This post is close, but this is in PostGIS, not QGIS Calculate distances between series of points in postgis This post got me part of the way there, but since I'm new to QGIS the answer doesn't provide enough detail for me. For example, once I've installed the GRASS plugin, I think I need to save/import the Shapefile with my series of GPS points into GRASS so I can use the v.distance module, but I don't know how to do that. QGIS calculate distance of point along a line Is the GRASS v.distance module the only way to go? Or is there a more straight forward way? If v.distance is the only way would someone be able to point me to or explain more step-by-step how to do this?


title: SQL 2012 - How can I convert MultiLineString to Geography?

text: I have a Geometry LineString that I wish to convert to Geography. Since the Geometry isn't valid , When I use MakeValid on it, SQL convert it to MultiLineString and than the Geography creation is failing with error :               The specified input does not represent a valid geography instance. Use MakeValid to convert the instance to a valid instance. Note that MakeValid may cause the points of a spatial instance to shift slightly.      What can I do in order to convert the Geometry to Geography ? Geometry text : http://pastebin.com/Ly8ex41U Geometry text after MakeValid : http://pastebin.com/HGNFKnS3


title: ModelBuilder doesn't add my output feature class to TOC

text: I'm trying to create a tool in ModelBuilder that will take a new point from a user, create a feature class for this point and display it on a map. What I've done is:   1. I've created a feature class for a point-like feature in ArcCatalog.   2. I've prepared a really simple model in ModelBuilder, that starts with a feature set variable (which takes the feature class as a schema and symbology file). The variable is marked as a parameter.   3. Next I'm using the Copy Features tool that takes the variable as input and saves it as a feature class in my scratch workspace.   4. The output feature class is marked as a parameter with the 'add to display' option. After the module is run I get a shapefile in my scratch workspace. I can add it manually to TOC and it contains the selected point. However no matter what I do, I can't display it in TOC automatically after the model completes. There is a known issue described in several places that doesn't allow to display newly created layers in TOC when model finishes with the 'add to display' option. However there's a quick workaround, which is to mark the output layer as parameter. In my case it doesn't work. I'm using ArcGIS for Desktop 10.1 on Windows 7 x86_64.


title: Align OpenLayer raster with OS Mapping Vector

text: I have an OS Vector map of a town in Ireland. When I load the Google Satellite Raster from the Open Layers plugin, the results are slightly off, by about 2m, probably due to ariel persepective. I'm trying to digitise information from the Satellite raster onto the OS Mapping and it would be easier if I could force the Satellite raster to move temporarily so they align in that persepective. Is this possible? Thanks in advance.


title: How to connect to MSSQL-Server using Geoserver?

text: I have Geoserver on Ubuntu 13.04. On another machine I have MS SQL-Server 2008 R2. Now I want to connect it together. Google says that its possible to connect to SQL-Server from Ubuntu using ODBC but Geoserver works through JDBC. In Windows its easy, you just to need to download MS JDBC driver. But what do you do in Ubuntu?


title: GRASS v.rast.stats missing in sextante toolbox

text: I am new in the QGIS world and also with GRASS. I am really interested in using v.rast.stats as is doing exactly what I need. My problem is that I am not able to find it inside sextante in QGIS and I don't now how to make it appear. Can anybody help me?? Thank you in advance!


title: ggmap, ggimage… - Saving and loading maps in R

text: In my current R project I'm working with the ggmap package to download and display maps. What bothers me is that you need internet access to load maps from google etc. I'd like to load a map from the internet and save it so I can work offline as well. An ideal solution should minimize the amount of code between get_map() (or a related command) and using the saved map within ggmap() (or a related command). I know already that get_map() has a filename attribute to save maps as an image, and I also know the ggimage() method (I've never actually worked with it, though.) Is there a simple way save such a map, or are these two the best tools I have, and if so, how do you use them effectively?


title: Model builder to call python script

text: So I'm looking to write a script that can be called by model builder. What I would like this script to do is to select three records from a feature class that have the highest field values ('AverageGasPrice'). Basically, I want the script to select the three most expensive gas prices. I am not sure where to start, however. Is there also a way to have the model launch and allow the user to enter which Feature Class they would like to specify as the input FC? If this is not possible, I will also require assistance with how to 'hardcode' the input file using relative path (everything will be stored in one file geodatabase). Thanks for the help!


title: How do I create a heatmap in ggplot?

text: I'm using R with ggplot2 to generate maps from shapefiles and from data in SQL Server. I'm trying to get to a stage where I can automate a large chunk of the mapping I do. So far this is working well and I'm able to generate some fairly nice looking maps (just need to work on labelling). I'm wondering if it's possible to create a geographical heatmap from a list of point co-ordinates (which are irregularly spaced)? Similar to using geom_density2d, but raster.