dimanche 30 novembre 2014

shapefile for mumbai city area or section wise


please don't kill me if this question is listed anywhere else, i really searched high and low to find answer to this question.


i am working on a project in which i need to plot data on Mumbai city map section wise or area wise but i am not able to get Shapefile for the same. can anyone please provide me the link where i can get this Shapefile or guide me on creating one using QGIS.I mentioned QGIS because i know a little about OGIS and hence it will save time.Also i tried certain links which are listed below. but i couldn't find what i need.Guys it would be a great help because i really need this and i am stuck on this from few days.Also i mentioned shapefile because these are easy to manipulate, later i will convert it into geoJSON and work on it so if geoJSON files are there then it would also help. An image is also provided below.


Note: only the important links are listed: 1. http://ift.tt/1ysNCJn 2.http://ift.tt/1jOi21j 3.http://ift.tt/1ysNAkP this is what i want





What QGIS line styles work best with Garmin maps?


I am using QGIS to create some maps to load on to a GPS using the plugin. When I display the maps on the GPS the lines are very faint and I can not distinguish colours and some colours vanish altogether. But even the black lines don't stand out.


I have tried increasing the thickness of the lines in QGIS but this does not seem to make any difference.


I am displaying the maps on an ETrex 30.





arcgis flex api queryTask keep loading


In my flash builder mxml file, I load a datagrid with data get from database through web service. In the datagrid I create a context menu "select and zoom" to select and zoom to a feature (attribute value from the datagrid selected row), but when I select a row, and call this "select and zoom" function, the busy cursor will always there.


From the map service page, I can successfully query the same feature as follows:


enter image description here


I can also confirm that the selected feature exist and loaded, but the situation is when I call the select function on datagrid selected row, it keep loading.


I am stuck at the problem for quite a long time but still cannot figure out what is wrong.


Don't know whether it is code issue or setting issue, my mxml code is listed below:



<viewer:BaseWidget xmlns:fx="http://ift.tt/rYUnH9"
xmlns:s="library://ns.adobe.com/flex/spark"
xmlns:mx="library://ns.adobe.com/flex/mx"
xmlns:viewer="com.esri.viewer.*"
xmlns:esri="http://ift.tt/1tsyegi"
right="420" bottom="220" width="850" height="449" layout="absolute"
widgetConfigLoaded="init()">
<fx:Script>
<![CDATA[
import com.esri.ags.FeatureSet;
import com.esri.ags.Graphic;
import com.esri.ags.geometry.Extent;
import com.esri.ags.geometry.Polyline;
import com.esri.ags.layers.GraphicsLayer;
import com.esri.ags.tasks.QueryTask;
import com.esri.ags.tasks.supportClasses.Query;
import com.esri.ags.layers.FeatureLayer;

import mx.collections.ArrayCollection;
import mx.controls.Alert;
import mx.rpc.AsyncResponder;
import mx.rpc.events.FaultEvent;
import mx.rpc.events.ResultEvent;
import mx.utils.StringUtil;


private var graphicsLayer:GraphicsLayer;
private var hashmapOfExistingGraphics:Object = new Object();

public function init():void {

graphicsLayer = new GraphicsLayer();
graphicsLayer.symbol = sfs;
map.addLayer(graphicsLayer);
ws.GetLoc();

var cm:ContextMenu = new ContextMenu();
cm.hideBuiltInItems();
var item:ContextMenuItem = new ContextMenuItem("Select and Zoom");
item.addEventListener(ContextMenuEvent.MENU_ITEM_SELECT, onSelectItem);
cm.customItems.push(item);
datagrid.contextMenu = cm;
}

private function onSelectItem(e:ContextMenuEvent):void{
//event.contextMenuOwner
if(datagrid.selectedIndex >= 0){
this.cursorManager.setBusyCursor();
var query:Query = new Query;
query.where = "SELECT * from dev_db.DBO.BQ_Segment_1 WHERE BQ_no_ = '" + datagrid.selectedItem.bqNo + "' and Segment_ID = '" + datagrid.selectedItem.segmentID + "'";
queryTask.execute(query, new AsyncResponder(onResult, onFault));

function onResult(featureSet:FeatureSet, token:Object = null):void
{
graphicsLayer.clear();
if (featureSet.features.length == 0)
{
Alert.show("No feature found. Please try again.");
}
else
{
Alert.show("feature exist....");
var unionExtent:Extent;
var myFirstGraphic:Graphic = featureSet.features[0];

unionExtent = Polyline(myFirstGraphic.geometry).extent; // problem line

for each (var myGraphic1:Graphic in featureSet.features)
{
graphicsLayer.add(myGraphic1);
unionExtent = unionExtent.union(Polyline(myGraphic1.geometry).extent);
}
map.extent = unionExtent;
}
}

function onFault(info:Object, token:Object = null):void
{
Alert.show(info.toString());
}
}else{
Alert.show("Select a record first","ATTENTION");
}
}



public function GetLoc(event:ResultEvent):void {
// Databind data from webservice to datagrid
datagrid.dataProvider = event.result;
}

public function fault(event:FaultEvent):void {
// Oppps some error occured
Alert.show(event.toString());
}
]]>
</fx:Script>
<fx:Declarations>
<mx:WebService id="ws" wsdl="http://localhost/TestWebService/Service.asmx?WSDL" fault="fault(event)">
<mx:operation
name="GetLoc"
resultFormat="object"
result="GetLoc(event)"
/>
</mx:WebService>
<esri:SimpleFillSymbol id="sfs" alpha="0.7" color="0xFF0000"/>
<esri:QueryTask id="queryTask"
showBusyCursor="true"
url="http://localhost:6080/arcgis/rest/services/MyService/MapServer/1"
useAMF="false"/>

<esri:PictureMarkerSymbol id="sps"/>

</fx:Declarations>

<viewer:WidgetTemplate id="wTemplate" x="0" y="2" width="834.5" height="447"
textDecoration="none">
<viewer:layout>
<s:BasicLayout/>
</viewer:layout>
<s:HGroup x="44" y="10" width="757" height="350">

<s:DataGrid id="datagrid" width="732" height="329">
<s:columns>
<s:ArrayList>
<s:GridColumn headerText="BQ No" dataField="bqNo"/>
<s:GridColumn headerText="BQ Segment" dataField="bqSegment"/>
<s:GridColumn headerText="Contract No" dataField="contractNo"/>
<s:GridColumn headerText="Project Title" dataField="projectTitle"/>
<s:GridColumn headerText="Completion Date" dataField="completionDate"/>
<s:GridColumn headerText="Drain Name" dataField="drainName"/>
<s:GridColumn headerText="Road Name" dataField="roadName"/>
<s:GridColumn headerText="Description" dataField="shape"/>
</s:ArrayList>
</s:columns>
</s:DataGrid>
</s:HGroup>

</viewer:WidgetTemplate>
</viewer:BaseWidget>


I am using ArcGIS API for Flex 3.6.





ArcGIS 10.2 "Drive time analysis" service area problem


Although reasonably experienced with ArcGIS I am a bit of a newbie with Network analyst. I have a situation in which things work Ok with the tutorial data but not with mine!


I am trying to determine an 8km road distance from a set of hospitals in East Africa. I have some roads data that I can build into a network. As far as I can see there is full connectivity. If I work using detailed polygons (which is what I want) I get this. net2


So the "service area does not go along roads A,B and C, which I would expect it to. and how on earth does it think I am going to get to Point D!!.


Using generalised polygons I get this: net2


Now I have lost the slightly anomalous point D, but the services area still does not traverse roads A and B


Any help gratefully appreciated


Rob Corner. Western Australia





ENVI 5.1 how to load two image?


I am using ENVI 5.1 and when I use the Data manager to load several image and if I want to load two adjacent image why would one load and other doesn't load. Am I missing something here ?





I am unable to georeference counties in North Carolina


I am unable to georeference counties in North Carolina. could you provide assistance?





Transfer Data from Database to Leaflet POP UP


Good Day too all! I just wanna ask a little question. Ahmm how do you show data from database to leaflet pop up?


Here is my code in php select command



<?php
$db = new PDO('mysql:host=localhost;dbname=poi', 'root', '');
$sql = "SELECT name,user_date,user_time,address,lat,lng,icon_name FROM tblmarker";

$rs = $db->query($sql);
if (!$rs) {
echo "An SQL error occured.\n";
exit;
}

$rows = array();
while($r = $rs->fetch(PDO::FETCH_ASSOC)) {
$rows[] = $r;
$name[] = $r['name'];
$user_date[] = $r['user_date'];
$user_time[] = $r['user_time'];
$address[] = $r['address'];
$icon_name[] = $r['icon_name'];
}
print json_encode($rows);
$db = NULL;
?>


and here is my code in showing the marker in the map from database



function getInfo() {
$.getJSON("get_info.php", function (data) {
for (var i = 0; i < data.length; i++) {
var location = new L.LatLng(data[i].lat, data[i].lng);
var marker = new L.Marker(location,{icon:Icon1});
var ll = marker.getLatLng();

marker.bindPopup("$name<br>$user_date<br>$user_time<br>$address<br>$icon_name").addTo(map);

}
});
}


I thought it works because when i refresh it the marker from my database shows in the map but when i click the marker for the pop up the pop up shows only is like this.



$name
$user_date
...


Whats wrong with my code? Am I missing something? TYFH and IA





Can someone help me do my 2nd year uni assignment?


I need help with an ArcGIS university assignment. I'm sure it is very basic for someone with any knowledge of ArcGIS but I'm struggling with it. Really I'd need someone to do the first part so I can use it as reference for the rest of my assignment. This is the first map that I have to make: A map showing Landform Profile dtm data appropriately symbolised with a hillshade. You should use 4 tiles and merge with the appropriate software. The data frame should be at a scale of approximately 1:50000 and a grid in British national grid units should be displayed. It should be centred on your home postcode. Is this something I can get help with here?





QGIS 2.4.0 issues with "Select By Location" for complex polygon and points and alternatives


I am attempting to use the "select by location" tool to select polygons from one layer that overlap with certain polygons in a different layer. Unfortunately this is causing QGIS to freeze. Doing a bit of research, it seems like this is a noted bug see here. However I have been able to successfully perform this operation with a simpler set of layers.


As such I am wondering if: 1.This seems like a hardware issue or a bug 2 if there is an alternative method for achieving something like this?


I have also attempted to use v.select from the processing menu tools but this seems to crash also (although I haven't tested this extensively).


I am using QGIS 2.4.0, Windows 32-bit version. I have better processing power and the windows 64 bit version at home and will be able to check tonight if this works on that computer.


Please let me know if there are any other details required, and thank you in advance.





Pavement survey implementation


I know I am supposed to be specific...BUT


I am working as part of a research team on pavement study. I need to visually display survey data in Arcmap.


The surveys consist of x transverse measurements every y longitudinal feet. I need to pick the longitudinal measurement location of the road polyline (figured that out already) and then create features to represent the traverse measurements. These features could be points or polygons.


Here is an example of what I am trying to do!


I have been stuck on this for some time and really need to start making progress! Thanks, Ryan





Overlay single, non georeferenced image in Leaflet?


I'm wondering whether it would be possible to somehow overlay a 3rd party, changing non-georeferenced image such as:


enter image description here


If the image was static, I'd use mapwarper.net. But it's not. So is there a way to add coordinates to the image's corners, then bend it on the fly, allowing it to be overlaid?





export data error message in ArcGIS desktop 10.2.2


When trying to do export a layer to create anther shapefile I receive this message: It's only giving me this message for this layer. I tried on another layer and it seemed to work.


Thanks in advance.


enter image description here





Highlighting features in qgis2leaf export


I have exported a QGIS map to Leaflet using the qgis2leaf plugin. I am looking for some assistance as to how to highlight features (polylines) when the mouse hovers over it. I just want the polyline to stand out (become bold); i only want popups to appear when the feature is clicked (which I have already customized).

I have seen plenty of examples such as the Leaflet Cloropleth tutorial but I have two geojson files in my map that users can click between (so i'm guessing that any 'highlighting' code would need to be entered twice.


The two geojson files are exp_AllRoutesAM and exp_AllRoutesPM. Here is the leaflet code that shows how the exp_AllRoutesPM file is handled.



function pop_AllRoutesPM(feature, layer) {
var popupContent = '<b>Survey: ' + feature.properties.Survey + '</b>' + '<br/>' + 'Section: ' + feature.properties.Section + '<br/>' + 'Speed In 2013: ' + feature.properties.OUT_2013 + ' km/h' +'<br/>' + 'Speed In 2012: ' + feature.properties.OUT_2012 + ' km/h' +'<br/>' + 'Speed In 2011: ' + feature.properties.OUT_2011 + ' km/h' + '<br/>' + 'Speed In 2010: ' + feature.properties.OUT_2010 + ' km/h' + '<br/>' + 'Speed In 2009: ' + feature.properties.OUT_2009 + ' km/h' ;
layer.bindPopup(popupContent);
}

var exp_AllRoutesPMJSON = new L.geoJson(exp_AllRoutesPM,{
onEachFeature: pop_AllRoutesPM,
style: function (feature) {
return {weight: feature.properties.radius_qgis2leaf,
color: feature.properties.color_qgis2leaf,
opacity: feature.properties.transp_qgis2leaf,
fillOpacity: feature.properties.transp_fill_qgis2leaf};
}
});
feature_group.addLayer(exp_AllRoutesPMJSON);
//add comment sign to hide this layer on the map in the initial view.
exp_AllRoutesPMJSON.addTo(map);


I hope that someone can help. Thank you.





How to display the slopes that only facing west using ArcGIS for Desktop?


Here is what I am struggling with ...


In ArcGIS for Desktop I have to Derive and symbolise areas considered susceptible to avalanche risk based on the following simple rules:



  1. areas with a slope between 20 and 45 degrees,

  2. facing towards the prevailing wind (west in this case)

  3. not in forested areas.


By far I manage to display the Slopes and Aspects using tools from Spatial Analyst.


But how I can display the slopes in 3 different categories where they all facing west?


Photo attached shows what I manage to come up with


Any help will be much appreciate! thxenter image description here





How to create this raster?


i got three layers:



  1. building data(polygon) contains park lot capacity in it

  2. on the road parking data(line) also contains parking space data

  3. car parks (polygon) includes parking space capacity


i would like to create a raster, the raster must bu in 50x50 m resolution, in every 50m square i want to know the parking capacity in it, how to i divide my city to pixels and calculate this parking space capacity in the pixels. I'm using arcgis 10.1





Buffer is not displaying on the map?


Im trying to create a buffer around 5 points on the map.


After input the file and selecting distance of 3km the buffer is successfully created and added to table of content, however there's no color circles displaying on the map.


Can anyone help me with it? thx





Import .gfs file into PostGIS?


Is it possible to import a .gfs file into PostGIS?


I can open the .gfs file in QGis, so I assume it is a geographic file of some kind: but I'm not sure how to get it into PostGIS.


I have Googled but can't find any examples of doing this with ogr2ogr. Do I need to convert it into some intermediate format first?





Spatial join by categories in QGIS


I am using QGIS 2.4.0 and I have a vector layer which contains a remoteness indicator which is broken up into 5 categories. I would like to turn this into a layer that has these 5 categories collapsed into 2. On the included map I would like to join the green areas together into one group and the yellow and red areas into another group. Any advice would be much appreciated enter image description here


In the attribute table included I would like to join together the categories in the RA_Name11 field.


enter image description here


My apologies if this is a simple task, I just haven't got a good starting point for identifying keywords.





ogr2ogr: Can I add a field when I import to PostGIS?


I have a bunch of GeoJSON files that I want to import to a PostGIS database.


I can import individual files perfectly as follows:



ogr2ogr -f PostgreSQL PG:"host=x user=x dbname=x password=x" myfile.geojson -nln my_table


But as I import them, I want to add an extra column to these rows in my_table for the filename, so that I know which source file each field came from.


Is this possible using ogr2ogr?





QGIS georeferencer - GCP projection problem


In short, the problem is that when Geoereferencer is reopened with a previously edited raster, all the previously saved control points have moved off the image to a location thousands of kilometres away. (QGIS 2.6, Windows 64 bit).


Here is what I did. The project CRS is set to GDA94/MGA zone 55 (EPSG 28355), OTF projection is ‘on’ (to project some GDA94 lat/lon layers to the project’s MGA zone 55 tranverse Mercator projection). I open Georeferencer, and choose a raster file to georeference. My default is to prompt for a projection, so I choose MGA zone 55.


In Georeferencer, I create a number of GCPs (i.e. control points) by referencing to objects in the map canvas. After adding points, I start the georeferencing calculation, choosing Thin Plate Spine (cubic) for full rubber-sheeting. Operation is successful, the raster appears on the map canvas, along with the control points. The control points are also visible in the Georeferencer window, still overlying the points chosen for georeferencing. I save the control points and close the Georeferencer. All good so far.


Now for the problem. Having closed the Georeferencer, I reopen it, and open the same raster image. But now, the control points in the Georeferencer window no longer overlay the raster in the Georeferencer window. However, they do appear correctly on the map canvas (the main QGIS screen). In the Georeferencer window, they have in fact displayed, but half a world away to the south. That is, if the raster covers a rectangle with coordinates 300,000, 5834300 ; 362100, 5782700, the control points occupy an area bounded by coordinates 0,0 ; 6800, -5900. These small coordinates seem to correspond to the number of pixels in the raster image. So, it is no longer practically possible to edit the GCPs.


How do I solve the problem?


See screenshots below.


Image 1: Georeferencer window, and map canvas. After reopening Georeferencer, GCPs do not show in window, but show on canvas (the small red dots).


Image 2: Zooming out in Georeferencer window shows the GCPs are there but far away.


Image1 image2





QGIS crash when I try print geometry


I do (It`s work):



from osgeo import ogr
inputds = ogr.Open(r'c:\temp.shp')
inputlyr = inputds .GetLayer()
for i,feature in enumerate(inputlyr):
print i, feature.GetGeometryRef()


Output:


0 POLYGON ((4.286549707602339 0.526315789473684,3.797114870553876 -2.563854154275787,2.376719651351817 -5.351536733451042,0.164402230527078 -7.563854154275783,-2.623280348648173 -8.984249373477848,-5.713450292397642 -9.473684210526315,-8.803620236147113 -8.984249373477859,-11.591302815322369 -7.563854154275806,-13.803620236147115 -5.351536733451073,-15.224015455349186 -2.563854154275825,-15.71345029239766 0.526315789473643,-15.22401545534921 3.616485733223115,-13.803620236147164 6.404168312398374,-11.591302815322436 8.616485733223128,-8.803620236147188 10.036880952425204,-5.713450292397716 10.526315789473685,-2.623280348648239 10.036880952425237,0.164402230527027 8.616485733223191,2.37671965135178 6.40416831239846,3.797114870553858 3.616485733223211,4.286549707602339 0.526315789473684))


1 POLYGON ((3.894736842105264 1.16374269005848,3.4053020050568 -1.926427253690991,1.984906785854742 -4.714109832866246,-0.227410634969997 -6.926427253690988,-3.015093214145249 -8.346822472893052,-6.105263157894718 -8.836257309941519,-9.195433101644188 -8.346822472893063,-11.983115680819445 -6.926427253691011,-14.195433101644191 -4.714109832866277,-15.61582832084626 -1.926427253691029,-16.105263157894736 1.163742690058438,-15.615828320846287 4.253912633807911,-14.195433101644239 7.04159521298317,-11.983115680819513 9.253912633807921,-9.195433101644262 10.674307853009999,-6.105263157894791 11.163742690058481,-3.015093214145315 10.674307853010031,-0.227410634970049 9.253912633807985,1.984906785854705 7.041595212983256,3.405302005056782 4.253912633808007,3.894736842105264 1.16374269005848)) ... ...


but when I try add:



print inputlyr[0].GetGeometryRef()


QGIS crash - Crash dumped. minidump written to C...


But my goal is not to print the geometry.. I want to Union geometry



poly1 = inputlyr[0].GetGeometryRef()
poly2 = inputlyr[1].GetGeometryRef()
union = poly1.Union(poly2)
print union.ExportToWkt()


But I receives the same error. I test this on 2 comp (QGIS 2.6 and QGIS 2.2)


I do something wrong?





Locked out of all instances of Geoserver due to password problem


I beginning to deploy Geoserver (2.6.1), Tomcat 7 on Windows Server 2008 and everything was going just fine until I changed the master password but I now find that I am completely locked out of my instance since my password is not now recognised (not sure why, possible initial mistyping). I decided to start again so uninstalled Java, Tomcat + Geoserver but with a fresh install of Java, Tomcat and Geoserver I remain locked out of any instance of Geoserver whether the default admin/geoserver or my password. I changed the user.xml to admin/geoserver as others have advised, read the security documentation, went through the whole uninstall process again but I’m still locked out. I can only assume that the master password must be buried very deep somewhere in the operating system and re-surfaces for any new instance. Hoping I don't have to re-install the server OS! Any help appreciated.





What causes ERROR 000210: Cannot create output in_memory\AND from ArcPy?


The ArcGIS for Desktop "ERROR 000210: Cannot create output" has been reported on this site a number of times but an in_memory workspace seems to have been implicated only once:



ERROR 000210: Cannot create output in_memory



which was in Why do I get 000210 Error trying to wirte output to in_memory workspace in ModelBuilder


I have just encountered the same error message from ArcPy, using ArcGIS 10.2.2 for Desktop both in a Python script tool and from IDLE.


I have determined the cause and so I am asking a question in order to provide an answer to anyone who may encounter it in future.


What causes the error message below from ArcPy?



ExecuteError: ERROR 000210: Cannot create output in_memory\AND
Failed to execute (CopyFeatures)




PyQGIS export to JSON: writeAsVectorFormat and precission


In my QGIS plugin qgis2leaf I am exporting all layers as JSON into a JSON directory. Therefore I would like to enable the user to define the precision in the GUI. The normal "save as" dialogue in QGIS has the possibility to define the precision:


GEOJSON export with the save as dialogue


Is there any possibility to use a precision in the writeAsVectorFormat? The found documentation gives me no clue:



QgsVectorFileWriter::WriterError QgsVectorFileWriter::writeAsVectorFormat ( QgsVectorLayer * layer,
const QString & fileName,
const QString & fileEncoding,
const QgsCoordinateReferenceSystem * destCRS,
const QString & driverName = "ESRI Shapefile",
bool onlySelected = false,
QString * errorMessage = 0,
const QStringList & datasourceOptions = QStringList(),
const QStringList & layerOptions = QStringList(),
bool skipAttributeCreation = false,
QString * newFilename = 0
) [static]




Using valueAsText, giving me unicode data instead of a list


I am learning Python to create toolbox right now. When I'm trying to get the parameter, using valueAsText, the date returned is unicode instead of a list that I was expecting.


This is the code that I was using:



inFeatures = parameters[0].valueAsText


I was hoping the input data would be stroed in a list in this way, put when I use for loop, this is what I have:



messages.addMessage(type(inFeatures))
for inFeature in inFeatures:
messages.addMessage(inFeature)


enter image description here


It is dividing the path into single characters.


Any way I could get the result in a list?


Thanks.





Using tool validation to create parameter based on fields of another parameters input?


I'm assuming I'm missing something and this is a pretty easy thing to do but I'm stumped...


I'm working on a script that currently takes and input excel file and uses a "hard coded" field name that I know is in the table. I would like to make the tool a little more user friendly in that I would like to have a drop-down of field values from within the table, letting the user select the appropriate field from the table. This would ensure that the field being used is indeed in the table.


So... If I have an input parameter for a file in which the user selects an excel workbook, how would I create another input parameter that would then be a drop-down menu of field headers from the selected excel workbook. Is there some kind of validation code for this?





Possible to show aggregate votes in a tooltip on a markercluster?


I'm building an interface that would take poll data in and show some results over a marker cluster using Mapbox but I don't have any experience with it.


The interface would look much like https://www.mapbox.com/mapbox.js/example/v1.0.0/leaflet-markercluster/


If you click on the cluster instead of zooming in a tooltip would show up with a breakdown of (yes/no) votes in a pie chart like this https://www.dropbox.com/s/hv4c8708wlhe4mt/Screenshot%202014-11-30%2015.27.53.png?dl=0


In the GeoJSON data I was thinking each data point would have a property to state the vote value and somehow I could look at them in aggregate at the cluster level


geoJSON structure



{
"type": "Feature",
"geometry": {
"type": "Point",
"coordinates": [125.6, 10.1]
},
"properties": {
"vote": "yes" //no
}
}


Is this possible using these libraries and geoJSON data structure? Am I on the right track? Thanks for your guidance





SQL Server ShortestLineTo doesn't return an intersecting line


I have a LineString describing a path, and a Point, which is not on the line. I am trying to find the closest point to my POI on the line.

I figured I'd get the ShortestLineTo between the two objects, and then intersect it with the line to find that point. It doesn't seem to work - the line does not intersect with the shortest line. I suspect it's because of different rounding along the way.

here is my code:



declare @section_a geography;
declare @poi geography;
declare @shortest_line geography;
declare @intersection geography;

set @section_a = geography::STGeomFromText('LINESTRING (-116.4703362 32.6030246 0 1880.19645062108, -116.4703371 32.603065 0 1884.67750272083, -116.4703262 32.6031078 0 1889.53294225379, -116.4702971 32.603152 0 1895.14435545006, -116.4702939 32.603177 0 1897.9330160081, -116.4703446 32.6035651 0 1941.23463376485)', 4326);
set @poi = geography::STGeomFromText('POINT (-116.470132 32.603062)', 4326);
set @shortest_line = @poi.ShortestLineTo(@section_a);
set @intersection = @shortest_line.STIntersection(@section_a);

select'section a' name, @section_a data, @section_a.ToString() toString
union all
select 'poi', @poi.STBuffer(1), @poi.ToString()
union all
select 'shortest line', @shortest_line, @shortest_line.ToString()
union all
select 'end point', @intersection.STBuffer(1), @intersection.ToString();


And my results show the intersection as an empty geometry collection. When I try and add an STBuffer around the shortest line, it only creates a LineString as the intersection, and I want to get the intersecting point. What am I missing?





TypeError/unbound local error - help!



I'm new to python and stuck with this error....help?

import graphics
from graphics import*

def house():
win=GraphWin(800,500)
win.setCoords(0.0,0.0,3.0,4.0)#reset coordinates
Text(Point(1.5,3.5),"click spot to designate 2 corners of house").draw(win)
p1=win.getMouse()
side1=(Point(p1.getX(),(p1.getY()))).draw(win)
p2=win.getMouse()
side2=(Point(p2.getX(),(p2.getY()))).draw
rectangle = Rectangle(Point(p1.getX()),(Point(p1.getY()))).draw(win)
house = ((side2.getX()-(side1.getX))).draw(win)
return house
house()

Python 3.4.1 (v3.4.1:c0e311e010fc, May 18 2014, 10:38:22) [MSC v.1600 32 bit (Intel)] on win32
Type "copyright", "credits" or "license()" for more information.





================================ RESTART ================================



Traceback (most recent call last):
File "E:\ICS 140\ass 8.py", line 23, in <module>
house()
File "E:\ICS 140\ass 8.py", line 20, in house
Rectangle = Rectangle(Point(p1.getX()),(Point(p1.getY()))).draw(win)
UnboundLocalError: local variable 'Rectangle' referenced before assignment


================================ RESTART ================================



Traceback (most recent call last):
File "E:\ICS 140\ass 8.py", line 23, in <module>
house()
File "E:\ICS 140\ass 8.py", line 20, in house
rectangle = Rectangle(Point(p1.getX()),(Point(p1.getY()))).draw(win)
TypeError: __init__() missing 1 required positional argument: 'y'







How to get the layer name with a arcpy code


Using the FeatureToPoint tool I've having difficulty calling just the name for file without the extension or path. Here is the part of the script that I'm creating that Im stuck on.


eg. C:\user\shapefiles\polygon.shp I just want the word polygon



import arcpy

file1 = arcpy.GetParameterAsText(1) # this is a .shp file
savepointlocation = arcpy.GetParameterAsText(2)

arcpy.FeatureToPoint_management(file1, save1, "INSIDE")


The problem is that its calling for the path of C:\user\shapefiles\polygon.shp and not just the name polygon.


I'm new and I'm sure there is some function or way to just call the name of the file without the path or the .shp in the file.





PyQGIS scripting: switch rule-based styles in loop


I'm a working on a QGIS script that could activate/deactivate rule-based styles of the vector layer in a loop (one active layer per iteration). On the every iteration this scenario should:



  1. Activate one of the prearranged style rules for a given vector layer and deactivate all other;

  2. Save current view as JPEG screenshot;


In my particular case the vector layer is a PostGIS table that contains the contours derived from some raster layer. Rules are based on the distinct values of the isolines. I take screenshot of the map for the every activated rule and finally put them all into the animated gif.


Animated gif


Since list of the rule-based styles can be very long, I'm looking torward the PyQGIS code snippet that could go through the one in a loop.





PostGIS doesn't use index for spatial query


I can't get PostGIS 2.1 running on PostgreSQL 9.3.5 to use a spatial index even for the simplest queries. The whole dataset is 8 million points (population count grid from here http://sedac.ciesin.columbia.edu/data/set/gpw-v3-population-count), the table is created as CREATE TABLE points (population DOUBLE PRECISION NOT NULL, location GEOGRAPHY(4326, POINT) NOT NULL), index is CREATE INDEX points_gix ON points USING GIST(location). The queries are as simple as they get SELECT SUM(population) FROM points WHERE ST_Distance(location, ST_GeographyFromText('SRID=4326; POINT(0 0)')) < 1000. PostgreSQL always uses Seq scan for it, I've tried a subset with 10000 points - still Seq scan. Any ideas?





How can I automatically make new computings after I have modified some value in a field in ArcGIS Desktop?


Hy all!


I'm a computer programmer and I have no idea about geographic things. I have got a long algorithm with a lot of computings, computed from a previous computed field. There are 7-8 computed fields and they refer to each other, so it's a bit complex algorithm. There is a chance, that the original values in the beginning of the long algorithms will change (some empty columns will have values and the algorithm needs to compute from them). How can I make this to be as automatic as it can be? I think there is other solution than run the python scripts on the fields again and again.


Thanks.





Generate random points in geometry near heatspots


I have an algorithm, that generates random points with in a given polygon (Generate random points in polygon. Algorithm question). The Polygons the points are generated in are communities or districts of germany. So instead of just randomly placing the points somewhere in the Polygon I would like to place them near cities, town etc. since a Point should represent a person.


So my Idea was to use OpenStreetMap data to query cities, town etc., but I yet have to find out which data from OSM I need to use. Under the assumption that the data is available how do I influence the algorithm to generate the Random points near the cities?





I am unable to remove the legend from a CartoDB layer. I have passed legends: false in createLayer options but with no luck


cartodb.createLayer(map, BikeTrails)



.on('done', function(layer) {legends:false

layerBikeTrails = layer;


L.control.layers(null, {'Dedicated Bike Lanes': layerBikeLanes, 'Bike Share Toronto': layerBikeShare, 'Bike Stores': layerBikeStores, 'Access to Multipurpose Trails': layerBikeTrails}) .addTo(map);





Generate dynamic XY coordinate of 4 corners of dataframe in Arcgis Layout!


enter image description here


I need to label the 4 corners of a dataframe in this way.


Is there a possibilities with dynamic text scripts?





Python Geolocation Software to Filter Cities By Population Within a Radius?


What I'm trying to achieve


I want to enter a US city/state or zip code value. Then enter a radius (ex: 30 - would stand for 30 miles). Then I enter a minimum population (ex: 100000 - would tell the program logic to only look for cities with populations greater than or equal to 100000). The program returns the names of cities within the radius that has populations greater than or equal to the minimum population.


I don't care about super-precise numbers. For example, if a major city barely touches the outside radius, it still counts so I want to include that city's entire population (even though it's only perhaps a small sliver of the actual population) in the output.


My Question


Is there an open-source Python software or solution out there that would allow me to achieve this? If yes, what do you recommend?


What I'm considering



  • Google geolocation but don't know if it's possible to do what I'm trying to do

  • Some Python radius software package I don't know about

  • Utilizing a giant database lookup table or an API but don't know which route is best from a strategic standpoint?





Replace Domains of ArcGIS 10.1


I need to Replace the Domain on Field level on my Geodatabase.


My Geodatabase having multiple features and tables where at each on field one Domain is assigned, I need to replace that Domain and Delete old domain.


To better way understand. (example)



  1. Table1 having Field XYZ on which Domain A1 is assigned


now I like to replace A1 with A2 and once it change for all the tables I would like to delete A1 from my Domains.


I search on internet but I could not find any helping materials.


Do I do this with performing some query on database ?


select * from sde.GDB_ITEMS


Above query get me all the domains which exists for my Geodatabase but I could not to find relation how this domains is assigned to my tables (how's the relationship between them?)


I do believe that with that relationship I can replace the domain.


I am not python developer.. but I am trying to write this script..



>>> import arcpy
... admin_workspace = "Database Connections/WGis01-SDE.sde"
... arcpy.env.workspace = admin_workspace
... domain_to_search = r"MM Work Function_1" #supply the CV domain name
... datasets = arcpy.ListDatasets()
... for dataset in datasets:
... if dataset == N.Electric:
... featureclasses = arcpy.ListFeatureClasses()
... for featureclass in featureclasses:
... print featureclass
...


Parsing error IndentationError: expected an indented block (line 9) >>>


Please Help to get it done.





How to create a map and save it to an image with GeoTools


I would like to create a map with GeoTools and save it to an image (e.g. JPEG). My requirements are simple:



  1. Create a map of the world with 2 layers: Political boundaries and a graticule. The layers are from different sources and different projections.

  2. Output the map to different projections (e.g. "EPSG:5070", "EPSG:4326", "EPSG:54012", "EPSG:54009", etc.)

  3. Clip the output to different AOIs (e.g. -124.79 to -66.9 lon, 24.4 to 49.4 lat).


I want to do this programmatically, via the API. So far, I have had limited success. I have learned to create a map and output in various projections using this approach:



//Step 1: Create map
MapContent map = new MapContent();
map.setTitle("World");

//Step 2: Set projection
CoordinateReferenceSystem crs = CRS.decode("EPSG:5070"); //Conic projection over US
MapViewport vp = map.getViewport();
vp.setCoordinateReferenceSystem(crs);

//Step 3: Add layers to map
CoordinateReferenceSystem mapCRS = map.getCoordinateReferenceSystem();
map.addLayer(reproject(getPoliticalBoundaries(), mapCRS));
map.addLayer(reproject(getGraticules(), mapCRS));

//Step 4: Save image
saveImage(map, "/temp/graticules.jpg", 800);


The save method is straight from the GeoTools website:



public void saveImage(final MapContent map, final String file, final int imageWidth) {

GTRenderer renderer = new StreamingRenderer();
renderer.setMapContent(map);

Rectangle imageBounds = null;
ReferencedEnvelope mapBounds = null;
try {
mapBounds = map.getMaxBounds();
double heightToWidth = mapBounds.getSpan(1) / mapBounds.getSpan(0);
imageBounds = new Rectangle(
0, 0, imageWidth, (int) Math.round(imageWidth * heightToWidth));

} catch (Exception e) {
// failed to access map layers
throw new RuntimeException(e);
}

BufferedImage image = new BufferedImage(imageBounds.width, imageBounds.height, BufferedImage.TYPE_INT_RGB);

Graphics2D gr = image.createGraphics();
gr.setPaint(Color.WHITE);
gr.fill(imageBounds);

try {
renderer.paint(gr, imageBounds, mapBounds);
File fileToSave = new File(file);
ImageIO.write(image, "jpeg", fileToSave);

} catch (IOException e) {
throw new RuntimeException(e);
}
}


The reproject method is my invention. It's a bit of a hack but its the only way I could find to output an image to a specific projection.



private static Layer reproject(Layer layer, CoordinateReferenceSystem mapCRS) throws Exception {

SimpleFeatureSource featureSource = (SimpleFeatureSource) layer.getFeatureSource();


//Define coordinate transformation
CoordinateReferenceSystem dataCRS = featureSource.getSchema().getCoordinateReferenceSystem();
boolean lenient = true; // allow for some error due to different datums
MathTransform transform = CRS.findMathTransform(dataCRS, mapCRS, lenient);


//Create new feature collection
SimpleFeatureCollection copy = FeatureCollections.newCollection("internal");
SimpleFeatureType featureType = SimpleFeatureTypeBuilder.retype(featureSource.getSchema(), mapCRS);
SimpleFeatureIterator iterator = featureSource.getFeatures().features();
try {

while (iterator.hasNext()) {

SimpleFeature feature = iterator.next();
Geometry geometry = (Geometry) feature.getDefaultGeometry();
Geometry geometry2 = JTS.transform(geometry, transform);
copy.add( SimpleFeatureBuilder.build( featureType, new Object[]{ geometry2 }, null) );
}

}
catch (Exception e) {
e.printStackTrace();
}
finally {
iterator.close();
}


//Return new layer
Style style = SLD.createLineStyle(Color.BLACK, 1);
layer = new FeatureLayer(copy, style);
layer.setTitle("Graticules");
return layer;
}


The output is really bad:


Output from reprojection


So, I guess I have a couple different questions:



  1. Is this the right approch? Do I really need to reproject layers manually or is the MapViewport supposed to do this for me?

  2. How do I clip the output to a specific AOI? I have tried setting the bounds using the MapViewport.setBounds(envelope) method but the saveImage method seems to ignore the bounds.

  3. How do I get my latitude lines to render as arcs? Is there a transform setting that I am missing?


I am using GeoTools 8.7.


Thanks in advance!





add points when joining two datasets together in ArcGIS


I am trying to join an excel file with a shapefile based on a zipcode. My problem is that I have multiple persons with the same zip code, but it joins only 1 of those persons. Is it possible to automatically add multiple points when trying to join several persons from the excel file with my shapefile (which only has one zip code point)?





Spatial Clustering of Point Data in shapfiles with R


Helle,


i want to divide a number of point observation into k clusters, based on their spatial attributes only.


stations.txt


Code :



cur=setwd("C:/Users")

Algerie <- shapefile("DZA_adm/DZA_adm0.shp")
# read in ascii file, and assign column names
x <- read.table('E:/cluster/stations.txt', sep=";")
names(x) <- c('long', 'lat')

# subset original object, return only x,y cols
y <- data.frame(x[,1:2])
row.names(y) <- x$cat

# simple plot of x,y data
plot(y, pch=3)


# load cluster package
library(cluster)
library(flexclust)

# figure out a good number of clusters use a range of 2 to 10 clusters, with 20 reps each
s <- stepFlexclust(y, k=2:10, nrep=20)
#plot(s)

# 5 clusters in a good number
y.pam <- pam(y, 5, stand=TRUE)

# add the clustering vector back to the original dataframe
y$cluster <- y.pam$clustering

# plot the clusters by color

plot(y$long, y$lat, col=y$cluster, main="Bugsites Spatial Clustering, 5 classes", cex=0.5, pch=16, xlab="long", ylab="lat",add=TRUE)

# add the cluster number to the original dataframe
y$cluster <- y.pam$clustering
y$orig_cat <- as.numeric(row.names(y))

# save as a text file and quit
write.table(y, file='bugsites.clust', row.names=FALSE)


enter image description here I started in R, looking for a way ploter the colorful dots Each one with its cluster on a shapefile (Algeria) and then unscrew the board (shapfiles) as cluster exists.


thank you





Why does my modulo expression return no results?


CORRECTION:


I've got a big set of contour lines with 5 m and 10 m stepping from an ArcGIS File Geodatabase (wasn't a shape file, as previously stated). From time to time I need to modulo out the 5 m lines in the viewer. Given a column "contour" this should work (shouldn't it?):


"contour" % 10 = 0


However, it returns 0 lines with both, double or integer type. The column contains just values of xx0 and xx5.


Thanks for any hint!





AttributeError: Help!


new to python(3.4.1) gui and can't fix this error...help?


import graphics from graphics import*


def house(): win=GraphWin(800,500) win.setCoords(0.0,0.0,3.0,4.0)#reset coordinates Text(Point(1.5,3.5),"click spot to designate 2 corners of house").draw(win) p1=win.getMouse() side1=(Point(p1.getX(),(p1.getY()))) side1.draw(win) p2=win.getMouse().draw(win) side2=(Point(p2.getX(),(p2.getY()))) side2.draw Rectangle=Rectangle(side1,side2).draw(win) house = ((side2.getX()-(side1.getX))).draw(win) return house house()


=ERROR<--------------------------------------------------------------------------------<





Traceback (most recent call last): File "E:/ICS 140/ass 8.py", line 28, in house() File "E:/ICS 140/ass 8.py", line 21, in house side2=(Point(p2.getX(),(p2.getY()))) AttributeError: 'NoneType' object has no attribute 'getX'








samedi 29 novembre 2014

Selected features with a join are identified improperly from a service


We couldn’t figure out why a service fails to provide correct identification for the selected feature. However, the scenarios below may highlight the issue:



  1. In the “UrbanMasterPlans_WithJoin.mxd” file, the attribute table of the “UrbanMasterPlans_PNA_Approved” layer is joined with the “Regulations” table, then the ArcGIS Server fails to identify all the selected features.


enter image description here


enter image description here



  1. In the “UrbanMasterPlans_WithoutJoin.mxd” file, the attribute table of the “UrbanMasterPlans_PNA_Approved” layer is not joined with the “Regulations” table, then the ArcGIS Server is able to identify all the selected features.


enter image description here


enter image description here


Here is data link to test the problem: http://www.mediafire.com/download/ev5wthdiijhp51n/PublishingIssue.zip


What might be the issue here? Why does the join affect the Identify results when applied in a service?





Looking for succesful implementaion of PPDM in form of geodatabase


I am have difficculties in iplementing spatial version of PPDM for our GIS database. Does anyone has a reference to a sucessful case of implementing PPDM in as a geodatabse. Or how to build a spatial PPDM database.


Thanks,





Aggregating a really big shapefile in ARC 10


I have a huge DOT shapefile I need to trim down. I have tried everything I know how to do. I have hit the wall.





How to split linestrings where start/end points don't touch nearby lines


I am trying to figure out some basic routing with PostgreSQL/pgRouting and I'm struggling a bit.


Our network of linestrings are not topologically correct, tracks that should be connected don't touch, or if they do, it's by chance. example linestrings


It is my understanding that to make a network topology for routing I will have to split all the tracks at each intersection or crossing. The problem is most the tracks end/start points don't touch, though they are close.


So I'm looking for a way to either split the tracks into segments using a tolerance for points that touch or to modify the linestrings themselves, moving the end/start point to snap to the nearby line.


But I don't know how to go about that?


I have this query for splitting the tracks, but it relies on points touching. It works for creating segments where the trail crosses another line, perfectly. But the ST_Touches() part is useless when there is no intersection.



CREATE TABLE trail_split_points as
SELECT DISTINCT ST_GeometryN(ST_Intersection(a.track, b.track), 1) as geom
FROM
trails as a, trails as b
WHERE ST_Touches(a.track, b.track)
OR ST_Crosses(a.track, b.track)
AND a.id != b.id
GROUP BY ST_Intersection(a.track, b.track);


UPDATE


I determined I had no choice but to fix my tracks so they are all topologically correct. It was a bit of a pain to figure out, but I've mostly done that.


The first steep was to get the start and end point of each line and then see if there was a nearby point on another line within <5 meters. If there was then I move that point to that nearby point.


I then use some query to populate a table full of all the points where these lines touch or cross.


Then using that new point table generate new lines for each segment between points. This new table is the one used for routing. The one I created a network map and indexes on for pgRouting.


UPDATE 2


I was snapping the nearby trail to the line, but not creating a new point on that line at the intersection. This is how I managed to do that.



ST_LineMerge(
ST_Union(
ST_Line_Substring(track, 0, ST_Line_Locate_Point(track, ST_GeomFromText('".$point['newpoint']."', 4326))),
ST_Line_Substring(track, ST_Line_Locate_Point(track, ST_GeomFromText('".$point['newpoint']."', 4326)), 1)
))




Web site to show GPX with elevation?


I notice that a GPX file contains not only lat + lon but also elevations:



<trk>
<name>Some route</name>
<trkseg>
<trkpt lat="50.64536" lon="3.05657">
<ele>23.6</ele>
</trkpt>

<trkpt lat="50.63836" lon="3.06491">
<ele>32.7</ele>
</trkpt>
</trkseg>
</trk>


I was wondering if someone knew of a web site where I could simply upload a GPX file and have it drawn on a map with the route colored differently based on the elevation, ie. where the flat parts are drawn in green while the steep parts are shown in red?


Thank you.





OpenLayers Only draw features in a vector layer within view


I am dealing with a kml that I am importing into openlayers as a vector layer. This kml has many features so it dramatically slows down the web app. I would like to make it so only the features within a certain extent are drawn (based around the user's position and also only at a higher zoom level). So far I set up an event listener to check for the zoom level. Once the appropriate zoom level is reached the vector layer is drawn however, it draws the whole layer. Is there any way to do this such that it only shows part of the KML? I can't find the answer or an example anywhere.


map.events.register("zoomend", this, function(e) { var zoom = map.getZoom(); if (zoom > 10){



kml = new OpenLayers.Layer.Vector("Deer Accident Locations", {visibility: true,
strategies: [new OpenLayers.Strategy.Fixed()],
protocol: new OpenLayers.Protocol.HTTP({
url: "deer.kml",
format: new OpenLayers.Format.KML({
extractStyles: false,
extractAttributes: true,
maxDepth: 0

}),
style: {externalGraphic: 'images/Deer_green.png' } })
});

map.addLayers([kml]);




using real time web for map api


I plan to do a master thesis on a map api showing performances of a simulated network of wind turbines located across a map(maybe using snmp polling mibs ,ems etc). I was planning on doing a case study on the different web services/apis like REST, soap, JSON-rpc, OSLC etc that could be tested and determine which would best suite this system. But i was wondering if i should instead be looking at Real time web protocols. the system only needs to update in near real time, but immediate updates wouldn't be a bad thing of course. Is using real time web alot more expensive/difficult to use, or should just using web API's/services that i have mentioned above be suffice for a college project?


Any help would be very much appreciated!!!





CartoDB named template events and info windows


I'm trying to apply what I'm learning from https://github.com/mhkeller/cartodb-templates/blob/master/basic/polygon-hover.html to this scenario except hover isn't triggering any info window (or any console logging for that matter). The map does render correctly. What is the problem here?:



var account = data.layergroupid.split('@')[0];
var url = 'http://'+account+'http://.cartodb.com/api/v1/map/'+data.layergroupid+'/{z}/{x}/{y}.png';

L.tileLayer(url).addTo(map)

.on('featureOver', function(e, pos, latlng, data) {
$('.leaflet-container').css('cursor','pointer');
if (data.cartodb_id != polygon.cartodb_id) {
drawHoverPolygon(data);
}
cartodb.log.log(e, pos, latlng, data);
})
.on('featureOut', function(e, pos, latlng, data) {
$('.leaflet-container').css('cursor','default');
removePolygon();
})
.on('error', function(err) {
cartodb.log.log('error: ' + err);
});

var polygon = {};

var polygon_style = {
color: "white",
weight: 2,
opacity: 1,
fillOpacity: .45,
fillColor: "white",
clickable: false
}
function drawHoverPolygon(data){
removePolygon();
polygon = new L.GeoJSON(JSON.parse(data.geometry), {
style: polygon_style
}).addTo(map);
polygon.cartodb_id = data.cartodb_id;
}

function removePolygon(){
map.removeLayer(polygon);
polygon.cartodb_id = null;
}




CartoDB named template programmatically set different colour for each polygon


Is there an easy way to set a different colour for each of an n number of polygons in a named template in CartoDB API?


If I'm currently creating my template this way, how could I go about having each polygon a different colour for instance by randomly generating hex codes?



{
"version": "0.0.1",
"name": "test",
"auth": {
"method": "open"
},
"layergroup": {
"layers": [{
"type": "cartodb",
"options": {
"cartocss_version": "2.1.1",
"cartocss": "#layer { polygon-fill: #FFF; }",
"sql": "select * from table"
}
}]
}
}




Network Analyst--join, iterate, change attributes and save [duplicate]



I would like to create multiple routes at the same time using the iterate function. Attached is an image of a model for three different networks. One location is the same (Bell 2) but the other location differs for each network (OID_533, OID_534).


Specifically how do I do the following?:


1) Create a model, if possible, that will generate individual routes at the same time like in the example below but will automatically recognize location Bell2 goes with location OID_533. I'm assuming a join or something, then iterate needs to occur here but not certain.


2) Change the attribute fields generated in the route file from "FirstStopID" and "LastStopID" to recognize names of the attributes of the stops--e.g. Bell2 and 533.


3) Save the routes as a shapefile or a feature class with an unique name in the model. I can do this in the ArcMap dataframe but can't seem to figure out how to do it in the model.


Thanks!


Here's a picture of the model: enter image description here





XML ElementTree in ArcPy consuming WMTS XML


How do I work with the xml Element Tree to obtain the ScaleDenominator parent where ows:Identifier = default028mm?


With a fair bit of swearing I have managed to hack my way into all of the WMTS Scale's, I noticed that the OGC ones list first, I've taken advantage of this. But I want to do it properly and read the OGC Scales only. Alas, I've hit a brick wall. Any ideas?



import arcpy
import urllib
import xml.etree.ElementTree as et

WMTS = urllib.urlopen("http://basemap.nationalmap.gov/arcgis/rest/services/USGSTopo/MapServer/WMTS/1.0.0/WMTSCapabilities.xml")
print "Started"

import xml.etree.ElementTree as et
tree = et.parse(WMTS,parser=None)
root = tree.getroot()
test = 0
MaxScale = 600000 #600K No Scales Above

for ScaleDenominator in root.iter('{http://www.opengis.net/wmts/1.0}ScaleDenominator'):
Scale = float(ScaleDenominator.text)
ScaleCon = ((Scale/90.7144671432)*95.9999999999)
print 'WMST: ', ScaleDenominator.text
print 'Esri: ', ScaleCon
if test == 0:
test = (ScaleCon+1)
if test > ScaleCon:
print 'OGC028:', ScaleCon
test = (ScaleCon+1)
if ScaleCon < MaxScale:
print '< OGC028:', ScaleCon


WMTS XML: National Map USGS Topo WMTS XML


Example





Copy Feature Class to Feature Dataset


I try to copy a feature class to a feature data set but i get this massage see the screen shot , what might be the issue here?enter image description here





Comparing unclassified super-objects with classified sub-objects


I have two shapefiles exported from ecognition. the first shapefile contains the land cover classification and the 2nd one contains unclassified superobjects of the 1st shapefile. i want to compute the percentage of each landcover classes within my superobjects. I also want to get the NDVI percentage of the superobjects based on the vegetation class. I found this link http://resources.arcgis.com/en/help/main/10.1/index.html#//000800000044000000 but i cant seem to find this in ArcMap. I dont know if I can use this to achieve what i want.





How to open a second window in QGIS plugin from a button inside the main window


If I already created a qgis plugin, and i need to create a button inside that plugin that opens a new window, how is that possible? Qt creator won't add a slot unless there are a header file. I looked inside the plugin directory and didn't find the .h file, so I compiled it using uic, but it still can't find the .h file it says : no document matching ui_test.h could be found, rebuilding the project might help.





Diverse values to one field for soil analysis


I have soil analysis data in an exel table. There are for one analysis point three depths. I want to transform these data into gis-data. What is the common way, to add three values in an field of one row or similar.


Example (exel table): point depth pH-value 1 0-20 4 1 20-40 5 1 40-60 7 2 ....


This could look like this in a shape file: point depth pH-value 1 0-20;20-40;40-60 4;5;7


I think this is not nice for doing queries. But when I make one column for each depth and soil value the table would have about 120 columns :(


Any common way to do so?





how to refresh the map automatically


I am using OL3 and geoserver and postgis as database. My database is connected to geoserver and OL3 conected to geoserver and reads the points in database as vector layer and shows on the map. my data are changed each 5 second (position of points). I want refresh the map automatically each 5 second to show the movement of the points. is there any way to do this?





How to make QGIS abbreviate labels with a dictionary?


I'm trying to use qgis "expression dialog" for labeling a layer. What I'm looking for is an expression that permits to apply some abbreviation rules on the layer. For now i'm using this



CASE
WHEN nome2 is not null THEN replace("nome1" || '$' || "nome2", 'Wasserfall', 'Wssf.')
WHEN nome2 is not null THEN replace("nome1" || '$' || "nome2", 'CascatA', 'Casc.')
WHEN nome2 is not null THEN replace("nome1" || '$' || "nome2", 'Cascate', 'Casc.e')
WHEN nome2 is null THEN replace("nome1", 'Wasserfall', 'Wssf.')
WHEN nome2 is null THEN replace("nome1", 'Cascata', 'Casc.')
WHEN nome2 is null THEN replace("nome1", 'Cascate', 'Casc.e')
END


Obviously it doesn't work because only the first condition (wasserfall) is used for the label.


It happens something like that


Wasserfall$Cascata -> Wssf.$Cascata


instead of


Wasserfall$Cascata -> Wssf.$Casc.


I can solve the problem directly in postgis (with a view or a simple update) but if it is possibile to do it in qgis it would be great. Most of all it would be great to use an abbreviation dictionary but I'm quite sure that qgis (for now) is not able to doing it. I'm on W7 qgis 2.4.





What does the Join button in the categorized style tab do?


I'm just studing gis, and I use QGIS. I have found in the layer properties some JOIN and don't understand how can I use It? Maybe someone can help me to undersan it, thank's! enter image description here





QGIS map composer - grid showing unreasonable lat/lon


I need to print a map with coordinates in the format of Degree, Minute and second. CRS is set to WGS 84 (EPSG 3857) for both layers and on the project properties it is set to enable transformation on the fly. But then, when in map composer I get the following coordinates, which make no sense to me.


enter image description here





SQL Server : find Intersecting geometries with a geometry using sql query


I asked this question on parallel site (h++p://stackoverflow.com)and think this is a better place to ask such questions I want to find all geometries in a table that intersects with spacial geometry using SQL Server, I tried to use this SQL query



DECLARE @g geometry;
set @g=Geometry::STGeomFromText('POLYGON ((578060.3439000080000000 3450265.1968000000000000, 578049.6538000090000000 3450236.4844000000000000, 578030.9746000080000000 3450193.9723000100000000, 577995.7130000080000000 3450133.1654000000000000, 577974.8281000070000000 3450090.7579000000000000, 577811.4680000080000000 3449745.6367000000000000, 577804.8125000080000000 3449733.9653000000000000, 577778.4473000080000000 3449675.1851000000000000, 577769.3170000080000000 3449669.8222000000000000, 577751.0450000080000000 3449674.6951000000000000, 577687.3179000100000000 3449718.1878000000000000, 577667.0708000070000000 3449740.7260000000000000, 577668.7438000080000000 3449750.0009000000000000, 577710.9384000080000000 3449836.2113000000000000, 577735.3070000080000000 3449889.5231000000000000, 577744.7248000080000000 3449902.4691000000000000, 577779.5568000080000000 3449976.3343000000000000, 577794.5004000070000000 3450004.6225000000000000, 577974.2512000080000000 3450373.9550000000000000, 577984.9893000070000000 3450399.4253000000000000, 577987.9717000080000000 3450404.3992000100000000, 578021.4391000080000000 3450432.3046000000000000, 578053.5285000080000000 3450441.3588000000000000, 578085.6470000080000000 3450452.9310000000000000, 578130.7084000080000000 3450465.0823000000000000, 578140.7284000080000000 3450462.3602000000000000, 578132.0483000100000000 3450431.9346000000000000, 578077.0528000090000000 3450304.0990000000000000, 578056.8904000090000000 3450265.6079000000000000, 578060.3439000080000000 3450265.1968000000000000))',32639).MakeValid();

SELECT
p.[owner_ID], p.[OBJECTID] as id,
[Shape], o.[OBJECTID], o.[Name], o.[Family], o.[Father],
o.[owner_id] as oid
FROM
[dbo].[Property] AS p
JOIN
[dbo].[OWNER] AS o ON p.owner_ID = o.OBJECTID
WHERE
[Shape] IS NOT NULL
AND o.idshahrestan = 13
AND [Shape].STIntersection(@g)


but it does not work...Is this even possible using SQL Server? I think for a better preformance If I buffer my shape with a distance around 50m and find intersecting shapes or other with some boundry boxes this worsks better but I don't know how to do that :( I need some hints about this. Thank you for your help





OpenLayers and GeoJSON


Recently I started learning Openlayers 2x and I've been struggling for a while now with a GeoJSON and it's projections.


So, here is the deal;


I successfully added wms layers from ( http://geoportal.dgu.hr/podaci-i-servisi/svi-servisi-i-aplikacije/ ) with a projection HTRS96/TM (EPSG:3765) and then when I try to add vector layer it doesn't render it. Vector layer was created using QGIS software with a projection code:


"crs": { "type": "name", "properties": { "name": "urn:ogc:def:crs:EPSG::3765" }


Here is the OpenLayers code:



var map;

function init() {
map = new OpenLayers.Map('map_element', {
projection: 'EPSG:3765',
maxExtent: new OpenLayers.Bounds(426906.034398, 5101783.64664,
468058.881494, 5123976.6968),
numZoomLevels:15,
maxResolution: 'auto',
units: 'm'
});
var wms = new OpenLayers.Layer.WMS(
'DOF 1:5000',
'http://geoportal.dgu.hr/wms?layers=DOF',
{}
);

var wms2 = new OpenLayers.Layer.WMS(
'TK 1:25 000',
'http://geoportal.dgu.hr/wms?layers=TK25',
{}
);

var wms3 = new OpenLayers.Layer.WMS(
'HOK',
'http://geoportal.dgu.hr/wms?layers=HOK',
{}
);

map.addControl(new OpenLayers.Control.LayerSwitcher({}));
map.addControl(new OpenLayers.Control.MousePosition());


var KZZ = new OpenLayers.Layer.Vector("KZZ", {
protocol: new OpenLayers.Protocol.HTTP({
url: "KZZ_JLS.json",
format: new OpenLayers.Format.GeoJSON()
}),
strategies: [
new OpenLayers.Strategy.Fixed()
]
});

map.addLayers([wms,wms2,wms3,KZZ]);
map.zoomToMaxExtent();
}


And here is GeoJSON layer (fragment of it, whole file can be downloaded from this link https://dl.dropboxusercontent.com/u/126827789/KZZ_JLS.json):



{
"type": "FeatureCollection",
"crs": { "type": "name", "properties": { "name": "urn:ogc:def:crs:EPSG::3765" } },
"features": [
{ "type": "Feature", "properties": { "POVRINA": 27.027842695, "NAZIV": "Stubièke Toplice" }, "geometry": { "type": "Polygon", "coordinates": [ [ [ 457191.01612696185, 5094147.6274271049 ], [ 457199.98992743361, 5094130.2420611409 ], [ 457231.34151837695, 5094112.0906051407 ], [ 457283.71295271069, 5094094.2316290718 ], [ 457310.09355065826,....


How can I render GeoJSON and add it to this file? I'm really struggling with this and any help would be much appreciated. (this problematic openlayers map can be viewed here htt()p://interaktiva.pregrada.hr/nc/--> to low reputation to add 2nd link on my post, pls remove ())


Thank you!





How to solve coordinate differences between packages [on hold]


I have two maps generated by two different methods. The first, I used maps package and the function map() and symbols(). On the other hand, I used the ggmap package with the ggmap() and geom_point() functions. This the structure of my data and my code:



> str(city_frequency)
'data.frame': 144 obs. of 4 variables:
$ city: Factor w/ 144 levels "A Coruna","Agen",..: 1 2 3 4 5 6 7 8 9 10 ...
$ freq: int 7 2 4 1 4 2 1 1 1 5 ...
$ lon : num -8.412 0.616 -8.446 -4.262 2.149 ...
$ lat : num 43.4 44.2 40.6 42.8 43.9 ...

xlim <- c(-13.08, 8.68)
ylim <- c(34.87, 49.50)
map("world", lwd=0.05, xlim=xlim, ylim=ylim)
symbols(city_frequency$lon, city_frequency$lat, circles=sqrt(city_frequency$freq), bg="#93ceef", fg="#ffffff", lwd=0.5, inches=0.15, add=TRUE)

m <- get_map("France", zoom = 5, maptype = "toner")
ggmap(m) + geom_point(data = city_frequency, aes(lon = lon, lat = lat, size = sqrt(freq)), alpha = 0.7) + labs(x = "longitude", y = "latitude", title = "City Frequency")


The thing is that in the second map all the points are located a little bit North from the right place. How can I solve this problem? Thanks in advance.


Note: in order to get the coordinates of my data frame I used the geocode() in a early step.


---> I have created a public repository in github in order that everyone can access my data frames. Feel free to play with them.





Fix a repaired issue (gpx bug in 2.6)


sorry for problaly a stupid question but how do i use a cpp file to repair my qgis?


im talking in this case of the issue in this link: http://hub.qgis.org/projects/quantum-gis/repository/revisions/11c2bae74079ab008d7e7fb9233d714409aff008


can anyone help?





How to create this raster?


i got three layer 1-building data(polygon) contains park lot capacity in it 2-on the road parking data(line) also contains parking space data 3-car parks (polygon) includes parking space capasity i would like to create a raster, the raster must bu in 50x50 m resolution, in every 50m square i want to know the parking capacity in it, how to i divide my city to pixels and calculate this parking space capacity in the pixels





Make ldouble list in qgis


I want to create a table whith two drop-down lists.


If I select somthing from first drop-down list in the second will change attributes that I can choos. For example, if in the firstlist I pick a house, in te second I can choose -stone or wood. And if I choose in the first- tree than in the second list I can choose from tree or palm.





Create Random Points with a specific table


That's my question.


I have 3 tables in a Database: The first table is joined with another table (TABLE 3) where are stored 4 or more pairs of coordinates (so for one record of TABLE 1 there are more correspondences of record in TABLE 3). So, the TABLE 1 stores records that are represented like polygons in QGIS [For example ID 1(TABLE 1) = x1,y1; x2,y2; x3,y3; x4,y4 (TABLE 3)]. In the TABLE 2 there are records that are represented in QGIS like points linked with the TABLE 1 through an unique ID. I need the coordinates of the "TABLE 2" to calculated in a range between the coordinates of the "Table 3" (even randomly). The result in QGIS should be to random points inside a polygon.


I know that there is a command in fTools that randomly distributes the points, but the problem is that those random points that creates the program I do not know how to join them with the "Table 2". I hope I have expressed the problem correctly, because my English sucks :)





vendredi 28 novembre 2014

Modeling SQL Query to add statistical data to a geometry


Hello I have a problem modeling a SQL Query for the following problem:


I have two tables with geometries from Shapefile der Verwaltungsgrenzen (WGS84). Of intrest for me are Kreise and Gemeinden, whereas a Gemeinde is always part of the bigger Kreis. Each Kreis has a unique id called rs which is exactly 5 chars long, where as the unique id of Gemeinden is 12 chars long and the first 5 chars are equal to a Kreis id.


For some Gemeinden in de_commuter_gemeinden I have statistical data as well as for most of the Kreise in de_commuter_kreise. These data describes how many commuter there are.


For each geometry with statistical data I have to create random points within the geometry. For the Gemeinden it is no problem, since the are a subset of a Kreis and I can select the geometry from the table de_shp_gemeinden.geom. But for the Kreise I have to substract the already generated points as well as the geometry they were created in. My current SQL for this is:



SELECT
k.rs, k.gen, st_difference(k.geom, (
SELECT ST_Union(geom) AS geom
FROM de_shp_gemeinden g, de_commuter_gemeinden c
WHERE c.rs ~ ('^' || k.rs) AND c.rs = g.rs)
) AS geom
FROM de_shp_kreise k


Which looks like this


Kreise without Gemeinden


Now I need a smart way on calculating the statistical data for each Kreis, meaning get the data of commuters for the Kreis and substracting all the data from it subset Gemeinden.


Could anyone point me in the right directions? If you need more info please ask. Thank you.





Coordinates on Allotment map


I have an allotment map with angles and distances but no North/East coordinates. instead it has at one corner some numbers wich I suspect are the coordinates but in some different form. Please help me to convert em in N/E coordinates. The numbers are;


U 19419.157 B 22999.288





QGIS georeferencer - GCP projection problem


In short, the problem is that when Geoereferencer is reopened with a previously edited raster, all the previously saved control points have moved off the image to a location thousands of kilometres away. (QGIS 2.6, Windows 64 bit).


Here is what I did. The project CRS is set to GDA94/MGA zone 55 (EPSG 28355), OTF projection is ‘on’ (to project some GDA94 lat/lon layers to the project’s MGA zone 55 tranverse Mercator projection). I open Georeferencer, and choose a raster file to georeference. My default is to prompt for a projection, so I choose MGA zone 55.


In Georeferencer, I create a number of GCPs (i.e. control points) by referencing to objects in the map canvas. After adding points, I start the georeferencing calculation, choosing Thin Plate Spine (cubic) for full rubber-sheeting. Operation is successful, the raster appears on the map canvas, along with the control points. The control points are also visible in the Georeferencer window, still overlying the points chosen for georeferencing. I save the control points and close the Georeferencer. All good so far.


Now for the problem. Having closed the Georeferencer, I reopen it, and open the same raster image. But now, the control points in the Georeferencer window no longer overlay the raster in the Georeferencer window. However, they do appear correctly on the map canvas (the main QGIS screen). In the Georeferencer window, they have in fact displayed, but half a world away to the south. That is, if the raster covers a rectangle with coordinates 300,000, 5834300 ; 362100, 5782700, the control points occupy an area bounded by coordinates 0,0 ; 6800, -5900. These small coordinates seem to correspond to the number of pixels in the raster image. So, it is no longer practically possible to edit the GCPs.


How do I solve the problem?


See screenshots below.


Further note that may or may not be relevant: even though the GCPs no longer overlie the features on the raster image, a further run of the Georeferencer will still produce a correctly positioned and rubber-sheeted image in the main map canvas.


Image 1: Georeferencer window, and map canvas. GCPs do not show in window, but show on canvas. Image 2: Zooming out in Georeferencer window shows the GCPs are there but far away


Image 1 Image 2





Generate delaunay triangulation from 3D(x,y,z) points


I have 3D Points array.each 3D point have x,y,z data. I want to generate delaunay triangulation.I am new in delaunay triangulation. My array as below.

-256.24 85.35 5.0

-262.24 84.35 15.0

-263.24 87.35 25.0

................

I am working in c++ with Visual studio IDE. I want to generate indexes of triangulation. As well as array of triangles points.





Labels with buffers creates bloated PDF file in Composer


I know this question has appeared before, but the earlier questions appear to be linked to the old labeling engine. I'm running qgis 2.0, which I thought uses the new-generation labeling engine. In my map there are three labels (not three layers, but just three labels) that are buffered. The pdf file grows from 1MB to 18MB by turning on buffering.


Am I wrong in thinking that 2.0 uses the new labeling engine? The rest of the team I'm working with is on 2.0 and moving everyone is not easy for a variety of reasons.





image size has been increased after i.atcorr


I have compiled i.atcorr and run for different data sets. The problem is, it is writing the output 4 times of input image. I checked for data types, they are fine. I am not able to catch where exactly problem is. Can anyone suggest me?


I also want to mention that image first processed for toar using grass gui 6.4.4 and then the image passed for i.atcorr (customized for our purpose) actually in toar itself it is getting doubled. Is there any problem in gui version of i.landsat.toar? I have used for landsat5 tm data path141row43. Please suggest me.





GGIS: what line styles work best with Garmin maps


I am using QGIS to create some maps to load on to a GPS using the plugin. When I display the maps on the GPS the lines are very faint and I can not distinguish colours and some colours vanish altogether. But even the black lines don't stand out.


I have tried increasing the thickness of the lines in QGIS but this does not seem to make any difference.


I am displaying the maps on an ETrex 30.





Draw line overlays and center map


(Apologies if this is a trivial question; I'm completely new to the GIS world.)


I have a dataset with start and end points of lines. Columns in the table are



unique_line_id, line_start_point_lat, line_start_point_lon, line_end_point_lat, line_end_point_lon


I need to use this in my website (python,mysql) to draw lines as map overlays connecting these points.



  • According to the user selection, the map will display one or more lines from the above table

  • I need to zoom the map to fit the line/lines in the view when the page is loaded, and center the map accordingly

  • I would like to add a small thumbnail image in the middle of each line


I am looking for open source technologies to use in my website which would not need a paid license. (I have a vague idea to use google maps as base map images and use google JS API or D3.js, I've also seen references to GeoDjango) Please let me know any demo sites too, with one or more of these implemented. Thanks In advance!





Cross information


I'm a total amateur at ARCGIS, this is the first time I use it for real. The thing is I have several shapefiles with information about:


a. Schools in my city (elementary, middle and high school) b. City districts


The schools' information has Points that represents the unit (the place itself) and also a Buffer (influence's radio of every unit). ONE LAYER BUT SEPARATE SHAPES (One for point one for buffers). So, I need to cross this information so I can generate data that tells me the districts that are best equipped with each kind of educational level. I mean, which districts have more coverage of elementary schools, more coverage of middle high schools and more coverage of high schools. Even if they don't have any school within, but just some coverage from any buffer.


I don't know if I explained myself very well, but I'll be so thankful for any advice. I think I must merge the schools and then join this shape with the districts... (From what I've learned in my two days ARCGIS experience)





How to convert proprietary ESRI geodatabase to open standards and consume data in a web service?


I have a ESRI file in geodatabase format (v9.3). I intend to use this in a non-commercial application. I think I will have to convert it into a open standard since I am confused if ArcGIS online would need a paid subscription even for a non-commercial 'production' application. (Please let me know if i am wrong about their developer license)


So I am here basically looking for a solution to convert the gdb file/directory to another open format and do below operation: I need to use the gdb data (available data is polygon map overlays) to derive connecting lines from mid point of one polygon to the other polygon. (there would be a huge number of sets of two polygons like this). The end dataset would be a single table having below columns



set_id, line_start_point_lat, line_start_point_lon, line_end_point_lat, line_end_point_lon


and I plan to use this data further in any open standard (google maps/ d3.js - still need to figure out); to use in my website to draw lines as map overlays connecting these points.


Could anyone tell me what open source technologies I should be looking at to convert and do above operations?





QGIS FOR IOS, IPhone or IPad


Why is there not an IOS version of QGIS? The overall file size is not excessive. Given the alternatives, looks like a niche QGIS could well prosper within.





osm2pgsql cannot connect to postgis database, running on a Bitnami project installer


Based on this guide I am trying to import an .osm.pbf file into my PostGIS database, which is running on my Postgres Bitnami installation here:



[bitnami@ip-xxx]/usr/share/proj$ which postgres
/opt/bitnami/postgresql/bin/postgres


but I am getting an error, where osm2pgsql seems unable to connect to my database:



[bitnami@xxx]/home/.../geoserver_data/OSM$ osm2pgsql -E 900913 -d geodb -U geouser -W -S /usr/share/osm2pgsql/default.style washington-latest.osm.pbf
osm2pgsql SVN version 0.82.0 (64bit id space)

Password:
Error: Connection to database failed: could not connect to server: No such file or directory
Is the server running locally and accepting
connections on Unix domain socket "/var/run/postgresql/.s.PGSQL.5432"?


Yet my Postgres install is listening to port 5432 under /opt/bitnami/postgresql/.s.PGSQL.5432 as follows:


[bitnami@ip-xxx]/home/.../postgresql/share$ sudo netstat -nlp | grep 5432 tcp 0 0 127.0.0.1:5432 0.0.0.0:* LISTEN 13969/postgres

unix 2 [ ACC ] STREAM LISTENING 1876491 13969/postgres /opt/bitnami/postgresql/.s.PGSQL.5432


It seems that osm2pgsql is having trouble finding my database. I know that I have my postgres server running, as shown here:



[bitnami@ip-xxx]/usr/share/proj$ ps aux | grep postgres
postgres 982 0.0 0.0 247712 468 ? S Oct21 0:17 /usr/lib/postgresql/9.3/bin/postgres -D /var/lib/postgresql/9.3/main -c config_file=/etc/postgresql/9.3/main/postgresql.conf
postgres 986 0.0 0.6 247856 7028 ? Ss Oct21 0:00 postgres: checkpointer process
postgres 987 0.0 0.0 247712 140 ? Ss Oct21 0:14 postgres: writer process
postgres 988 0.0 0.0 247712 104 ? Ss Oct21 0:12 postgres: wal writer process
postgres 989 0.0 0.1 248568 1088 ? Ss Oct21 0:12 postgres: autovacuum launcher process
postgres 990 0.0 0.0 103520 536 ? Ss Oct21 0:22 postgres: stats collector process
postgres 11111 0.0 1.1 171996 11876 ? S 04:28 0:00 /opt/bitnami/postgresql/bin/postgres -D /opt/bitnami/postgresql/data
postgres 11113 0.0 0.0 171996 876 ? Ss 04:28 0:00 postgres: checkpointer process
postgres 11114 0.0 0.1 171996 1696 ? Ss 04:28 0:00 postgres: writer process
postgres 11115 0.0 0.0 171996 932 ? Ss 04:28 0:00 postgres: wal writer process
postgres 11116 0.0 0.2 172844 2084 ? Ss 04:28 0:00 postgres: autovacuum launcher process
postgres 11117 0.0 0.1 27812 1120 ? Ss 04:28 0:00 postgres: stats collector process
bitnami 13358 0.0 0.0 11752 924 pts/6 S+ 09:10 0:00 grep --color=auto postgres
root 26079 0.0 0.0 64536 32 pts/7 S Nov26 0:00 su postgres
postgres 26080 0.0 0.0 21100 200 pts/7 S Nov26 0:00 bash


I also have the following under my /etc/postgresql/9.3/main/pg_hba.conf, /opt/bitnami/postgres/share/pg_hba.conf and my /opt/bitnami/postgres/data/pg_hba.conf, which seems like is should allow access for my geouser:



# TYPE DATABASE USER ADDRESS METHOD
local geodb geouser md5


Any suggestions about how I can import my OSM data?





Spatialite in GeoDjango - "...they do not properly define db_type..."


I've started a new project and choose GeoDjango for it. Because I already developed with Django and sqlite as a DBMS and all worked perfectly fine for me, I found it a great option to use spatialite to make the first development steps and perhaps switch to PostGIS later on, when the app needs scaling.


Although the docs tell different I faced a couple of problems now bringing the spatialite datastore up to life.


First I created a new empty sqlite database file on the sqlite command line (because let django do that for me with syncdb resulted in an error where syncdb didn't hae the geometry_columns table to register my spatial table...so it's about creating the spatial repository first).


Second I ran spatialite database.db "SELECT InitSpatialMetaData();" to create the repository tables.


At least I ran python manage.py syncdb to create my model tables. It does a couple of thing, but ends up in the following error:



Operations to perform:
Synchronize unmigrated apps: gis
Apply all migrations: admin, contenttypes, walmap, auth, sessions
Synchronizing apps without migrations:
Creating tables...
Installing custom SQL...
Installing indexes...
Running migrations:
Applying contenttypes.0001_initial... OK
Applying auth.0001_initial... OK
Applying admin.0001_initial... OK
Applying sessions.0001_initial... OK
Applying walmap.0001_initial... OK
Applying walmap.0002_campsite_capacity...AddGeometryColumn: "duplicate column name: point"
OK
Applying walmap.0003_auto_20141126_1434...Traceback (most recent call last):
File "manage.py", line 10, in <module>
execute_from_command_line(sys.argv)
:
:
File "/home/juergen/Development/wal/env/local/lib/python2.7/site-packages/django/db/backends/schema.py", line 460, in alter_field
new_field,
ValueError: Cannot alter field walmap.Campsite.point into walmap.Campsite.point - they do not properly define db_type (are you using PostGIS 1.5 or badly-written custom fields?)


The model class just consists of an id field, a name, a slug, a Point Geometry and an int field for capacity, so nothing special here.


Any ideas about this issue?


EDIT: SpatialLite is version 4.1.1


PROJ.4 is version 4.8.0


GEOS is version 3.4.2


Django version 1.7


python version 2.7.6