Darb Dabney's 3d GIS – OpenSim

Thanks to the astronauts aboard the International Space Station, their time-lapse photography at very high ISO that helps to share some of what their eyes may well see, and of course Michael Koenig for his care and smoothing of the HD video, with some loungy score, too.
Take five (minutes) and watch it on HD in a darkened room. You might find yourself pausing, reviewing, and spending 20 minutes enjoying.

Earth | Time Lapse View from Space, Fly Over | NASA, ISS from Michael König on Vimeo.

I was fascinated by an orange wiggle, that turned out to be the astoundingly well-lit India-Pakistan border, around 1000 km long.

Meanwhile, I’m forming some plans for next semester’s course, and have realized that it may well be possible to offer students training in multi-user virtual environments without hacking one of the lab workstations to image it as an Opensim server. Thanks to the incessant business analytics of Maria Korolov over the past few years, it was possible for me to quickly get caught up in the new and improved options for cloud hosting of Opensim regions.

Right away it became clear that the business model of Kitely was quite compatible with my modest but area-expansive needs for real-life terrain simulations.  I’ve found it quite easy to get set up with a single region, and that’s a really big start.  I was able to use the latest beta Second Life 3.2 viewer to connect to the latest Opensim 0.7.2 stable release, tweak terrain and set up a few flexi-prims to test the weather.  Nice work technically, and a very nice pricing scheme for my sort of use.  I’m also very sympathetic to Ilan Tochner’s philosophy of “just keep building new regions”—it’s a consistent theme with cloud solutions, and refreshing to see it in connection with Opensim.

The Marin Community Map cache has been draft-built, and can be explored as the highest four zoom settings here

I’m settling in to a new GIS workstation, custom-configured HP Z400 with twelve execution threads.  Windows 7 configures with a  Windows Experience Index of 7.3 and running very new browsers is a treat.  The preview version of Internet Explorer 10 checked out thusly with its Windows-loving HTML5 Fish Bowl Speed Test.  Two of them running, each with 2000 fish, both hitting 31 FPS.

Fast HTML5 Browsing, according to Microsoft

I I think that this will meet my needs for now. ;^)

The CPU cycles are burning as I set up some test cache building for our local edition of ESRI Community Maps – Large-Scale Topographic Base Map.  It’s rather a treat to see a bunch of execution threads and trying my best to keep them all busy!

As I work through glitches here and there, I’m catching up on some of the good information that is online to help folks in my situation, like Advanced Map Caching from the 2010 ESRI Developer Conference.

The enormously tragic megathrust event east of Honshu island in Japan has been in my thoughts throughout the weekend.   It is natural to focus on the human loss and images of civic destruction; I’ve got little to add to that story.

When I saw the first news items, it was already Friday morning and eight hours past the event.  My concerns were for Hawaii and the arrival of a tsunami.  I was shocked by the energy released, as reported at NEIC:  Mo=3.9×10^22 Nm,  Mw=9.0 — for as many earthquakes as have been recorded and studied in Japan, this was a much larger event.  For disaster planning, that is Not A Good Thing.

The beach ball diagram shows compression back up the thrust. Depth of event is 24km

While there was some concern for Hawaii, there is something to be said for the diffraction caused by the long line of seamounts extending northwesterly from the big island.  While the energy won’t be cancelled out, the height of the crest would be widened some.  By the time the tsunami reached California, it was a real concern and a literal wake-up call for emergency workers.  Our damage here was mostly an economic annoyance, with Del Norte county once again taking the brunt of the damage in Crescent City harbor.

With a check-in call to relatives in Hawaii saying everything OK, their eyes and mine turned towards Japan.  As I watched shocking videos of tsunami damage, I struggled to calibrate what I was seeing with my mental model of how fragile most human habitation structures are.  I’ve thought much about the effects of shaking, liquefaction, and occasionally about wildland fire.  I’ve read about losses in Marin county during 1982 flooding where debris flows destroyed houses or literally rolled them end-over-end down a slope.  I’d seen some videos from near Sumatra of debris flowing up streets after waves climbed up the beach.  Yet images from Japan recorded damage unfolding at an entirely worse scale.

My concern became much more engaged when I saw helicopter video of the south end of nuclear power reactor facility Fukushima Daiichi (No. 1).  You see, in the video clip I recognized not the reactor, but the cut slope behind it.  In a very strange sense of model-based deja vu, my memories were unequivocal: “I know that place!”.

Why?  Because about 17 or 18 years ago, I was fortunate to be part of a site design project for the (still yet to be constructed) Unit 7 and Unit 8 reactors.  Fortunate for me, because it was a first opportunity to learn not just ESRI Arc/INFO, but how to work with Bentley CAD software to create a 3D cut slope design.  My small contribution was to create a realization of slope that was  not simply faceted as a buffer surrounding the building footprints, but instead create a semblance of the natural cliff slopes  adjacent to the plant, while meeting the engineering requirements for not-to-exceed slope steepness, and a more natural-looking accumulation of drainage from the slope.

So apparently, after spending a couple of weeks learning to use 3D design software and creating a design that extended the existing cut farther southward, I had an image of the plant, its cliffs, and the breakwater that guided cooling water that stuck with me.  After having the flash of recognition with the video, I opened up Google Maps and found the site on my first inward zoom.  It was a bit spooky.

So now when I watch coverage of hydrogen venting that leads to building explosions, I feel a curious terrain-based sense of connection with the site.  I wish them well, and the safest possible return to production.  The TEPCO power is needed by many people.

And more ominous for the Pacific Northwest, I can’t help but reflect on what a similar megathrust event would mean for Cascadia.  Both Portland and Seattle would be in some sort of peril, although I don’t have a clear understanding of how tsunamis are modeled for either lower Columbia River or for Puget Sound.  But the possibility of a 9.0 megathrust event along the Cascadia subduction zone was a whole lot more abstract for me until last Friday.  It may not be the best time to reflect on it now while Japan is suffering—but the risk to the northwest has been a matter of public record for several years now and it not the time to forget that, either.  An event of that size and location would have tsunami implications both locally, and perhaps a far greater risk for windward Hawaii.

—————-<>—————-

Closer to home,  this past week we’ve created a prototype 3D product that provides a facsimile of our 45cm terrain model—imported to Google SketchUp 8 with a georeferenced orthophoto texture on the terrain.  Things are not fully tuned up yet, but we are able to take TIN and decimated TIN surfaces out of ArcGIS 10, by way of clipping polygons that are interpolated to multipatch (3D multi-polygon features) and then exported to Collada with a KML point for georeference.  The Collada can be imported to SketchUp, at least up so some limit of detail.  Once there, I’m trying now to find the right way to smooth the facets and improve the rendering of the surface without having high-contrast dark facets—because the orthophoto textures are arriving intact and draping over the surface.  Only at the moment, many black facets are covering up almost half of the orthophoto in an aggravating triangular patchwork.

Next step: smoothing within SketchUp.

And inspiration from others who have gotten gridded terrain into SketchUp.

Each year for the past 40 or so years, the American Geophysical Union has met in San Francisco around December for what is now the Fall Meeting.  I’ve been an AGU member for about 28 years, and for a time was attending each and every Fall meeting—but these days it’s about once every three years.

This was one of those years, and it was a great pleasure today running into a good handful of friends and former school colleagues from years past!  Also, it was much fun to present a poster that summarized an analysis of synthetic flow lines built on the integrated topographic-bathymetric surface model.  Basically, with a very detailed 3D surface grid that runs continuously from mountaintop to offshore out to the 3-nautical-mile legal boundary of California counties, it is possible to draw streams as they would have flowed when sea level was lower, like 7000 years ago.

Much of the interesting topography from those streams got clobbered by sea level rise.  As the Ice Age retreated and continental glaciers melted out, the waves from the Pacific Ocean pounded the coast back to where it is today, and planed off much of where the streams once ran.

With ArcHydro-style drainage analysis on our terrain model that has fused detailed multibeam bathymetry from the California Seafloor Mapping Project, it is possible to identify extremely subtle signatures in the portion of the offshore platform that is Santa Cruz mudstone formation, a harder Miocene formation that expresses bedding in its surface.

With the analysis, when synthetic drainage paths are symbolized to emphasize flow lines with greater catchment area one can observe suggestions of right-lateral offset.  In California, this is a signature pattern for tectonic offset of drainages that cross strike-slip faults with right-lateral offset.  Because the formation where the analysis has detected possible offset is older (Miocene is more than 5 million years old, but the offset is perhaps only in the last 1 million years), this result should not cause much excitement with regard to modern seismic hazard.  It could however prove helpful to those who would decode the geologic structure of the Point Reyes peninsula.

(comments on the Cr-48 have been moved to its own blog at http://cr.bq.bz

Sounds great!  Less work!  And yet, frequently, it often doesn’t work because valid SDE geometry must still be made OGC-valid.  Three days ago I found myself with a working syntax to apply the (.NET CLR invocation of) MakeValid() method to my existing GIS data, and so I had valid OGC spatial tables available.  However, these OGC tables were missing a couple of things that needed to earn their registration back into the SDE fold.  Things needed, or just good-to-have when invoking sdelayer to register an OGC table include:

1. A defined primary key field with its own clustered index
2. OGC Spatial indexes (up to four levels available in OGC tables, and three in ArcSDE), with a name to pass along to SDE
3. A numeric expression of the geometry’s extent, as bounding rectangle limits
4. A pointer to the geometry’s projection “srid”  in the SDE_spatial_references table, linking it to real SRID as “auth_srid”
5. A friendly 63-character name to describe the spatial features
6. A known instance of ArcSDE, such as port ‘5151‘ or DirectConnect ‘sde:sqlserver:servername‘
7. A resolvable name or IP address for the ArcSDE server
8. The name of the ArcSDE database on the server
9. Login credentials for the database under the ArcSDE instance, most useful with CREATE privileges

One last detail: you’ll likely need to have the administrative parts of ArcSDE installed on your workstation.  This means installing ArcSDE 10, but not configuring it to run as a server locally.  Mercifully, this type of install does not require access to an ArcGIS Server license.  Of course, you’ll use it to perform actions a server instance where that license was required.  But after that install, your workstation should be able to run command line administrative tasks on a remote SDE server, often with either the ArcSDE monitor ‘sdemon’ or the Feature Class utilities of ‘sdelayer’, each of which has very many options.

1) Update the OGC table to ensure you have the correct SRID on every Geometry feature

UPDATE [SpatialTesting].user.PARCEL_OGC_LI
	SET Shape.STSrid = 2872

2) Define the primary key in your OGC table using SQL Server Management Studio (SSMS), and doing so will create a cluster index for you. In SSMS 2008 R2 this could mean a right-click on the OGC table name, selecting the Design view (second choice), then highlighting the OBJECTID attribute row by clicking its box along the left, right-click and choose Set Primary Key.  Save the change, or close the view and accept to save the change.  This should provide the “gold key” icon for OBJECTID as a Primary key, add an entry in the table’s Keys folder, and list an index named like “PK_your_table name  (Clustered)” in the Indexes folder.
Please note that you can’t skip this step and succeed in building a spatial index.   The spatial index appears to require a clustered index on a defined Primary Key.

3) Create the spatial indexes for the OGC table by right-clicking the table’s Indexes folder and choosing “New Index…” option.
With the General page chosen in “Select a page”, give the index a name (e.g. “spx_50″) and choose Index type: Spatial, then click “Add…”
In the “Select Columns from” popup, hopefully you will only  have one spatial column’s row to choose from (e.g. Shape).
Check its box and click OK.
Then in “Select a page”, choose Spatial, and here you will enter bounding coordinates for the OGC Geometry extents.
(e.g. 5816131, 2125631, 6028424, 2312384 as X-min, Y-min, X-max, Y-max)
For Geometry, ensure that the “Geometry grid” Tessellation Scheme remains chosen, and 16 Cells Per Object seems OK thus far.
While SQL Server 2008 offers four levels of spatial index grids, it only has qualitative labels for its options; by contrast ESRI file or SDE geodatabases have three grid levels, but one can specify precise index grid scales.
Given the prompts in the interface, I have typically understood Level 1 / top-level grid to be the coarsest, and chosen “Low”, then at Level 2 chosen “Medium”, and for both Level 3 and Level 4, chosen “High”.   Once you’ve chosen, click OK and wait while the index is built.

4) Once the spatial index has been built on the OGC spatial table, it should be ready to register with SDE and return it to ArcGIS service. The example below (less credentials) is what worked for us with line feature class with over 275,000 features.
Just open up a ‘CMD’ command-prompt window  and invoke the ‘sdelayer’ command with options along these lines:

sdelayer -o register -l PARCEL_OGC_LI,SHAPE -e l -t GEOMETRY -C OBJECTID,SDE
 -E 5816131,2125631,6028424,2312384 -R 2 -S "OGC Version of 2010 10 Parcels"
 -i sde:sqlserver:sqlgisdev -s SQLGISDEV -D SpatialTesting -u <user> -p <password>

And that should do it.

Another bit of syntax that I feel more comfortable with is making use of “UNION ALL” to chain queries together, so that they are easier to visualize in SSMS as overlays.

select Shape.STCentroid().STBuffer(10000) from SpatialTesting.bquinn.COUNTYBNDY_PG_OGC
UNION ALL
select Shape.STExteriorRing() from SpatialTesting.bquinn.COUNTYBNDY_PG_OGC

Takes our county boundary polygon, finds its centroid, buffers that centroid point out 10k feet, and then overlays that centroid blob on the outline of the county polygon, garnered from the poly’s exterior ring. Sounds simple enough, and the pasted code worked.  But 10 days ago I would have been very hard pressed to come up with that syntax!

The rules for geometry validation are just different between the ESRI world and the OGC world as instantiated by Microsoft in SQL Server 2008 R2.
I’ve been pounding my head against the wall over what I thought were the challenges of taking fine examples from Alastair Aitchison (a.k.a. tanoshimi on the MSDN forums) that are presented in Beginning Spatial with SQL Server 2008, and trying to make these happen in useful ways on real GIS data.  That is to say, data that were not created by hand-entry or even well-known-text (WKT) entry into SSMS, which are often  simple single features.

I was particularly concerned reading the first paragraph in Chapter 11: “All of the techniques discussed in this chapter apply only to a single item of geography or geometry data, and generally do not require any parameters.”  As I understand now, Mr. Aitchison meant “All of the techniques as discussed…” and all the various Spatio-Temporal methods like STIntersects() work perfectly well on perfectly valid spatial tables.

I thought that these methods just weren’t going to work on my GIS data,  typically large tables with lots of objects, not the scalar objects of the book examples.  But as it turns out, I had a different problem.  My geometry, as transcribed via Shapefiles and the Shape2SQL tool from SharpGIS.net, was arriving in a form that was occasionally invalid for certain features.  One invalid feature in a spatial table makes the entire SQL statement crash.  Errors were reported with reasons that did not help me find my actual problem very easily.  Bad geometry meant that I needed to clean it up by running the MakeValid() method.

After a couple of hours goofing around with T-SQL syntax, I’ve made my first baby steps toward cleaning up my stuff.  Here’s how that is happening:

(select ID, PARCEL, Prop_ID, Shape_area, Shape_len,
(CASE Shape.STIsValid() WHEN 1 THEN Shape ELSE Shape.MakeValid() END) as Shape
INTO spatial_01.dbo.Parcel2 FROM spatial_01.dbo.parcel)

As so often with code, it looks stupid-simple to see it written out in the way that actually works; a major “duh” factor.
The pay-off is that with 100% valid geometry, table vs. table spatial queries are now working.  I’ve not yet broken any speed records in query performance vs. ArcGIS 10, and yet I have finally proven to myself that I can use  Spatial SQL to do powerful geoprocessing that, in some cases, will be easier to write and maintain in SQL than in ArcGIS / Model Builder / Python.

For context, I’ve struggled through the past four days trying to figure out where my Spatial SQL syntax was messed up.  As it turned out, exactly ZERO of my uploaded spatial data tables arrived with 100% valid GEOMETRY data.  Anyway, once I made everything OGC-valid, a query like this

SELECT NAME, b.PType INTO spatial_01.dbo.School_Geology
FROM spatial_01.dbo.school2 a, spatial_01.dbo.geology2 b
WHERE a.Shape.STIntersects(b.Shape) = 1

"NAME","PType"
"TOMALES HIGH SCHOOL","Twg"
"BOLINAS-STINSON SCHOOL (STINSON)","Kfs"
"LAGUNA JOINT ELEM. SCHOOL","Qal"
"SAN GERONIMO VALLEY ELEM. SCHOOL","fsr"
"LAGUNITAS ELEM. SCHOOL","fsr"
"MARIN HORIZON SCHOOL","fsr"
"COUNTY COMMUNITY SCHOOL","Qal"
"ST RITA SCHOOL","Qal"
"ST PATRICK SCHOOL","fsr"
...

On a last note, this week I have validated the experience of using ESRI ArcSDE 9.3.1 settings to configure SDE (when running atop MS SQL Server 2008) to store its geometry in native OGC GEOMETRY binary form rather than the classic SDEBINARY form.  To do this, someone with db owner privileges can update the SDE_dbtune table with these directions to cause the GEOMETRY format to be used by default.  In principle, this should allow one to exploit the live feature classes in SDE for Spatial SQL queries—without the muss and fuss of exporting through Shapefiles and Shape2SQL (which is not intended for production-environment use, anyway), or need to code something in .NET to export a copy.

Wouldn’t it be nice to have ArcSDE running just like it always has, and then also be able to use Spatial SQL queries on the exact same feature classes that are being used throughout the enterprise?  Sounds great, in principle.  After tonight’s insight into the need to MakeValid() some perfectly functional feature classes, I do have some concerns that SDE feature classes could get broken by MakeValid().  On the other hand, it appears that one can take a Spatial SQL GEOMETRY-based table that meets SDE limitations: (1) only one geometry column, (2) all rows of the table must have the same geometry type, (3) all rows of the table must have the same Spatial Reference ID (SRID) defined—and any SDE feature class that was cleaned through MakeValid() should retain those characteristics—and then register the resulting spatial SQL table once again as an SDE feature class.  Hey, it sounds like it’s worth a try!

More working through Aitchison’s book has yielded many successful uses of OGC methods on our GIS data:

-- SELECT Shape.STGeometryType() FROM [spatial_01].[dbo].[road]
-- SELECT Shape.STGeometryType() FROM [spatial_01].[dbo].[Parcel]
-- SELECT Shape.STGeometryType() FROM [spatial_01].[dbo].[school_pt]
-- SELECT Shape.STGeometryType() FROM [spatial_01].[dbo].[geology]

-- SELECT Shape.STDimension() FROM [spatial_01].[dbo].[geology]
-- SELECT Shape.STDimension() FROM [spatial_01].[dbo].[road]
-- SELECT Shape.STDimension() FROM [spatial_01].[dbo].[school_pt]

-- SELECT Shape.InstanceOf('Surface') FROM [spatial_01].[dbo].[geology]
-- SELECT Shape.InstanceOf('Polygon') FROM [spatial_01].[dbo].[geology]
-- SELECT Shape.InstanceOf('Curve') FROM [spatial_01].[dbo].[road]
-- SELECT Shape.InstanceOf('LineString') FROM [spatial_01].[dbo].[road]
-- SELECT Shape.InstanceOf('Point') FROM [spatial_01].[dbo].[school_pt]

-- SELECT Shape.STIsSimple() FROM [spatial_01].[dbo].[road]
-- SELECT Shape.STIsSimple() FROM [spatial_01].[dbo].[BldgFoot]

-- SELECT Shape.STIsClosed() FROM [spatial_01].[dbo].[BldgFoot]

-- SELECT Shape.STIsRing() FROM [spatial_01].[dbo].[road]

-- SELECT Shape.STNumPoints() FROM [spatial_01].[dbo].[road]
-- SELECT Shape.STNumPoints() FROM [spatial_01].[dbo].[geology]
-- SELECT Shape.STNumPoints() FROM [spatial_01].[dbo].[school_pt]

-- Find Centroids and blob them out big
-- SELECT Shape.STCentroid().STBuffer(3500) FROM [spatial_01].[dbo].[City]

-- Calculate Perimeter Lengths and sort them out
-- SELECT Shape.STLength() as Shape_len FROM [spatial_01].[dbo].[City] ORDER BY Shape_len

-- Calculate City areas and sort them out
-- SELECT NAME, CAST (Shape.STArea()/1000000 AS INTEGER) as Shape_km2
--   FROM [spatial_01].[dbo].[City] ORDER BY Shape_len

-- Reading the SRID
-- SELECT Shape.STSrid FROM [spatial_01].[dbo].[City]

-- Setting the SRID
-- UPDATE [spatial_01].[dbo].[City] SET Shape.STSrid = 32610
-- SELECT Shape.STSrid FROM [spatial_01].[dbo].[City]

-- oops, that's wrong, let's set it back to 2768
-- UPDATE [spatial_01].[dbo].[City] SET Shape.STSrid = 2768
-- SELECT Shape.STSrid FROM [spatial_01].[dbo].[City]

-- counting elements
-- SELECT SUM(Shape.STNumGeometries()) FROM [spatial_01].[dbo].[City]
-- SELECT SUM(Shape.STNumGeometries()) FROM [spatial_01].[dbo].[school_pt]
-- SELECT SUM(Shape.STNumGeometries()) FROM [spatial_01].[dbo].[geology]
-- SELECT SUM(Shape.STNumGeometries()) FROM [spatial_01].[dbo].[road]
-- SELECT SUM(Shape.STNumGeometries()) FROM [spatial_01].[dbo].[Parcel]
-- SELECT SUM(Shape.STNumGeometries()) FROM [spatial_01].[dbo].[BldgFoot]

-- Geometric Difference between cities and their own 3500-meter centroid blobs
-- SELECT Shape.STDifference(Shape.STCentroid().STBuffer(3500)) FROM [spatial_01].[dbo].[City]

-- Geometric Symmetric Difference between cities and their own 3500-meter centroid blobs
-- SELECT Shape.STSymDifference(Shape.STCentroid().STBuffer(3500)) FROM [spatial_01].[dbo].[City]

-- Creating 50-foot Buffer around roads
-- SELECT Shape.STBuffer(50) FROM [spatial_01].[dbo].[road]

-- Simplified Buffer around roads, tolerance 8 feet
-- SELECT Shape.BufferWithTolerance(50, 8, 'false') FROM [spatial_01].[dbo].[road]

-- Create convex hull around schools
--SELECT Shape.STConvexHull() FROM spatial_01.dbo.school_pt

-- buffer parcels
-- SELECT Shape.STBuffer(300) AS Shape INTO pcl_seed FROM parcel where Prop_ID = '001-001-01'
-- select * from pcl_seed  -- to verify the selection
-- SELECT Prop_ID, a.Shape FROM spatial_01.dbo.parcel a, spatial_01.dbo.pcl_seed b
--   WHERE a.Shape.STIntersects(b.Shape) = 1

Many interesting projects have been happening, so many that projects have been backed up while workstations grind through their days-long work flows. This creates opportunities to update systems while we wait, and so the past week has seen massive updating of Windows system stuff and deleting old development tools so that the decks were cleared for MS SQL Server 2008 R2 Management Studio. With a bit of shape2sql from SharpGIS and today we’ve had our first spatial tables loaded and fledgling spatial SQL queries made.  Oh, at the end of the day we got the third of five workstations updated to ArcGIS Desktop 10 and the full complement of seats moved over to ArcGIS 10 licensing.

What’s been holding things up has been some final processing of cartographic-grade contours, one of the key new feature classes developed for support of the ESRI Community Maps Program, where local jurisdictions grind their own cache tiles for large-scale topographic mapping.  Since it’s part of a worldwide seamless map, of course the contours and spot elevations must be metric—which meant generating new contours.  The bathymetry was done on one-meter interval, with some half-meter supplemental contours in shallower waters.  The topography was done at quarter-meter intervals up through 25 meters elevation, where we’d included all available LiDAR data, the half-meter intervals through 50 meters, and one meter up to the summits.  All non-integer meter contours were flagged as supplemental, and contour index attributes were also calculated for 2, 5, 10, 20, and 50-meter intervals.   Spot elevations were derived from VARGIS photogrammetric spot elevations that were screened down from about 75,000 points to just 440 points, using neighborhood focal statistics.

For performance’s sake, after index attributes were calculated, the contours were chopped into segments not-to-exceed 500 meters shape length.  When all the chopped segments from all the various elevation ranges had been merged to a single polyline feature class in a feature dataset in a file geodatabase, spatial indices were built for the cache tile scales about 1:1200, 1:4800, and 1:19000.  The purpose is to make the contours as fast as possible not just for web app interaction, but for rendering the large-scale cache tiles.  In the end, the contours in this one feature class have about 1.4 million polyline features that fill the file geodatabase to about 6.8 GB of data.

Meanwhile, some great brainstorming has taken place with regard to hydro-enforcement of the terrain model, so that accurate synthetic flow lines can be generated county-wide.  Based on prototype work performed by Evan K. Babb before summertime, we knew how to deal with our terrain features and manually generate hydro-enforcement features.  Now, with a funded project to do the work county-wide, we’ve needed to devise a spatial analysis technique using the approach described by Poppenga and others in 2009.  Based on some correspondence and personal communication with Poppenga, we have sketched a workflow that should automate the creation of most simple (single road-crossing) hydrologic enforcement features, especially in the absence of accurate or complete spatial features for culverts.

So with the contours completed, it should be possible to have some progress made on automation of hydro-enforcement, completing more of the data development for large-scale topographic base mapping, help my colleagues continue the ArcGIS 10 migration, and follow up with some demonstration spatial SQL queries.  It’s an exciting couple of weeks!

There’s also been a thread of thought over the past six years or so that supports the workflow of using Visual Studio to compose SQL queries, in preference to the direct SSMS (SQL Server Management Studio, Tool Formerly Known As “Enterprise Manager”) ways of composing a query, because the typed pattern matching (MS: Intellisense) pulls up all the various ST__  spatial methods that are available in SQL Server 2008 spatial, but saves one needing to look-it-up and type-it-in correctly.  Such productivity gains really help boot up the learning experience a bit faster!

Opensim has slipped to a back burner for these past few weeks, as I button up a new edition of county terrain–this one fused topography and bathymetry.
Another focus has been to compile the map data layers necessary to publish a local edition of the large-scale topographic map, using the ESRI Community Maps Program template. This (Opensim-distracting) work forms a synergy, as the terrain drives revised hillshade and topo-bathy contour lines, as well as forming the basis for refined surface flow and stream mapping, all of which show up in the large-scale topo map. Also, I’ve compiled a great deal of web research into a single feature layer: classified building footprints, many of which have use categories [commercial, medical, educational, religious, etc.] and specific organization names included, for every structure on the associated parcel.

Why generate this large-scale topographic base map? It’s purpose is to provide local control over local detail, in a cloud-based service suitable for desktop, tablet, and mobile phone users who may consume applications that mash or overlay the base map with relatively sparse point or line features. While we’ll be able to serve up the base map through our own MarinMap servers, one of the key benefits of the conforming large-scale topographic map style (as made possible by the template) is that the caches of the four largest-scale views will be rendered and contributed for inclusion in the ArcGIS Online map services, as viewed at arcgis.com and consumable by many interactive map sites. The cached map tiles are essentially raster prints of the map template (which has styles customized at four local zoom levels) that can be served extremely fast.

Those fast maps can also be efficient, for mobile phones and tablets consuming base maps with no special rendering. Along that path, I’ve rebuilt my mobile phone this past week, taking the T-Mobile G1 and rooting it, then installing various new images until I found one that was up-to-date, stable, and attractive to use. As of this weekend, I’m now running Chromatic 4.5 with SetCPU 3.0.2 overclocking. Chromatic 4.5 is a packaging of CyanogenMod 6.0.0 with the ADW launcher and many other tunes. CyanogenMod 6.0.0 is an Android Open Source Project-based build of Android version 2.2 “Froyo”. So now, when I view maps and browser content, I can pinch to zoom out, spread to zoom in, and still push around the map to pan it on the screen. As part of sharing content with mobile phones, I’ve learned how to generate QR tags to provide access to map URL or geolocation tags.

For fun today, we visited the Silicon Valley area to show a guest some of the sights. I hadn’t visited the Googleplex in a couple of years, and many of the adjacent buildings that I recall as having different companies in them are now part of the Google campus. We found one interesting sight, too–and it appears that today may have been Google Inc.’s 12th birthday. We didn’t see any partying in the mostly-empty parking lots, which wasn’t surprising for a Sunday!

Blemishes notwithstanding, nearly six months of back-burner work has reached a threshold of readiness and is outward bound to some engineering firms, flood mappers at FEMA, and interested parties within the county. A handful of known issues remain unresolved. Proper name is “tbsm45cm_20100823″, proper edition is “2010.08″.

This is the third version of the terrain. Second version was “2010.01″ and included multiple LiDAR data sets, but fewer than presently used, and was a topographic model only. First version was “2009.09″ and was mainly photogrammetry and FEMA LiDAR, and was the last version to be developed in California Coordinates. Once the massive NCALM LiDAR data sets were processed, it became easier to move everything into WGS84 UTM zone 10 north meters projection, WGS84 NAVD88 CONTUS Geoid 2003 vertical position.

The NOAA utility program VDatum, a brilliant Java-based application able to stream-process data sets of near-infinite size, brought the NCALM data to heel, and opened up decades of NOAA depth surveys to our use in integrated topographic-bathymetric surface modeling.

First-return NOAA ALACE LiDAR swaths were fused along the outer coast, as bare-earth filtered versions were not produced in 1997–2002; the benefits of LiDAR detail along the rocky coast do seem to outweigh the distracting appearance of structures near Rodeo Lagoon, Stinson Beach, and outboard Bolinas.

When ArcGIS 9.4 beta 2 reached its limit in ability to render the terrain dataset into 45cm grid over the full extent, the clipping quadrants created to resolve this problem ended up chopping a very small portion of Sonoma county that drains into Estero Americano; the full watershed remains intact in the 1-meter version of the terrain grid under analysis for county-wide hydrology. Likewise, the tighter clipping quadrants lost a few hundred meters of San Pablo Bay bathymetry just west of where Marin, Sonoma, Contra Costa, and Solano counties meet. Also, tighter clipping quadrants snipped a portion of the San Francisco Bar southerly of San Francisco’s Seal Rock that was intended to be part of the model. All of these areas exist in the 100cm grid, and will be part of drainage analysis.

Happily, we have updated the workstation to ArcGIS 10, and have been enjoying such great speed gains with Spatial Analyst that our ERDAS use has been noticably reduced. Finally, Spatial Analyst is often showing performance nearly on par with ERDAS. Thank goodness that the Raster Calculator survived the transition to version 10 ArcGIS!

Painfully, the existence of unutilized bathymetric data sets for upper broad-channel Corte Madera Creek and Bolinas Lagoon have been revealed this week. Hey, there’s already something to look forward to for the next build!

The new terrain is getting some immediate use in support of an effort to participate in ESRI’s Community Maps Program for large-scale topographic mapping. The Program provides a template geodatabase with 36 vector feature classes and two raster, into which local agencies may pour their data. Once tucked into a conforming schema, a template multi-scale map document is provided with 120 layers—30 at each of four large scales that correspond to Google Maps and Bing Maps projection and cache tiling schema. The difference is that the template document makes use of ESRI tools to allow much more local detail to be packed into a map designed with notably more sophisticated cartography than either Google or Bing maps now have. The Community Maps Program concept is that local agencies may publish their local detailed content in a fairly uniform style, while retaining a world-wide seamless context for their surrounding area.

Qualitatively, the effect is that, when viewing the ArcGIS.com topo map alongside either Google or Bing maps (on two monitors, with comparison made at the same scale), the ArcGIS.com map looks to be a larger scale. It isn’t, and I’ve measured the size of features to convince myself, but my mind insists that I’m zoomed farther in on the ArcGIS.com map for some reason. My guess is that it is a perceptual effect of the much greater amount of information that is cleanly displayed in the ArcGIS.com map versus the much sparser Google and Bing content at these large zoom levels. Try it out—it’s like a carto version of an optical illusion!

The 120 layers in the template large-scale topographic base map from the ESRI Community Maps Program are arranged to provide four precise cartographic designs for Google/Bing map cache levels 16 through 19, which correspond to these display scales
1:15000–1:6001 (level 16, a.k.a. ~9k)
1:6000–1:3501 (level 17, a.k.a. ~4.5k)
1:3000–1:1501 (level 18, a.k.a. ~2k)
1:1500–1:501 (level 19, a.k.a. ~1k)
One of the most attractive areas currently online is Toronto, ON where at levels 18 and 19, individual building outlines are graced with street addresses.

Anyhow, the new tbsm45cm model will serve County of Marin’s effort at large-scale topographic mapping several ways. First, it has made possible a very detailed hillshade that helps emphasize the grading around each hillside structure in the county. Second, it helps us to create the required metric topographic contours. These are necessary to meet world-wide mapping standards, and throughout this weekend, contours are being generated from a related (smoothed version) of the terrain on 50cm vertical interval. Needless to say, most of these won’t get used in the map renderings, but the ESRI cartographers have shared a very clever indexing scheme that will help us use this single set of metric contours to support the requirements for all four of our topographic map scales.