title: How to make a freehand polyline option from drawing toolbar enabled throughout in arcmap

text: Following code selects freehand polyline as default draw tool on first go: final ICommandBars commandBars = app.getDocument().getCommandBars(); final UID windowId = new UID(); windowId.setValue("{57E77868-D015-11D2-9F31-00C04F6BC8DD}"); final ICommandItem commandItem = commandBars.find(windowId, false, false); app.setCurrentToolByRef(commandItem); But once we finish with drawing first element, it goes back to its original settings i.e. it starts showing 'Polygon' as default draw option. Is there a way to set drawing tool to freehand polyline option throughout application execution in Arc map Regards, Viji


title: pg_escape_string() equivalent?

text: A client is hosting one of my maps on a server that does NOT have a postgres module installed, and will not install it (nor anything else... ) . We have been told also by the site owner that we can't move it due to wanting to retain the IP etc for seo (etc etc). so obviously my first call to escape a string for CartoDB results in "Call to undefined function pg_escape_string() in /usr/ etc etc... Is there a preferred method for circumventing this? it looks like there may be some remedy in PEAR , but it would be very preferable to not have to depend on being able to do /anything/ on this particular server.


title: NDVI time series in R -- how do I create a time series object from a data frame?

text: I am doing a time series analysis of NDVI using the bfast package in R. Bfast requires a time series (TS) object. I have successfully created a data frame with the NDVI at various point locations for tiffs in a given directory (code for this is at the bottom of the post). Problem 1: the code does not extract the date of the images. Landsat file names give year and Julian day starting at the 10th character (YYYYDDD) so this should be possible. But how? And will the Julian day be recognized as a date? Problem 2: once the date is extracted, how do I create TS object from the data frame? Thanks for any tips! I'm new at R, so be gentle. Code for extracting values at point locations:               require(raster)     require(rdgal)     require(sp)     setwd("C:\\Users\\Work\\LSTest")     pts <- readOGR(getwd(), "MyPoints")     LS.files <- list.files(getwd(),pattern="\\.tif$", full.names=FALSE)     India.LS <- stack(LS.files)     pts@data <- data.frame(pts@data, extract(India.LS, pts))


title: Convert an .xlsx datasheet with point feature data into a POINT feature class

text: I am trying to create a script tool to automate this process, eventually it will be used for batch processing but at the moment I can't even get test data to work. This is pretty basic at the moment, and I know it can be done through ArcGIS, but I would like to develop a python script for processing. My code is as below: It completes without error, but I am expecting the creation of a feature layer (even better would be a feature class with POINT geometry). I can't find or access the output of the SaveToLayerFile step. I apologise for the messiness of the code, I've been trying to beat this for 2 days now.. Any help would be appreciated.               #Import libraries     import arcpy     import os     import sys     arcpy.env.overwriteOutput = True          #Set workspace     home = arcpy.env.workspace     arcpy.env.workspace = arcpy.GetParameterAsText(0)     mxd = arcpy.mapping.MapDocument("Current")     df = arcpy.mapping.ListDataFrames(mxd, "*")[0]               #Create spatial reference object     prjFile = arcpy.GetParameterAsText(2)     spatialRef = arcpy.SpatialReference(prjFile)     arcpy.AddMessage('Setting Spatial Reference')          #Set Excel Datasheet for Migration     TEST = arcpy.GetParameterAsText(1)          try:     # Process: Table to Table Conversion with destination variable         FarmerTable = 'Farmers1'         arcpy.AddMessage('Migrating Data...')              arcpy.TableToTable_conversion (TEST, home, FarmerTable)     except:         arcpy.GetMessages()               # Set local variables     in_layer = "FarmerLyr"     out_layer = "FarmerLyr.lyr"     out_layer0 = "FarmerLyr"               try:        #Make XY Events of each Farmer         arcpy.AddMessage('Creating Features...')              result = arcpy.MakeXYEventLayer_management(FarmerTable, "X", "Y", out_layer0, spatialRef, "")         result.getOutput(0) = out_layer0         arcpy.AddMessage(arcpy.GetCount_management(out_layer0))         arcpy.AddMessage('Successfully Created Features')              # Execute SaveToLayerFile         arcpy.SaveToLayerFile_management(out_layer0, out_layer, "ABSOLUTE")     except:         print arcpy.GetMessages()


title: Is there demand for a higher resolution alternative to GPS tracking in the GIS market?

text: If there was an alternative to GPS that provided higher resolution tracking, what applications would there be for it in GIS? Note that this means higher resolution, but not necessarily higher accuracy. For instance, the absolute position (lat/long) would not be more accurate than GPS, but the relative movement during tracking could be within a few inches. This means you could get a very detailed track recording of a device's movement, allowing you to see small changes in position or elevation. If you were to carry this device while hiking a trail, you could see elevation changes with 1-5 inch resolution. I imagine this might be useful for GIS, but I am not familiar with the market. Any insight would be appreciated.


title: Thermal Map Data

text: I am curious about making a thermal map to display where the hot spots and where the cool spots in a community are. Obviously the roads will be hot and the shaded areas with trees will be cool but this is more to brush up on my skills. My main question is how would I go about collecting this data for an area where none currently exists? It will be some form of raster data but I am lost from there on out. Does this relate more to remote sensing? If anyone could help me get the ball rolling I would greatly appreciate it.


title: Iterating using search cursor and layer selection performance?

text: I'm trying to use the results of a search cursor to iteratively make selections on a feature class and then create new feature classes (so you get a new feature class or shp for each unique value in a field... kind of an explode by attribute tool). _As an example - if you had a feature class containing 14763 records and one of the fields had 18 unique values in it, after running the script you would have 18 new feature classes in your gdb._ So far I have this:               def unique_values(table, field):     with arcpy.da.SearchCursor(table, [field]) as cursor:         dict = sorted({x[0] for x in cursor})         print type(dict)         for i in range(len(dict)):                  whereD = dict[i]             fname = whereD[0:5] + "_" + whereD[-4:]             arcpy.MakeFeatureLayer_management(table, fname)             query = "relatedBirds_Species = '" + whereD + "'"             arcpy.SelectLayerByAttribute_management(fname, "NEW_SELECTION", query)             arcpy.CopyFeatures_management(fname,fname+ "_c")      This works but it's incredibly slow when it runs. Does anyone know of a quicker way to get the same results? **Requested edits below** _The Data_ Geodatabase workspace The feature class I'm running this with has 18 unique values in the field in question. Total number of records in the original feature class is 14,673 _Timings_ As the code is above = 2' 35.04" With the {where clause} added to MakeFeatureLayer (w/out copy feature either)= 1' 27.75" (thanks Stephen) With Select_analysis = 2' 10.58" I think most of the time was being taken to actually draw these in Arc as I was running this directly from the python window so I ran it again from cmd with arc closed. Stephen wins at 57.93" (with copy feature put back in) running from cmd.


title: How to parse xsd:schema file using openlayers?

text: I got layer column names using this thread. But I am stuck parsing it. How to parse this xml using openlayers?               <xsd:schema elementFormDefault="qualified" targetNamespace="http://www.openplans.org/topp">         <xsd:import namespace="http://www.opengis.net/gml" schemaLocation="http://localhost:11032/geoserver/schemas/gml/3.1.1/base/gml.xsd"/>         <xsd:complexType name="statesType">             <xsd:complexContent>                 <xsd:extension base="gml:AbstractFeatureType">                     <xsd:sequence>                         <xsd:element maxOccurs="1" minOccurs="0" name="the_geom" nillable="true" type="gml:MultiSurfacePropertyType"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="STATE_NAME" nillable="true" type="xsd:string"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="STATE_FIPS" nillable="true" type="xsd:string"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="SUB_REGION" nillable="true" type="xsd:string"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="STATE_ABBR" nillable="true" type="xsd:string"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="LAND_KM" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="WATER_KM" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="PERSONS" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="FAMILIES" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="HOUSHOLD" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="MALE" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="FEMALE" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="WORKERS" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="DRVALONE" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="CARPOOL" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="PUBTRANS" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="EMPLOYED" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="UNEMPLOY" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="SERVICE" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="MANUAL" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="P_MALE" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="P_FEMALE" nillable="true" type="xsd:double"/>                         <xsd:element maxOccurs="1" minOccurs="0" name="SAMP_POP" nillable="true" type="xsd:double"/>                     </xsd:sequence>                 </xsd:extension>             </xsd:complexContent>         </xsd:complexType>         <xsd:element name="states" substitutionGroup="gml:_Feature" type="topp:statesType"/>     </xsd:schema>


title: Problem in Spatialite database creation

text: When I try create a new DataBase, Spatialite plugin return the dialog error box:               The SQL query seems to be invalid. no suck table: views_goemetry_columns      When I click in OK button, a second dialog box (Transform DataBase into SpatiaLite) is open with the text :               This DataBase seems to be a valid SQLite DataBase but not a valid SpatiaLite One     Would you like QSpatiaLite to automatically transform it to SpatiaLite (information won't be lost) ?      When I click in the YES button, a new dialog error box open with the text:               Unable to convert Database to SpatiaLite      When I click in OK button, another dialog error box open with the text:               The SQL query seems to be invalid.      The DataBase was set up, but I cant import spatial layers. How resolve this problem?


title: How to convert a point grid to a raster?

text: If I load the original .csv file, I then have a map of the African continent represented by 14069 data points covering it at 40x40km resolution. Each data point shows solar radiation means across the year, (jan, feb, etc). What I want to be able to do is convert that vector into a raster using GRASS. What I have managed to do already, is import the .csv (tab delimited) file into the GRASS mapset using v.in.ascii. However, when I try use v.to.rast.constant to convert that vector into a raster, I get a solid square covering Africa, that cannot be edited. I.e. there is no data from the vector within the created raster. What am I doing wrong? This work is for my masters thesis and I'm constantly working on it. Thanks!


title: Repaint seems not working when called multiple times from the same method

text: I have the following scenario:   `   displaySelectedFeatures(IDs);   then some code... and then:   displaySelectedFeatures(features);    Display selected features looks like this:               public void displaySelectedFeatures(Map<Color, Set<FeatureId>> features) {              Style style;              style = createSelectedStyle(features);              Layer layer = mapFrame.getMapContent().layers().get(0);         ((FeatureLayer) layer).setStyle(style);              mapFrame.getMapPane().repaint();         }         What I need is the following. First IDs object should be drawn and then after some time e.g. 5 seconds, object - features should be repainted. What happens it that only the IDs is drawn while features(second repaint) is omitted. How to solve this?


title: How to install PostGIS 1.5 on Mac OS 10.5?

text: I am running PostgreSQL 9.0 on Mac OS 10.5. I've installed postgis- pg90-1.5.3-1-osx through "Stack Builder" but still don't see postgis in my list of databases when I open PgAdmin. When I open my Qgis 1.7, I can connect to my postgres database and display all tables with no spatial geometry. However, when I try to query a non-spatial table, I get the following message "Your database has no working PostGIS support". I appreciate your ideas to get postgis to integrate into postgresql. Thanks! * * * Thanks for your ideas! Unfortunately, the first link (http://www.kyngchaos.com/software/postgres) is limited to "...Lion and Snow Leopard, 64bit only." I'm using Leopard (Mac OS 10.5). The second link (http://trac.osgeo.org/postgis/wiki/UsersWikiPackages) provides two options. Both of these (http://yum.postgresql.org/repopackages.php) & (http://www.postgresonline.com/journal/archives/204-postgis15-install- yum.html) point to Linux installations. So, I'm still stuck. I've downloaded PostGIS 1.5 using Application Stack Builder 3.0.0 from PostgreSQL 9.0. It installs, but does not run. I appreciate any other suggestions.


title: How to convert WMS to SHP?

text: I'm looking for tool that will be convert WMS to vector layer. Are they ready tool?(for example in python)


title: Which parameter values should I use with QGIS Heatmap plugin?

text: I've read all the posts and docs about the heatmap plugin but can't seem to get it to work. The map is the size of the United States and has a little over 2000 points. I've tried small radi (1 meter) and large radi (5 million meters), but still get a red hot map. I'm using QGIS 1.8 and WGS84. The data was initially in a spreadsheet, but I saved it as a .shp vector layer and deleted the old csv layer. The goal is to convert that to contour lines to map density of these points like in this Mapbox blog post.


title: looking for an esri compatible javascript engine

text: I work on an esri heavy team, i'm looking for a mapping engine that will give me flexibility to build on top of and leverage only the pieces I need and or add more as the app grows. I've been using ESRI's API to great success for my indiviual code but leveraging their codebase into a reusable solution grows too large. We have arcgis server 10.1 and create featureserver and mapserver services. Along with those I will consume data thats created from various tabular locations and sent to my front-end via JSON/Ajax. I've started on leaflet, and have began testing on esri-leaflet to some success, dynamic layers are drawing slowly, but we are currently upgrading to 10.1 so that may wash out. Goal: find a lightweight javascript based mapping engine that will handle ArcGIS feature&mapserver services. Along with the ArcGIS services the engine needs to handle JSON data. Geocoding with the service of my choice would also be great(bing, google, esri, or internally created) Thx


title: Is there a way to reference a selection in model builder or python?

text: Is there a way to reference a selection in model builder or python so that I can edit my python script to only run if there is a selection. And even better if there is only one row selected? I made a model to make a selection with python script to ask for variables, make a selection based on those variables, and then zoom to the selection. When the tool does not make a selection for whatever reason it takes forever to finish and zooms to the entire extent of the map or crashes ARCmap entirely.


title: How to create a text box with inverted rounded corners

text: I just saw this post on the tex SE site and wondered how I would do the same in ArcGIS. I created a text box in a layout and drilled down in all of the dialog boxes but couldn't see anything similar. Sure, it's easy to create a rounded box, but how do you invert the corners? One solution I had was to manually draw a rectangle with circles on the corners. ![enter image description here](http://i.stack.imgur.com/td87c.png) Then, I used the _subtract_ tool on the Graphics Operations context menu for each of the corners. ![enter image description here](http://i.stack.imgur.com/eUREn.png) The result was satisfactory: ![enter image description here](http://i.stack.imgur.com/5icEm.png) I could easily put some text in there and it would look fine. However, if the amount of text changes and I need to stretch the box, the corners stretch as well: ![enter image description here](http://i.stack.imgur.com/BB8Ly.png) Does anyone know of a way to create such a frame around text in the layout in ArcGIS? I'd be happy to know if it's possible to do this using Python.


title: In ArcGis 10, spatial join of a hot street value to a fishnet cell is counting multiple values for same street segment

text: So, I am working on a fishnet, so that we can attribute multiple values to points within different cells (i.e., ranking light poles as optimum for traffic cameras based on the number of speeding tickets, accidents, etc) so that all poles in a high value cell would be ranked as prime candidates. To do this, we have assigned numerical counts of accidents to street segments (sections between intersections) and are attempting to join those values to the fishnet through a spatial join. The issue is, if a street segment runs through the cell, and has a value of 5, it is counted 15 times (or whatever for each cell) and the cell then has a value of 75 instead of 5). Any way to fix this and get one value from each unique street segment in the cell? (averaging as a run around doesn't work, as some cells have different street segments and if one street has a 0, we don't want that to discount the value of the other street that may have a higher value.) We've tried dissolving, interection, multi to single, and none of them seem to make each street segment unique or with only the requisite value display once. Thanks. Chris NYC as requested, here is a shot of the image. the issue is, and this example actually accidentally illustrated it clearer than I'd seen, is that for the horizontal street near the middle, with a "3" (there are only two dots because one has 2 values) looks correct, but if one tried to join that "3" value to a cell intersecting that street segment, it would likely join the value "45." This is demonstrated a street over where one street has a "3" a "0" and a "7." This seems to demonstrate the single street is actually 3 segments to arcgis. I can't figure out how to make, for example, that street equal "10" and then join it's value to a cell so the cell would only have "10" as the joined value, not "30" or "150" or however many times it joined. ![map example](http://i.stack.imgur.com/GOIX9.png)


title: Custom python scripts do not appear in SEXTANTE QGIS

text: Custom python scripts do not appear in Scripts menu in SEXTANTE plugin unless they are either empty or are as simple (and useless) as:               ##foo=group     ##input=raster      OS is openSUSE 12.3; QGIS 1.8 and Master both demonstrate the same behaviour. I was told that the issue may be caused by insufficient permission rights for the scripts folder. The path to the script's (for Master in this case) is that is defined by default: `/home/username/.qgis2/processing/scripts`. For both 'processing' and 'scripts' folders I have the following permissions: `drwxr- xr-x`. For the scripts in 'scripts' folder the permissions by default are `-rw-r--r--`, but even if I enable execution rights manually to be `-rwxr- xr-x` this doesn't affect their appearance in SEXTANTE scripts. So I suppose that permissions doesn't seem to cause trouble here... Please advise me how to overcome this issue.


title: ArcObjects SDK: Textured MultiPatch: problems with large texture images

text: I am using IGeneralMultiPatchCreator to create a textured MultiPatch. This works fine for small texture images, but fails for larger images (15 MB .png). In ArcScene the MultiPatch is displayed in grey, the same way it would be if the "Disable material textures" checkbox on the Rendering tab of the Layer Properties was checked. I have tried to increase the max texture memory for the layer:               I3DProperties2Ptr prop3D = layer;     long maxTextureMemory = 128*1024*1024;     prop3D->put_MaxTextureMemory(maxTextureMemory);     prop3D->Apply3DProperties(layer);      but this does not seem to have any effect :-(. I have also tried to split the MultiPatch into parts where each part have a smaller texture image, but this does not seem to work either.


title: Creating layer from in-memory shape data using ArcObjects?

text: In my current project, a feature layer is created based on shape data generated after some complicated surface analysis procedure. The sequence right now is like [Surface Analysis Algorithm] -> [array of xyz data] -> [shp file written to disk] -> [create feature class from shape] -> [create feature layer from feature class] I'm on the process of optimizing the project, so I wondered if there was any way to display directly the in memory shape data in a layer, bypassing the read/write to hard disk?


title: Is it possible to show Pie charts in QGIS Legends?

text: I would like to display some pie charts on my map. They show in the map frame perfectly but in legend on the print composer I am unable to show the legend for the pie charts. QGIS just shows the style of the point underneath the pie. Is it possible to get these shown? Using QGIS v1.8 from OSGEO4W on Windows XP


title: Geometric Network error in File Geodatabase<Network Dataset

text: I need to create a geometric network in my network dataset which resides within my file geodatabase. I have been receiving the following error: Cannot create geometric network. _The selected feature dataset does not contain any feature classes which can participate in a new geometric network._ This error should not come up because there are tons of functional feature classes in the network data set. Any workarounds or solutions?


title: How to extract vegetation polygons from Imagery using Python?

text: I have been trying to roughly polygonize vegetation areas from Imagery (RGB). I tried using Iso Cluster Unsupervised Classification (Spatial Analyst) but it is still too detailed for my purpose. I just need to roughly draw polygons wherever I have dense vegetation as can be seen here:![enter image description here](http://i.stack.imgur.com/6ac4u.jpg) Any idea/ thought to put me on the right track is much appreciated.


title: Plot the nodes of the shortest path found by dijkstra algorithm

text: I am using qgis and pgrouting plugin to find the shortest path between two nodes and display the route. Currently, the plot shows the shortest path. However, I also want to know which nodes the path passes through. So is it possible to actually get this information?


title: QGIS 2.4 plugin PyQT qWebview javascript interface problem

text: I have a QGIS plugin which contains a `qtwebview` and Javascript communication between the page and the plugin. It is working perfectly with QGIS 2.2 but not with 2.4 (both QGIS versions share the same plugin directory), The javascript side never gets called with `evaluateJavaScript()` and the python side never gets invoked from Javascript either. What could explain that it is working in 2.2 but not in 2.4?               def on_loadFinished(self):          ...               def init_Clicked(self,iface):         self.ui.webViewResult.page().mainFrame().addToJavaScriptWindowObject("QGISInterface", self.moTAZJsInterface)         self.ui.webViewResult.loadFinished.connect(self.on_loadFinished)         self.ui.webViewResult.load( QUrl(config["TAZ_URL_INTERFACE"]))         self.ui.webViewResult.show();


title: How to use pgsql2shp when PostGIS extensions installed in a separate schema as geometry

text: I have the PostGIS extension (postgis functions and spatial_ref_sys table) added to a specific schema called postgis. I'm wanting to run the pgsql2shp utility with a select query on geometries stored in other schemas that do not have postgis functions available within their respective schemas. Normally from psql I would first set the search path to contain both the postgis schema and the schema with the tables I'm querying on, but I'm having trouble figuring out how to do the same thing with the pgsql2shp utility. How would I go about setting a search path or having it look in multiple schemas so it picks up the postgis functions? Example Query (field_polygons is in different_schema):               pgsql2shp -f shapefile -h localhost -u user -P password gisdb "SET search_path TO postgis, different_schema; SELECT the_geom FROM field_polygons WHERE id = 1"      Results in error:               ERROR:  function postgis_version() does not exist     LINE 1: SELECT postgis_version()


title: How to query a PostGIS database through Geoserver and OpenLayers

text: I've built a web application using OpenLayers, which is connected to my PostGIS database through Geoserver. Now that I'm able to successfully display all of the layers, I'd like to add querying/interactivity capabilities to the map. For instance, I'd like to create a dropdown menu of checkboxes so the user could "Choose X who have done Y", another dropdown to "Choose X who haven't done Y", etc. I'm having trouble figuring out how exactly to add in this functionality. Does anybody know any good resources/examples for how I might go about doing this? Thanks for your help.


title: How to download hospitals and restaurants from OSM in shp format?

text: I need information on hospitals and restaurants in Guatemala. How can I download this information from OpenStreetMap? I want to use the data in Arcgis 10.1


title: Value of a Python programming course

text: I am considering upping my programming skills by taking a python geoprocessing course through ESRI. Before I part with the hefty registration fee, I would like to get your opinion on the value of 24 hrs of instructor led python programming. I am primarily interested in learning the following:   * For those of you who have taken a course such as this, did you walk away with applicable skills or more theoretical knowledge?   * For those of you who have recently been in the job market, did you find employers looked favorably on having a class such as this on your resume?   * For employers, does a specialty course like this really matter, or do most GIS analysts/specialist learn python programming on their own? Thanks for the advice.


title: Mapping zipcode to legislative district

text: Did anyone work with SunLight foundation data to map the zipcodes with legislative district data. I have to work with Illinois zipcodes only.


title: Do Different ArcGIS for Desktop License Levels Process Faster?

text: At my job, I have access to all three levels of ArcGIS licensing on my computer. My question is: Does the license you are using in ArcMap make a difference on the speed of ArcMap? I'm pretty sure that the drawing of the map isn't affected, based on tests I've done myself. But, it really does seem like some of the ArcToolbox tools do work faster in ArcInfo. I'm just trying to either validate my suspicious or prove I am going a bit crazy. Thanks for any replies!


title: QGIS systematically guesses wrong CRS

text: QGIS (1.7) systematically applies the wrong CRS when adding layers (or shapefiles from a given lot). They should be NAD83__MTM9 (EPSG:32189), as can be inferred from the prj-file, but the properties/general/options field show NAD27(CGQ77)/SCPQ_zone9. (EPSG:2015) Granted from there, I can specify/apply the correct CRS and the layers will line up. But it is a drag to do, and redo again, for many shapefiles. Is there anyway to correct that behavior or a better approach? Thanks for any guidance, cheers! **Update:** indeed corrected with the installation of 1.8, thanks to all.


title: Calculate distortion outside a UTM zone

text: One of my colleagues is working with data that is spread across two UTM zone. The majority of the data is in one zone, with a few outliers in another zone. He would like to know what the area distortion of those outliers would be if they were in the main UTM zone. Is there a formula to calculate the aerial distortion knowing how far into the other UTM zone the features were?


title: How to clean gap between two polygon?

text: I have a problem two layers, parcel and polygon, when I do create new building along side parcel it produces a gaps between parcel and building along the side, how can I clean this gab or how to avoid it. thank you.


title: How to create WKT file readable by QGIS 2.4.0

text: I am trying to create a well-known text (wkt) file that can be imported into QGIS ver. 2.4.0 Here's a minimal working example of the geometry I wish to represent.               POLYGON((-125 52, -125 21, -98 21, -98 52, -125 52),(-120 50, -120 35, -100 35, -100 50, -120 50))      I have created a text file named polygon.wkt containing only the above wkt representation of a polygon. When I try to import it into QGIS (Layer/Add Delimited Text Layer) the import dialog responds "no data found in file." How do I need to modify the file (a header, perhaps) to make it readable by QGIS? Solution: I needed to think of the file as delimited text with a header line. For some reason, I needed more than a single field to make it work. See the code block below for the working solution to the minimal example I provided.               name, geom     "my polygon", "POLYGON((-125 52, -125 21, -98 21, -98 52, -125 52),(-120 50, -120 35, -100 35, -100 50, -120 50))"


title: Adding (ir?)regular vector points

text: I run 2.0 qgis. (I am quite the rookie) I would like to add a point layer which is regularly spaced every 5 metres along the x-axis, and regularly spaced every 1m along the y-axis. I cannot seem to get this to work. I either get a 1x1m point layer or 5x5m point layer, but not the 5x1m required. ![Visual image of what I need:](http://i.stack.imgur.com/kAp8R.png) Any ideas much appreciated? Thank you!!!


title: Compute number of days with snow

text: I have a project and I want to map the effective number of days with snow cover for a year in some mountain area. Until now, I have processed the MOD10A2 (8 days snow cover) product derived from MODIS and I want to count the occurrence of 200 pixel value (this value correspond to the snow extent). I have tried to aggregate raster into the temporal gis offered by GRASS but the results are unsatisfactory. Any ideas? I can handle ArcInfo, GRASS GIS and R pretty well. Thanks!


title: Update OSM data already imported in qGIS

text: qGIS makes it really straightforward to import data from Open Street Map (there's a good tutorial here). However, despite this process is relatively straightforward, it takes some time to download the data, convert it, select the desired attributes to import, and then make queries to save a new layer. Since I am working on a project that will take some time to develop, I am concerned about having to repeat steps or losing data if I have to make new queries from OSM to add new data. Do you know if there's any way to automate this process as much as possible. Let's imagine that I have already downloaded OSM data, created several layers (with their corresponding attributes) and I spend some more time doing started doing some extra work (for example a heatmap). Now Imagine that I want to do the same with a more up-to-date version of OSM which will have more features. If I am not wrong I will have to repeat everything and it will be like starting from scratch. Do you know any way to just update all the layers without having to delete them and start over again? (I was wondering if there's any way to automate the process or create "smart" layers which query an specific layer or any way to update the original database)


title: Calculations in if statement

text: I want to perform a simple calculation with the Field Calculator using this piece of script (note that all the elif statements are attempts for another syntax and non of the worked):               def change(OldFor_Com):         if (OldFor_Com) == "OHWH":             return shape.area@Hectares / 40         elif (OldFor_Com) == "OHWH":             return !shape.area@Hectares / 30!         elif (OldFor_Com) == "OHWH":             return (!shape.area@Hectares / 30!)         elif (OldFor_Com) == "OHWH":             return !shape.area@Hectares! / 30         elif (OldFor_Com) == "OHWH":             return (!shape.area@Hectares! / 30)         if (OldFor_Com == "OTHH"):             return !shape.area@hectares! / 40      If the field is populated with `OTHH` I want another field to be populated with the area of the field, divided by 40, or 30, or... and so on. I modified my syntax a couple of times (as you can see) but it still says `Invalid Syntax`. Is a calculation possible with these if statement or do I need a workaround? **Note** that the shape.area@Hectares is an expression to calculate the area of the field, and not another field in my attributetable! ![screen1](http://s1.postimg.org/j9ry3wcsv/Untitled.png) ![screen2](http://s2.postimg.org/y4xdaf57t/Untitled.png) ![screen3](http://s17.postimg.org/y0ked1gu7/Untitled.png) ![screen3](http://s27.postimg.org/y2hnz4xo3/Untitled.png) **EDIT:** 7/15/2014 around 11pm I modified the script. I added a field to the table and calculate the area of the shape first (with !shape.area@Hectares!) and then I try to run this               def change(OldFor_Com):         if OldFor_Com == "OTHH":             return Area_Ha / 40         if OldFor_Com == "OHWH":             return Area_Ha / 30         if OldFor_Com == "OSFH":             return Area_Ha / 375         if OldFor_Com == "PINE":             return Area_Ha / 10         if OldFor_Com == "OOFH":             return Area_Ha / 375         if OldFor_Com == "None":             return Area_Ha         else:             return 0      It bounces back with the following:               ERROR 000539: Error running expression: change(u"OTHH", 7)      Traceback (most recent call last):     File "<expression>", line 1, in <module>     TypeError: change() takes exactly 1 argument (2 given)          Failed to execute (Calculate Field (10)).      Is it not possible to do any calculations in the if statement? **EDIT:** 7/15/14 12:30pm This how the attributestable looks like before I add the field `SizeScore` ![before](http://s29.postimg.org/h3rv6xmo7/Untitled.png) After I add the field `SizeScore` (FLOAT). ![after](http://s23.postimg.org/wlfa5r0sr/Untitled.png) And the following screen comes up if I use this code:               def change(OldFor_Com, Area_Ha):         if OldFor_Com == "OTHH":             x = Area_Ha / 40             return x         if OldFor_Com == "OHWH":             x =  Area_Ha / 30             return x         if OldFor_Com == "OOFH":             x = Area_Ha / 375             return x         if OldFor_Com == "OSFH":             x = Area_Ha / 375             return x         if OldFor_Com == "PINE":             x =  Area_Ha / 10             return x         if OldFor_Com == "None":             return Area_Ha         else:             return 0      ![screen](http://s30.postimg.org/lp0jzjvlt/Untitled.png) It only populates the rows where it does not need to calculate. The only statement that is working is that one.               if OldFor_Com == "None":          return Area_Ha      But for the others are no error messages coming back! **SOLVED, SEE ACCEPTED ANSWERE** Thanks a lot to everybody who went through this with me!


title: pyQGIS - Data Source is Invalid

text: Why I am getting `Data Source is Invalid` in python console of QGIS when running:               vlayer = QgsVectorLayer("home/sanjayr/Desktop/india-latest.shp/buildings.shp","sanjay","ogr")      I am using QGIS 2.4 on Ubuntu 12.04 LTS.


title: How to subsribe Web Map Service from GeoServer?

text: I have my `shapefile` inside `GeoServer` and I already published it as a `layer`. I am able to preview the layer under `GeoServer` web interface. Now, I want to subscribe / import/ retrieve the `shapefile` as image to my third party online map API. So that, the polygon can overlay on my third party map. The online map API I am using now is a non-popular API and they do not provide any sample or source for subscribing / importing / retrieving WMS layer. They just provide the following code to demonstrate the use of their map to import WMS layer:                       Xmap.ImportWMSLayer('layer title',             {             url:'http://WMS',             layers:'Layer_Name',             version: '1.3.0', // WMS version             crs:'EPSG:900913',              srs:'EPSG:900913',              format:'image/png' // format             });      I tried to key in the respective parameters but it just does not seems to retrieve. Hence, now I would like to diagnosis the problems               Problem 1: I set the WMS in GeoServer wrongly and hence can't retrieve the map          Problem 2: I set the WMS in GeoServer correctly but I key in the parameter name wrongly, and hence can't retrieve the map      So here is my question:               1) I am able to preview my layer in GeoServer web interface, does this mean, the layer is ready for retrieve to third party map API?          2) If answer for question (1) is a NO. Then how can I make sure my GeoServer is providing the map service for my layer?      **EDIT:** update based on the request of more info ![My Layer Name](http://i.stack.imgur.com/YNXQZ.jpg) Since the layer name in my `GeoServer` is `korea`, So I also key in '`korea`' as my `layer name` in as the parameter name. **EDIT:** Update the parameters information which I set                       map.ImportWMSLayer('korea',             {             url:'http://localhost:8080/geoserver/web/',             layers:'korea:korea',              version: '1.3.0',              crs:'EPSG:4326',              srs:'EPSG:4326',              format:'image/png'              });      **EDIT** : I update the url to `http://localhost:8080/geoserver/korea/wms` **EDIT** Add the preview layer display screen in my GeoServer Inspect Element in OpenLayer preview               <img id="OpenLayersDiv44" class="olTileImage" src="http://localhost:8080/geoserver/korea/wms?LAYERS=korea%3Akorea&STYLES=&FORMAT=image%2Fpng&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap&SRS=EPSG%3A4326&BBOX=121.54456007603,28.652528651739,136.26831368157,39.637829193375&WIDTH=512&HEIGHT=382" style="width: 512px; height: 382px; position: relative;">      ![preview Layer in GeoServer](http://i.stack.imgur.com/lxPbL.jpg)


title: How to calculate 90% confidence interval in a raster stack

text: I have a stack of grid files (23) in ASC format from a MODIS sensor and I need to calculate a 90% confidence interval, mean, min, etc, for each one in R, how can I do that? specially the interval.               ## set working directory     FilesPath <- "C:/Users/Grettel/ "          setwd(FilesPath)     FileList <- list.files()          # Loop for all grids by julian days and create stats     for (i in  unique(substring(FileList,40,42))){            print(as.character(paste(i,sep="")))            JulianDay<-as.character(paste(i,sep=""))            # Create output directory       WriteFolder<-paste(FilesPath,"StatsGrids",sep="")       dir.create(paste(WriteFolder))            # list files by pattern of julian day       FilesToRead<-list.files(pattern=paste("*",JulianDay,"*", sep=""))       print(FilesToRead)            # create raster stack       FilesStack<-stack(FilesToRead)              # mean files       meanRaster<-mean(FilesStack)          # min files       minRaster<-min(FilesStack)       # max files       maxRaster<-max(FilesStack)       # sd files       sdRaster <-calc(FilesStack, fun = sd)             # intervalo confianza       #CIRaster <- ci(meanRaster, conf.level = 0.9)            # write new asci into output folder       writeRaster(meanRaster,                    filename=paste(WriteFolder,                                  "/Modis_2000_2014_mean_JulDay_",                                  JulianDay, sep=""),                   overwrite=T, "ascii")      }


title: 3D Bulk Modelling tools in GIS

text: I work in an Architectural practice concentrating on large scale urban planning and Urban Design Schemes. Our current software is not really geared for this sort of work. We need to be able to accurately plot building block models onto accurate cadastral information. Is there any type of software that would allow us to do this sort of work?


title: Tool (Simplify Polygon) parameter is not visible anymore

text: I am having a strange problem in ArcGIS Toolbox. I am trying to use the Simplify Polygon (Cartography) tool from ArcToolbox. When I open the tool, there is no possibility to input Simplification Tolerance. Please see the attached image. I can run the tool as a batch and from the Python window without any problem. I am using ArcGIS 10(Editor) with Service Pack 5 on Windows 7 64Bit (English) and Windows Server 2008 64Bit (German). On both machines, I have the same problem. Any ideas? ![enter image description here](http://i.stack.imgur.com/u04Cm.jpg)


title: Problem overlapping shapefile and vector in R

text: I want to plot a shapefile over a raster file in R but I can't make them overlap perfectly: the raster appears to be rotated of few degrees counter- clockwise. Is it a problem with the projection? ![enter image description here](http://i.stack.imgur.com/RypsS.png) Please consider the following MWE               library(raster)     library(rgdal)           # Download from http://biogeo.ucdavis.edu/data/gadm2/shp/ITA_adm.zip     shape_file = "ITA_adm1.shp"     # Download from http://sedac.ciesin.columbia.edu/data/set/gpw-v3-population-density/data-download     # Setting Geography: Country, Italy; Data Attributes: Grid     pop_density_file ="w001001.adf"          italy_map <- readOGR(dsn = shape_file, layer = "ITA_adm1")     italy_map_dens <- raster(pop_density_file)          colPal <- colorRampPalette(c("white", "red"))( 500 )          par(mar=c(0,0,0,0))     plot(italy_map_dens, xlim = c(6.70, 18.32), ylim = c(35.2, 47.6), axes=FALSE, box=FALSE, legend=FALSE, col=colPal)     plot(italy_map, xlim = c(6.70, 18.32), ylim = c(35.2, 47.6), border="grey", add=TRUE)


title: Calculate surface elevation from DEM of pit and known volume with which it is being filled

text: I am using ArcGIS 10.1. I have a DEM of a pit. I have a known volume of material that I want to put in the pit. I want to know what the elevation of the flat surface of this material will be once it is all in the pit. There are a lot of tools for volume calculation that ArcGIS offers, but none seem to relate to this. What springs to mind is to iterate through a set of constant raster values and see which most closely matches my known volume when I use each constant raster with the DEM in the "cut/fill" tool. Any other ideas appreciated. Wil


title: openlayers clickevent for undercoat features polygon on multiple vector layers

text: I'm using two timers that give me coordinates and some data for my features (Sensors(de Sensoreinheit) & receiving unit (de Empfangseinheit)), deleting all features on the Layers and placing new ones. What I need is a `clickevent` that gives me all features under the `mouseclick` on specific layers (all would be okay too).                 var layerListeners = {         featureclick: function(e) {             console.log("click");             if (e.feature.data.type == 'empfangseinheit') {             } else if (e.feature.data.type == 'sensoreinheit')             {                 updateTabView(toFormParameter(e.feature.data.data.person));             }             return false;         },         featureover: function(e) {              }     };          var empfangseinheitenLayer = new OpenLayers.Layer.Vector("Empfangseinheiten", {         styleMap: empfangseinheintenStyle         , eventListeners: layerListeners     });     var sensoreinheitenLayer = new OpenLayers.Layer.Vector("Sensoreinheiten", {         styleMap: sensoreinheitenStyle         , eventListeners: layerListeners     });          map = new OpenLayers.Map("map", {         controls: [],         allOverlays: true,         units: "m",         scales: [100, 90, 80, 70, 60, 50, 40, 30, 20, 10]     });     map.addLayer(empfangseinheitenLayer);     map.addLayer(sensoreinheitenLayer);          var sensoreinheitenTimer = setInterval(handleSensoreinheiten, 1000);     var empfangseinheitenTimer = setInterval(handleEmpfangseinheiten, 1000);          function handleSensoreinheiten() {         $.ajax({url: "#{request.contextPath}/rest/Sensordaten/aktuelle",             cache: false,             dataType: 'json',             success: function(result) {                 sensoreinheitenLayer.removeFeatures(sensoreinheitenLayer.features);                 var features = [];                 $.each(result, function(i, iter) {                     var point = new OpenLayers.Geometry.Point(iter.sensordatum.x, iter.sensordatum.y);                     var circle = OpenLayers.Geometry.Polygon.createRegularPolygon(point, 0.25, 20, 0);                     var feature = new OpenLayers.Feature.Vector(circle);                     feature.data = {                         "type": "sensoreinheit",                         "data": iter                     };                     features.push(feature);                 });                 sensoreinheitenLayer.addFeatures(features);             }});     }         The problem is that the listener does not work.      ### **Additional Information** I am using Primefaces 5.0 on a TomEE.


title: Spatial Interpolation - defining a minimum value for agroup of points

text: I am doing spatial interpolation with a group of points. Some points have complete observation of events:   * for a period of 24 hours we have observed N events in each point. These points have N events. While others are incomplete:   * for a period of 6 hours we have observed M events in each point. This means that these points have >= M events, since we don't have 24h of observation we just know what is the value for the minimum amount of events. Although the second group of points have incomplete observations, they are still important for the interpolation. I would like to use them in the spatial interpolation (Krigging or co-Krigging) with the condition that these are not complete observations. I would like to interpolate new points and use both complete (24h) and incomplete (6h) points. How can I take into consideration in a spatial interpolation points that are incomplete? How can this be done?


title: How to integrate DEM in Kriging and other interpolation methods?

text: I would like to see the difference of including DEM in interpolation from without. Some studies on climate data interpolation suggest that interpolation performs better with DEM. How do I do this in ArcGIS 10.2? Thanks.


title: How to set up an offline maps viewer application on linux?

text: I'm having an 8-years old notebook (pentium 4 3Ghz 1Gb Ram) with Debian (installed recently), I just want to be able to navigate through the map around 1-2 cities, with ability to seek adresses and lay bike, car and public transport courses offline, preferably using openstreetmaps. I have already managed to install marble and download some rastered maps, but I can nither to seek adresses nor to lay courses offline. I don't even see house numbers offline. Maybe there are some instructions how to acomplish that? I have also tried to install other software - Gosmore and Navit but didn't manage to make work either of them (I couldn't see a map), so I gave up with this soft and came up with Marble. Any ideas appreciated. **Updated:** By the way, it seems that marble can not make public transport routes and monav can't seek for house numbers. And I didn't succeed to use an adress search online function in marble as well - I have no idea where to get and how to install such data (preinstalled only about 12000 cities as I read). I'm totally confused.


title: How to calculate the average distance among set of points as measure of closeness

text: I have three series (or tables) of points in PostGIS which I want to compare. I would like to know how "close" each series is when compared to the other two. In other words I would need some measure such as the average of the average distance of each point with all the others (that is, the average of the distance matrix)... Sounds complicated, but I am sure some common statistic exist to compare series of points. How my query would look like?


title: Creation point layer from linear layer with only ends of lines

text: I have ArcGIS 10 with ArcView licence and I need to make point layer from a linear layer, with only start and end points for each line. In ArcInfo I could use Feature Vertices To Points with "both_ends" as option but I don`t have such licence. I found "Polyline Vertices to Points" script for creation point layer from linear layer with only ArcView licence, but it returns all vertices not only ends. Please share your suggestions about how can I delete all vertices from lines but ends or maybe other variant how can I make point layer from line layer with only ends of lines


I have two feature classes, one ploylines, one a buffer. How can I populate a field of the polylines to indicate if they intersect?


I have my roads and locations added into my PostGis DB, and have installed pgRouting's SQL files. However, when I attempt to run:               SELECT * FROM shortest_path('        SELECT gid AS id,                start_id::int4 AS source,                end_id::int4 AS target,                shape_leng::float8 AS cost        FROM myr_roads_point',     1,     5110,     false,     false);      It gives me the following error:               ERROR:  function shortest_path(unknown, integer, integer, boolean, boolean) does not exist     LINE 1: SELECT * FROM shortest_path('                           ^     HINT:  No function matches the given name and argument types. You might need to add explicit type casts.          ********** Error **********          ERROR: function shortest_path(unknown, integer, integer, boolean, boolean) does not exist     SQL state: 42883     Hint: No function matches the given name and argument types. You might need to add explicit type casts.     Character: 15      Any idea?


title: Calculating from a location to another location (calculation must be on a line/road path)

text: I'm using QGIS. I'm new. I have a CSV file, with the source point (the_geom) and destination point (the_geom). I have a road network already imported into PostGIS with pgRouting, and I need to have the ability to calculate the distance from a source to a destination, by using the road path as the road. After calculating, it should save the distance back into the DB in a new column. Is there any quick and easy way to solve this? Thanks for all the help! :)


title: Query Layer not found ERROR 000229

text: My code executes in three steps:   1. create FeatureLayer,   2. create QueryLayer and   3. join QueryLayer to FeatureLayer:              arcpy.MakeFeatureLayer_management(**{         "in_features":features.get_dataset_fullpath(),         "out_layer": features.get_layer_name(),     })          arcpy.MakeQueryLayer_management(**{         "input_database":   features.get_attrib_QLD().get_workspace(),         "out_layer_name":   features.get_attrib_QLD().get_layer_name(),         "query"         :   features.get_attrib_QLD().get_spatial_query_definition(),         "oid_fields"    :   features.get_attrib_QLD().get_OID_field()     })          arcpy.AddJoin_management (**{          "in_layer_or_view" : features.get_layer_name(),          "in_field"         : 'RDUWI',          "join_table"       : features.get_attrib_QLD().get_layer_name(),          "join_field"       : 'RDUWI',          "join_type"        : 'KEEP_ALL'     })      However, while executing the step three (join), the code returns: arcgisscripting.ExecuteError: ERROR 000229: Cannot open permits_querylayer Failed to execute (AddJoin), even though the layer in question has been created in step 2. FYI: Step two is time consuming however, I thought the execution of the next statement would wait until the previous one has returned?


title: QGIS cannot open output layers

text: I have problem with QGIS 2.0.1. Everytime I run tool from SAGA or GRASS, I get exactly the same error, for example:               Oooops! The following output layers could not be open     Channel Network: /tmp/processing/7efeeab66b224094b23a2701a4fe30ca/CHNLNTWRK.tif     Channel Direction: /tmp/processing/34adb5dd62cb4273981f023aab6b6ded/CHNLROUTE.tif     Channel Network: /tmp/processing/eadce7c1d82e47f1a50dc1c91fd53125/SHAPES.shp     The above files could not be opened, which probably indicates that they were not correctly produced by the executed algorithm     Checking the log information might help you see why those layers were not created as expected     This algorithm requires SAGA to be run. A test to check if SAGA is correctly installed and configured in your system has been performed, with the following result:     SAGA seems to be correctly installed and configured      Im running QGIS on Linux Mint 64 bit, Windows Vista 32 bit and Windows 7 32 bit, but the problem is exactly the same (even with QGIS 1.8). I have unchecked "Enable SAGA 2.0.8. compatibility" and set the right route. I have tried completly reinstall everything (including cleaning libraries) and I've tried to install QGIS on brand new installation of Linux Mint with same result. It seems, that compute is running (full load of CPU core for a while - the bigger the layers are, the longer it takes till it stop with error) and log is running, but after stucking I can only close whole window. Temporary folders are created, but there is nothing in then. I have tried to set "hard" route for the result, but it didnt work either. For me it looks like everything is right except writing down an output, but probably I'm wrong. I've found plenty of similar problems solved by unchecking "Enable SAGA compatibility" or reinstalling SAGA/QGIS, but nothing worked for me. Can anybody help me? * * * EDIT 03/04/2014 21:24 UTC               Během vykonávání Python kódu se objevila chyba:          Traceback (most recent call last):       File "/usr/share/qgis/python/plugins/processing/ProcessingPlugin.py", line 142, in openHistory         dlg = HistoryDialog()       File "/usr/share/qgis/python/plugins/processing/gui/HistoryDialog.py", line 57, in __init__         self.fillTree()       File "/usr/share/qgis/python/plugins/processing/gui/HistoryDialog.py", line 65, in fillTree         elements = ProcessingLog.getLogEntries()       File "/usr/share/qgis/python/plugins/processing/core/ProcessingLog.py", line 90, in getLogEntries         lines = tail(ProcessingLog.logFilename())       File "/usr/share/qgis/python/plugins/processing/core/ProcessingLog.py", line 256, in tail         return Tailer(file).tail(lines)       File "/usr/share/qgis/python/plugins/processing/core/ProcessingLog.py", line 240, in tail         if not self.seek_line():       File "/usr/share/qgis/python/plugins/processing/core/ProcessingLog.py", line 201, in seek_line         bytes_read, read_str = self.read(read_size)       File "/usr/share/qgis/python/plugins/processing/core/ProcessingLog.py", line 179, in read         read_str = self.file.read(read_size)       File "/usr/lib/python2.7/codecs.py", line 671, in read         return self.reader.read(size)       File "/usr/lib/python2.7/codecs.py", line 477, in read         newchars, decodedbytes = self.decode(data, self.errors)     UnicodeDecodeError: 'utf8' codec can't decode byte 0x9a in position 0: invalid start byte          Verze Pythonu:     2.7.4 (default, Sep 26 2013, 03:23:21)      [GCC 4.7.3]               QGIS verze:     2.0.1-Dufour Dufour, exported          Python cesta: ['/usr/share/qgis/python/plugins/processing', '/usr/share/qgis/python', '/home/hanbanaroda/.qgis2/python', '/home/hanbanaroda/.qgis2/python/plugins', '/usr/share/qgis/python/plugins', '/usr/lib/python2.7', '/usr/lib/python2.7/plat-x86_64-linux-gnu', '/usr/lib/python2.7/lib-tk', '/usr/lib/python2.7/lib-old', '/usr/lib/python2.7/lib-dynload', '/usr/local/lib/python2.7/dist-packages', '/usr/lib/python2.7/dist-packages', '/usr/lib/python2.7/dist-packages/PILcompat', '/usr/lib/python2.7/dist-packages/gtk-2.0', '/usr/lib/pymodules/python2.7', '/usr/lib/python2.7/dist-packages/ubuntu-sso-client', '/usr/lib/python2.7/dist-packages/wx-2.8-gtk2-unicode', '/usr/share/qgis/python/plugins/fTools/tools', '/media/hanbanaroda/A80C5D840C5D4F06/\xc5\xa0kola/Ostatn\xc3\xad/QGIS']      "Během vykonávání Python kódu se objevila chyba:" = During processing Python code an error appeared, "verze" = version


title: Alternative to Farm Service Agency (FSA) CLU dataset?

text: Is there an alternative dataset that can I find that is similar to the Farm Service Agency (FSA) Common Land Unit (CLU) ?


title: Can I verify that attribute indexes are being used in a query from ArcObjects?

text: We have some software as part of a larger project written using ArcObjects via .NET. This software migrates and merges data from multiple file geodatabases in one schema to a single file geodatabase in a similar but different schema. Some steps in the migration involve looking up an existing row and updating it with data from a different table. This can get very slow when the destination geodatabase is large, so I have added some multi-field attribute indexes to improve the query/update performance. Is there a way that I can confirm that the queries are using the new attribute indexes, and perhaps how the query is being executed? I guess I'm really looking for the equivalent of a query plan - something similar to `EXPLAIN PLAN` in Oracle. The queries are on a file geodatabase using ArcObjects directly (no RDBMS backend or SQL).


title: How to fix problems opening JPEG2000 in QGIS Mac version 1.8?

text: I am running QGIS v. 1.8 on both my Windows and Mac partitions; Every time I try to open a jpeg2000 (jp2) raster on the mac version the program locks up and I need to force quit. Does anyone have any ideas how to fix this problem?


title: From RGB raster to grayscale

text: I have been searching for a way to rasterize shp files (Unexpected output when using gdal_rasterize) and it worked with a little help. Now, my code is first rasterizing the largest sip file and than using those dimensions for further layers. For each layer I am inputting a value (COST) which is used as value for each pixel but it is using a RGB format but I actually only need the files to be in grayscale. The RGB is giving a problem when I enter values higher than 255 => than it subtracts n*255 until the value is in the range of 0-255. Can anybody help me with changing this code to the one that just uses grayscale and where the cost values will be retained like I want it to?               import os     from osgeo import gdal, ogr     import config3          RASTERIZE_COLOR_FIELD = "__color__"     xmi = 0     xma = 0     ymi = 0     yma = 0     syspath = str()     newmap = str()     cost = 0     outputfilep = str()     source_srs = 0          def one():         while True:             try:                 global syspath                 syspath = str(input("Systempath: "))                 break             except (NameError, SyntaxError, ValueError, TypeError):                 print ""                 print("This is not a string. Please do not forget to type your input between quotation marks.")                 print ""             except (UnicodeEncodeError):                 print("Please do not use special characters in the systempath.")                 print ""          def three():         while True:             try:                 global cost                 cost = int(input("Cost per squared meter for this type of land-use: "))                 break             except (NameError, ValueError, SyntaxError, UnicodeEncodeError):                 print ""                 print("This is not a numerical value, please re-enter.")                 print ""          def rasterize(pixel_size=1):              one()              (shapefilefilepath, shapefilename) = os.path.split(syspath)         (shapefileshortname, extension) = os.path.splitext(shapefilename)               outmap = shapefilefilepath + "/"              three()                orig_data_source = ogr.Open(syspath)         source_ds = ogr.GetDriverByName("Memory").CopyDataSource(orig_data_source, "")         source_layer = source_ds.GetLayer(0)         global source_srs         source_srs = source_layer.GetSpatialRef()         x_min, x_max, y_min, y_max = source_layer.GetExtent()         global xmi         xmi = x_min         global ymi         ymi = y_min         global xma         xma = x_max         global yma         yma = y_max         field_def = ogr.FieldDefn(RASTERIZE_COLOR_FIELD, ogr.OFTReal)         source_layer.CreateField(field_def)         source_layer_def = source_layer.GetLayerDefn()         field_index = source_layer_def.GetFieldIndex(RASTERIZE_COLOR_FIELD)         for feature in source_layer:             feature.SetField(field_index, cost)             source_layer.SetFeature(feature)         x_res = int((x_max - x_min) / pixel_size)         y_res = int((y_max - y_min) / pixel_size)         global outputfilep         outputfilep = outmap + shapefileshortname + ".tif"              target_ds = gdal.GetDriverByName('GTiff').Create(outputfilep, x_res,                 y_res, 3, gdal.GDT_Byte)         target_ds.SetGeoTransform((                 x_min, pixel_size, 0,                 y_max, 0, -pixel_size,             ))         if source_srs:             target_ds.SetProjection(source_srs.ExportToWkt())         else:             target_ds.SetProjection('LOCAL_CS["arbitrary"]')         err = gdal.RasterizeLayer(target_ds, (1,2,3), source_layer,                 burn_values=(0,0,0),                 options=["ATTRIBUTE=%s" % RASTERIZE_COLOR_FIELD])          def rasterize2(pixel_size=1):         continue1 = 0              while continue1 != "STOP":                  while True:                 try:                     continue1 = str(input("Systempath: "))                     break                 except (NameError, SyntaxError, ValueError, TypeError):                     print ""                     print("This is not a string. Please do not forget to type your input between quotation marks.")                     print ""                 except (UnicodeEncodeError):                     print("Please do not use special characters in the systempath.")                     print ""                  while True:                 try:                     cost = int(input("Cost per squared meter for this type of land-use: "))                     break                 except (NameError, ValueError, SyntaxError, UnicodeEncodeError):                     print ""                     print("This is not a numerical value, please re-enter.")                     print ""                  (shapefilefilepath, shapefilename) = os.path.split(continue1)             (shapefileshortname, extension) = os.path.splitext(shapefilename)             outmap = shapefilefilepath + "/"                       outputfilep = outmap + shapefileshortname + ".tif"                  orig_data_source = ogr.Open(continue1)             source_ds = ogr.GetDriverByName("Memory").CopyDataSource(orig_data_source, "")             source_layer = source_ds.GetLayer(0)             field_def = ogr.FieldDefn(RASTERIZE_COLOR_FIELD, ogr.OFTReal)             source_layer.CreateField(field_def)             source_layer_def = source_layer.GetLayerDefn()             field_index = source_layer_def.GetFieldIndex(RASTERIZE_COLOR_FIELD)             # Generate random values for the color field (it's here that the value             # of the attribute should be used, but you get the idea)             for feature in source_layer:                 feature.SetField(field_index, cost) # DIT IS AANGEPAST                 source_layer.SetFeature(feature)             x_res = int((xma - xmi) / pixel_size)             y_res = int((yma - ymi) / pixel_size)             target_ds = gdal.GetDriverByName('GTiff').Create(outputfilep, x_res,                     y_res, 3, gdal.GDT_Byte)             target_ds.SetGeoTransform((                     xmi, pixel_size, 0,                     yma, 0, -pixel_size,                 ))             if source_srs:                 target_ds.SetProjection(source_srs.ExportToWkt())             else:                 target_ds.SetProjection('LOCAL_CS["arbitrary"]')             err = gdal.RasterizeLayer(target_ds, (1,2,3), source_layer, burn_values=(0,0,0), options=["ATTRIBUTE=%s" % RASTERIZE_COLOR_FIELD])                  continue1 = input("Press enter if you want to continue inputting shapefiles. Otherwise type STOP: ")              if err != 0:             raise Exception("error rasterizing layer: %s" % err)                     if __name__ == '__main__':         rasterize()         rasterize2()


title: How can I put coded domains into a list to be created in python script?

text: I have a large list of coded domains. They are numbers like: 01.6530.01, 05.7991.08, 05.7991.01, 05.7991.06 The code will be the exact same as the description. I've tried to make a list using:               FRClist = [01.6530.**01**, 05.7991.08, 05.7991.01, 05.7991.06]      but it gave me a syntax error on the 01 that is bolded. So I tried to change all the values to a string, using:               FRCstr = str(01.6530.**01**, 05.7991.08, 05.7991.01, 05.7991.06)      but it still gave me a syntax error on the same bolded 01. In the end I hope to put the list into this:               arcpy.AddCodedValueToDomain_management(GDBpath, "FRC", FRClist[x], FRClistp[x] )      I have quite a few domains which require a lot of coded values, so any help is appreciated.


title: Can FileGDB SDK be made to support custom coordinate systems?

text: In May 2012, a question was asked about File Geodatabase (*.gdb) support in QGIS? In the comments there, it was pointed out by @MattWilkie that the FileGDB SDK does not support custom coordinate systems. Has there been a workaround developed for custom coordinate systems?


title: PyQGIS: Editable Combobox for Value Map widget isses

text: I am having issues with an editable combo box in QGIS 2.4. The edit widget is set to Value Map with an uploaded CSV. However, if a new values (not included in the key value pair list) is added, the attribute value is not saved, thus defeating the purpose of the editable combo box. Currently, I have it set editable like this in my python form (set under python init in layer properties).               Attributefield.setEditable(True)      I saw a checkbox "Editable" option for value map as well, but nothing happens when it is checked (combo box is not editable). My original goals are mentioned in this post, using setLineEdit PyQt method. But I decided to try this editable combo box option, which still does not work.


title: Mask Bridge with Road

text: Is there a tool to create mask between bridge line/point on road? The requirement to have gaps 0.2mm between road and bridge and the offset fill with white color. For now, I add a white line attach to bridge symbol but the result not that great. Looking for improvement and better result. Thanks! ![enter image description here](http://i.stack.imgur.com/gMnHU.png) ![enter image description here](http://i.stack.imgur.com/6VW7h.png) ![enter image description here](http://i.stack.imgur.com/CtmQK.png)


title: QGIS: Need to remove unwanted white bands when raster layer are overlapping

text: In my project when raster layers are overlapping a white band appears the edge of the top layer. How can I remove it and obtain a seamless overlapping of the raster ? (I have try to built a vrt with those layers but the with bands remained) I am working with QGIS 2.4 _Here is a example of the phenomena (each layer has 10% transparency)_ ![white band on overlapping raster](http://i.stack.imgur.com/xCOZT.png) _Here is the gdalinfo output of one of the raster_ ` Driver: GTiff/GeoTIFF Files: /vdata/scot_009h_georef.tif Size is 5985, 2349 Coordinate System is: GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433], AUTHORITY["EPSG","4326"]] Origin = (-3.873394354081865,58.119237085214102) Pixel Size = (0.000601740178426,-0.000601740178426) Metadata: AREA_OR_POINT=Area Image Structure Metadata: COMPRESSION=LZW INTERLEAVE=PIXEL Corner Coordinates: Upper Left ( -3.8733944, 58.1192371) ( 3d52'24.22"W, 58d 7' 9.25"N) Lower Left ( -3.8733944, 56.7057494) ( 3d52'24.22"W, 56d42'20.70"N) Upper Right ( -0.2719794, 58.1192371) ( 0d16'19.13"W, 58d 7' 9.25"N) Lower Right ( -0.2719794, 56.7057494) ( 0d16'19.13"W, 56d42'20.70"N) Center ( -2.0726869, 57.4124932) ( 2d 4'21.67"W, 57d24'44.98"N) Band 1 Block=5985x1 Type=Byte, ColorInterp=Red NoData Value=0 Band 2 Block=5985x1 Type=Byte, ColorInterp=Green NoData Value=0 Band 3 Block=5985x1 Type=Byte, ColorInterp=Blue NoData Value=0 `


title: How to hide overlays in Leaflet layerControl if geojson var not found?

text: i would like to perform a specific task: _i have a set of geojson layers, loading features from local files. i have the layercontrol on the map showing available layers to show/hide_ now, regarding various set of local files, i would like to show ONLY the available data (means "if geojson.js file found). By this way, i can add/remove set of geojson files and see the layercontrol showing them accordingly. i suppose a If statement...but which one?


title: Displaying large number of polygons on web map?

text: I've started working on a project with a huge number of polygons, and I'm unsure of how to display/host it. The map I'd like to end up with will have three layers: 1) Shaded polygons representing each of the 70K+ U.S. census tracts. I currently have this as both a shapefile and a spqllite file in Tilemill. 2) Another 15K points representing locations on roadways (this layer will be interactive) 3) A street layer, such as OpenStreetMap I would set the opacity of the first layer so that the street layer remains visible. I've looked at possibly using UTFGrids, but haven't found much documentation. Any ideas on how I might want to implement this one? UPDATE: I'd prefer to use TileMill, but am open to different platforms.


title: ArcPy add point to a FeatureServer Layer

text: What is the best practice to add a point to a FeatureServer Layer using arcpy? using rest services, write to db, etc? Where can I find a code sample?


title: Why does arcpy.SelectLayerByLocation_management give ERROR 000368: Invalid input data?

text: When i trying to do Select by location i got : arcgis scripting. ExecuteError: Failed to execute. Parameters are not valid. ERROR 000368: Invalid input data. Failed to execute (SelectLayerByLocation) I am trying to excecute select by location on current selected feature after select by attribute.               env.workspace = "C:\temp"          input = 'C:\\temp\M0TASHT\M0TASHT.gdb\sLine'     input1= 'C:\\temp\M0TASHT\M0TASHT.gdb\sManhole'     arcpy.MakeFeatureLayer_management('sManhole', 'sManhole_lyr')     in_features = 'sManhole_lyr'     where_clause = '"ManholeSubtype" = ('"4515" ')'     arcpy.SelectLayerByAttribute_management(in_features, "NEW_SELECTION", where_clause)     arcpy.SelectLayerByLocation_management("C:\\temp\M0TASHT\M0TASHT.gdb\sLine", "INTERSECT", "sManhole_lyr", "", "NEW_SELECTION")


title: Are there templates for adding a search filter on a CartoDB embed?

text: I'm looking for any examples of an embedded CartoDB map with a search filter built in. So, a map of schools, say, with a pulldown that will let you select one school and then zoom to it. Or a map of schools and a pull down that will let you toggle between Elementary, Middle and High schools. I know I can use filters or a custom query but I'm looking for something viewers of the final map can play with themselves. Here's the deal: I'm trying to figure out what it would take to make a map that a colleague wants. If I have to do a lot of JS work, that's a bigger project. If there's a cut-and-paste solution that makes it easier. And if there's a GUI, they can make it themselves.


title: GMT raster (*.tiff + *.pgw) in epsg4326 shifts few kms (4326-900913-3857 problem?)

text: This seems to be an old problem with an already known solution (ie. "...stop using 900913 and start using 3857..."), but not sure whether this problem is about software versions, the raster itself or just projections. The thing is that I did a quite big raster (TIFF image + its PGW projection text file based on lon/lat values) with the software GMT (Generic Mapping Tools). It covers a big area (range: -83ºW, -75ºW, -17ºS,-5ºS). Its center is around -79ºS and -11ºS. To produce it, GMT requires both the lon/lat ranges and the projection. For the former I used the previous lon/lat range and for the latter I used a Mercator projection and a arbitrary width (in GMT language: -JM10i). So, the output raster was beautifully produced (as usual), and I loaded into my HTML with OL (of course) and served by Mapserver (as WMS). It overlays google, bing and esri base maps (so in MAP file I added the corresponding EPSG:4326 and 900913 as usual, and the latter was added to the epsg file in the /usr/share/proj directory, plus the usual "projection" and "displayProjection" calls in OL). So, it loaded normally as the other ones served by Mapserver (I have some more TIFF like this one, ie. GMT raster, but much more small ones, ie. less area). However, I noticed that around -80ºS and -8ºS, it should line up with other WMS layers (eg. seismic lines) but it didn't (it has an offset of ~3 km to the north in regard to other WMS layers), so I started changing things (900913 by 3857) in my MAP file and the HTML but nothing, the problem didn't go away. Interestingly, I did some tests and merged all the TIFFs and got only a unique TIFF (ie. new range: -83 W,-70 W,-19 S,-5 S), and this one produced an offset of ~10 km (instead of 3 km at the same location). Also, I loaded all the WMS layers and used the unique TIFF as base layer (because all the WMS layers served by Mapserver are in epsg4326), but the erroneous offset is still there. So, I was wondering if this problem could be caused by the OL 2.11, or Mapserver 5.6.7 and Proj 4.8.0 (these two are RPMs from ELGIS for RHEL5) I'm using. Another possibility is that the -JM10i used in GMT to generate the TIFF was the problem, or it is because the TIFF area (in lon/lat) covers two UTM zones (UTM zone 17 South and UTM zone 18 South). I think the "EPSG:900913 or EPSG:3857" thing is not the problem. After testing both in Mapserver and OL, I checked in Firebug and saw that the request (SRS=EPSG) were the same for both cases, ie. when 900913 was used it appeared SRS=EPSG=900913, and when 3857 was used: SRS=EPSG=3857. I also read that using either one of them it should work all in the same way (even though this may imply the other way: http://alastaira.wordpress.com/2011/01/23/the-google- maps-bing-maps-spherical-mercator-projection/). In my case I prefer using 900913 because google v3 uses it as far as I know. I hope someone could enlighten me here and provide me a solution, Thanks in advance, Gery


title: QGIS 2.4.0 Chugiak "reduced" real-type numbers. Precision problem? Bug?

text: I have a SHP point file created from CSV with QGIS 2.4.0 Chugiak. It contains some data on noise in a certain point. When I import from CSV, it seems to be fine, but when I export to SHP and open the attribute table - the real-type numbers are wrong ("reduced" number). The funny thing is, that when I open this SHP with older QGIS NextGIS 1.8 - everything is fine. Is it some kind of bug? I attach the printscreen: ![Left is Nextgis, right QGIS Chugiak. The suspicious fields are marked red.](http://i.stack.imgur.com/c9ws4.png)


title: PostGIS change column type varchar => geometry

text: I have a Postgres database with PostGIS extension on it. Also, I have a table with column lat with type `(varchar 40)`. Column values are like this:               52.1432514     52.1431855     52.1431586      I want to convert the whole column to `geometry` type. When I do `ALTER TABLE m_copy ALTER COLUMN latitude TYPE geometry USING (latitude::geometry);` I see an error:               ERROR: parse error - invalid geometry     SQL state: XX000     Hint: "52.1432514" <-- parse error at position 11 within geometry      What am I doing wrong and what is the right query? P.S Also I have a lon column. Maybe it is better to combine that columns?


title: interpreting the negative and positive values of profile curvature map in GIS

text: I would like to have help in interpreting the values returned by ArcGIS for profile curvature. By reading this Esri's webpage (http://blogs.esri.com/esri/arcgis/2010/10/27/understanding-curvature- rasters/), it seems that negative values should correspond to profiles upwardly convex (fig. 16.20A). Now, it seems counter intuitive to me, since I would imagine that negative values would best convey a sense of "concavity". I was just wondering if that reference is correct or if some error (typo?) would account for that seemingly oddity.


title: The best way to convert points to polygons using ArcObjects

text: I have a set of points from GPS I want to convert them to polygons..I have two approaches to do this.1- convert points to polygon..I do this with below code:                       private IPolygon PolygonFromPoint()         {                  try             {                 IGeometryBridge2 pGeoBrg = new GeometryEnvironment() as IGeometryBridge2;                 IPointCollection4 pPointColl = new PolygonClass();                      int numPoints = pList.Count;                 WKSPoint[] aWKSPointBuffer = new WKSPoint[numPoints];                 for (int i = 0; i < pList.Count; i++)                 {                     WKSPoint A = new WKSPoint();                     A.X = pList[i].X;                     A.Y = pList[i].Y;                     aWKSPointBuffer[i] = A;                 }                 pGeoBrg.SetWKSPoints(pPointColl, ref aWKSPointBuffer);                 IPolygon pPointColl2 = pPointColl as IPolygon;                 pPointColl2.SpatialReference = spatialrefrence;                 OutputPolygon = pPointColl2;                 status = true;                 return pPointColl2;             }             catch (Exception e)             {                 MessageBox.Show("PolygonFromPoint: " + e.Message);                 status = false;                 return null;             }         }      This works some times but some times I get the error that says the polygon is not colosed.I could not manage to fix this problem Another way is converting to closed polyline and then convert them to a polygon.I managed to convert them to polyline using geoprosseor but can not convert it to polygons.(I can not use feature to polygon either because of licence problam). So what ways do you suggest me to do this task?I can not find any samples or something like that.Can you please help me do this? I am using ArcObjects 10.2.1 and C# thanks a lot


title: Best ArcGISTiledMapServiceLayer Caching technique

text: 1) I have small portion of data from that I have to create **ArcGISTiledMapServiceLayer** .I want to cache only that portion of the map where data is present rather than caching for the entire earth surface of blank tiles. 2) There are some basemap where in they show data is not present in area on the tile map what is what is the technique behind that. 3)Best practice to cache the map service. I am using ArcGIS server 10.1.


title: What is the best tool to measure latitude the most accurately?

text: I'm thinking it's the sextant but not sure. I tried looking it up on Google and came across this site but have found nothing useful. I am doing this for a research paper and need help. Thanks to anyone who may reply.


title: Change Defination query of a layer in ArcGISDynamicMapServiceLayer

text: Is there any solid example to change the definition query at run time of the **ArcGISDynamicMapServiceLayer** using javascript API 3.5 or higher version and ArcGIS server 10.1. I found API documentation setdynamiclayerinfos but was not able to find example to find how it exactly works.


title: Calculating countries flown over on flight path

text: Say I have two airports given. I can then get their coordinates. Now that I have the beginning and end coordinates I can draw flight-path from the two locations. Now my question is that given two end points is there a "preferred" way to calculate which nations the flight path flies over? My two possibilities: Approach 1: You "ping" the flight path at increments to see where the flight currently is using google's Geocoder. However I am under the impression this would cost a large amount of data. Approach 2: Find some sort of polygonal data to represent the nations and see when countries' borders the flightpath intersect. Now I'd prefer approach two; but I am not sure where to find good resources/data to use. Any suggestions on a best approach? Thanks guys. EDIT 1: So far information I've found are as follows;   1. Link to http://www.naturalearthdata.com/features/   2. Link below from the first responder   3. http://thematicmapping.org/downloads/world_borders.php I notice that world_borders uses `shapefile`s; is there a way to read/import those into a java environment? Any open-source engines that can read the data? I am assuming that the `shapefile` is what holds the polygonal country data? **Edit 2** For future readers: I found that for programming Java; the link in the comment of the answer below had good information. namely using GeoTools to complete the project they have a step-by-step guide for Netbeans as well as Eclipse that I will report back with if it solved my problem. Furthermore, I am using their step-by-step guide in conjuction with data from NaturalEarth.


title: Many to Many Efficient Routing

text: Not sure if this has been posted yet. But I couldn't find any software that can do this yet. I am given a large set of random locations and a set of set office location. I would like to map the random locations to the set of office locations based on travel time and efficiency od the route and mapping. So basically, given a time or distance that each office can cover(which can change base on the office), is there anyway to create a most efficient routing? Any help is greatly appreciated.


title: WPS Process call

text: I'm simply trying to call a WPS process as shown on the following page: http://suite.opengeo.org/docs/processing/wpsclient/index.html The JTS:splitPolygon example they give works just fine for me, but if I try any other process I get nothing but errors. For example:               var myline = OpenLayers.Geometry.fromWKT('LINESTRING(117 22,112 18,118 13, 115 8)');     client.execute({         server: "local",         process: "JTS:isSimple",         inputs: {             geom: myline         },         success: function (outputs) {             console.log(outputs);         }     });      This throws an "Uncaught TypeError: Cannot read property 'supported' of undefined" error. What am I missing?


title: Geoserver Style Error

text: I was applying a style and there was no validation error. But when I am about to press submit I am getting the following error org.apache.wicket.WicketRuntimeException: java.io.FileNotFoundException: /var/lib/tomcat7/webapps/geoserver/data/styles/bangalore_poi_others.sld (Permission denied) I am using geoserver 2.1.2 and apache tomcat 7.


title: How do I break up a shape file into contiguous area?

text: I've got a shapefile which is a union of non-contiguous shapes. I'd like to be able to break them up into individual groups. For example, suppose I have a polygon of the entire US, currently only queryable as one large shape in PostGIS. I want to be able to query Alaska, the Aleutian Islands, each of the Hawaiian Islands, Puerto Rico, Guam, and the continental states separately. How would I break up that polygon into those different groups to do this? Assume that only land masses are specified as a shape and that I have to be able to generalize to N different non-contiguous shapes–I can't do the groupings by hand in software like ArcGIS or QGIS.


title: Can'tinstall arcgis for server 10.2

text: i have a problem installing arcgis server 10.2, i installed arcgis desktop 10.2.1, But when trying to install arcgis server it said that it wont be installed until i upgrated this list of product to 10.2, but i get empty list. i have a 64bits machine, and my arcgis for desktop is installed on program files (*86) could it be the source of the problem?


title: How to Get a Vector Layer to add to a QGIS Project

text: I used QGIS a little bit a few years ago and used some layers for my city that I had. I have a simple project that I would like to use it for and am having a little problem. I would like to simply get an orthophoto of a parcel that we recently purchased. I found this thread. It is practically what I want to do. I just want the data for another township in the same state. I can set up the coordinate system and know the WMS layer I would like to use for the orthos, but, the instructions say to add a vector layer of the specific area I'm interested in. The example gives Acme Township as an example. The instructions don't tell me where I can find the layer so I can add it. In my case I am looking for Mecosta Township, MI. Can anyone tell me where I might be able to find some sort of shapefile for Mecosta Township. If the layer includes parcels that would be ideal.


title: GeoRSS reprojection

text: Sorry I'm newbie in ESRI API javascript, I'm trying to to change the GeoRSS projection from from 102100 to 54030. could anyone give me a hand?! please! my project script is:               var map;     require([     "esri/map", "esri/SpatialReference", "esri/layers/ArcGISTiledMapServiceLayer",   "esri/layers/GeoRSSLayer", "esri/InfoTemplate",     "dojo/parser", "dojo/_base/array", "dojo/dom-style",      "dijit/layout/BorderContainer", "dijit/layout/ContentPane", "dojo/domReady!"     ],      function(     Map, Tiled, GeoRSSLayer, InfoTemplate,     parser, arrayUtils, domStyle     ) {     map = new Map ("map");     var tiled = new Tiled ("//.../arcgis/rest/services/.../MapServer");     //==> ROBINSON     map.addLayer (tiled);               // create layout dijits          parser.parse();          var georssUrl = "...gml";      //==> MERCATOR      var georss = new GeoRSSLayer(georssUrl, SpatialReference);      georss.on("load", function() {       domStyle.set("loading", "display", "none");       // create an info template       var template = new InfoTemplate("${name}", "${description}");       // set the info template for the feature layers that make up the GeoRSS layer       // the GeoRSS layer contains one feature layer for each geometry type       var SpatialReference=new SpatialReference (54030);       var layers = georss.getFeatureLayers();       arrayUtils.forEach(layers, function(l) {         l.setInfoTemplate(template);               });     });      map.addLayer(georss);     });


title: Why openlayers show wrong measures of line from ECW geoserver?

text: I have loaded the following raster into geoserver and added measurement control. The thing is that the measure gives me wrong result (it shows the car length is 49m) Any idea what can cause this ?! For the code , I used the generated code from "Preview Layer" in geoserver and just added the dynamic measure control to it. ![enter image description here](http://i.stack.imgur.com/LVuka.png)


title: get wkt from geometry in lat lon format

text: Using the **plane geometry editor** plugin in qgis 2.0, WKT representation of a feature is shown. The result is in the `lon1 lat1,lon2 lat2` format. e.g. `POLYGON((85.43523762804622379 27.69592296476042037,85.43592529210572195 27.69598802249519665,85.43609697201434017 27.69598624204515147,85.43523762804622379 27.69592296476042037))` I want result in this format `lat1 lon1 lat2 lon2`. Is there a simple easy way of doing this?


title: Style Manager Symbols and Ramps Missing in QGIS 2.0

text: This is very frustrating. I already posted a question about the missing color ramps in QGIS 2.0. I´ve even installed the development version without success. No color ramps available, and I have to create a new one everytime (or save a style which is not the right way to do it). Also, i found that the style manager is completely empty. Even if I save my symbols in the Style manager. After I start a new session, the symbol is lost... Same happens with the color ramps I create. Any help would be much appreciated. Andres


title: Failure to install QGIS 2 - held broken packages

text: I pasted the issue about not being able to install QGIS 2 on Ubuntu Raring on an ubuntu forum: launchpad.net apt questions and the answer seemed to indicate the problem is with QGIS? Thread: https://answers.launchpad.net/ubuntu/+source/apt/+question/236907 A copy of the thread is attached via the preceding link hoping there's someone who knows what is going on with this program. Unfortunately I can't use Linux or a version of it at the moment as I can't get QGIS to work on it so your help would be much appreciated.


title: Census.gov TIGER down, suggestions?

text: Census.gov is down for now. I have a friend who is a teaching a GIS class in the San Francisco bay area and only has so many datasets downloaded for specific projects. Does anyone have a large copy of this area or know of any mirror that the data could be downloaded from? The class is using QGIS and is going to stall or have to start using data from another source for now to continue. Any suggestions on other data sets to explore?


title: Local geolocation to FIPS county using Latitude and Longitude

text: I millions of lat,long pairs and I want to convert to FIPS county codes that match with the US census data. In the past I used the census block api but that is too slow for the amount of numbers that I have. Is there a better library / tool that I can load locally ? If not would be possible to use the TIGER data from the census website to create a service like that ? What knowledge is required (projections, minimum distances, ...) ?


title: Export topology tables to txt files using ArcObjects, generating directed graph, and finding intersections

text: I want to export topology tables to text files:               for points:        1  id  x   y          for lines:     line id   fompoint.x  from point.y  topoint.x  topoint.y          for polygons:     id   line id 1, line id 2...          2   arc list     3   node list      Using ArcObjects VBA.


title: Distribution of Human Population with Altitude

text: I recently found Joel E. Cohen and Christopher Small's work on hypsographic demography, i.e. the distribution of human population with altitude. Does anyone know of any other research in this area? I am looking for ideas on how to visualize the relationship between elevation and population in a map and would like to see what has been done already. I could not find any maps beyond those created by Cohen and Small.


title: raster not displaying below certain scale

text: i've clipped a raster in 20 smaller images and i can zoom in until 1:150 in the original raster but, in the smaller ones, when i go below 1:4000, the workspace goes blank. i'm working with arcgis 9.3 and can't find any configuration/solution to this problem.


title: QGIS GDAL polygonize - how to set scale

text: I've got some old contour plots, paper copy only, that I'd like to be able to put into QGIS. I get as far as a successful polygonize into shp and it looks fine, but the scale is wrong. How might I either: 1. tell the raster routine the correct scale or 2. change the scale of a layer to the correct value in QGIS?


title: GeoJSON vector layer created using OpenLayers is broken

text: What i want to do here is to create a vector layer (using openlayers) from a GeoServer JSON Response. Actually i'm doing this:               var url = "http://"+ip+":"+puerto+"/geoserver/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=tiger:poi&outputFormat=text/javascript&format_options=callback:getJson";          jQuery.ajax({                         jsonp: false,                         jsonpCallback: 'getJson',                         type: 'GET',                         url: url,                         async: false,                         dataType: 'jsonp',                         success: function(jsonObject) {                             alert(jsonObject);                             //handler(jsonObject);                             create_GeoJSON(jsonObject);                         }     });      With the code above (after a long time of headaches with cross domain issues) we can retrieve data from GeoServer (through WFS) using JSONP. The function below handles the response, tries to create the vector layer and put that vector over the map, using GeoJSON.               function create_GeoJSON(json){              var in_options = {                 'internalProjection': window.map.projection,                 'externalProjection': window.map.projection                  };              var jsonString = JSON.stringify(json);              var geojson_format = new OpenLayers.Format.GeoJSON(in_options);         var vector_layer = new OpenLayers.Layer.Vector("Consulados");         vector_layer.addFeatures(geojson_format.read(jsonString, "FeatureCollection"));         window.map.addLayer(vector_layer); }      Now, the problem is that the map isn't showing the layer components properly. As you can see i only see a point (or maybe a lot of superposed) drawed over the (0,0). This seems to be an issue related with the coordinates, projections or something like that. By the way, here is the map init:               var options = {        projection: "EPSG:900913",        displayProjection: "EPSG:4326",        controls: [],     };          window.map = new OpenLayers.Map('map_canvas', options);      What you think about that? I "waste" a lot of time over this right now and still don't have a solution. Thanks. PS: Feel free to edit the post, my english is a little poor. ![How the map looks like](http://i.stack.imgur.com/Sndnt.png)