PageSourceSearch

https://jimkeener.com/

html jimkeener.com collected 2026-10-08 16:56:39 UTC 37,884 bytes, 660 lines download raw bytes

1<!DOCTYPE html>
2<html lang="en">
3  <head>
4    <meta property="fb:pages" content="520981908088655" />
5    <meta charset="utf-8">
6    <title>My Musings</title>
7    <meta name="viewport" content="width=device-width, initial-scale=1.0">
8    <meta name="description" content="">
9    <meta name="author" content="">
10
11    <!-- Le styles -->
12    <link href="/css/bootstrap-spacelab.min.css" rel="stylesheet">
13    <style>
14      body {
15        padding-top: 60px; /* 60px to make the container go all the way to the bottom of the topbar */
16      }
17    </style>
18    <link href="/css/bootstrap-responsive.css" rel="stylesheet">
19    <link href="/css/custom.css" rel="stylesheet">
20
21    <!-- HTML5 shim, for IE6-8 support of HTML5 elements -->
22    <!--[if lt IE 9]>
23      
23<script src="http://html5shim.googlecode.com/svn/trunk/html5.js"></script>
23
24    <![endif]-->
25
26    <!-- Fav and touch icons -->
27    <link rel="apple-touch-icon-precomposed" sizes="144x144" href="/ico/apple-touch-icon-144-precomposed.png">
28    <link rel="apple-touch-icon-precomposed" sizes="114x114" href="/ico/apple-touch-icon-114-precomposed.png">
29    <link rel="apple-touch-icon-precomposed" sizes="72x72" href="/ico/apple-touch-icon-72-precomposed.png">
30    <link rel="apple-touch-icon-precomposed" href="/ico/apple-touch-icon-57-precomposed.png">
31    <link rel="shortcut icon" href="/ico/favicon.png">
32    <link rel="stylesheet" href="/css/prism.css">
33    
33<script type="text/javascript" async src="/js/prism.js"></script>
vendor: 1 bytes, line 33
33
34<script type="text/javascript" async
35  src="https://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-MML-AM_CHTML">
36</script>
36
37        
37<script src="/js/code.jquery.com/jquery-1.12.3.min.js"></script>
37
38        
38<script src="/js/cdn.datatables.net/1.10.12/js/jquery.dataTables.min.js"></script>
38
39        <link rel="stylesheet" href="/css/cdn.datatables.net/1.10.12/css/jquery.dataTables.min.css" />
40        
40<script>
41          $(function(){
42            $(".datatable").DataTable({"pageLength": 50});
43            $(".layerck").change(function(e) {
44              layers_dict[$(this).attr('id')].setVisible($(this).is(":checked"));
45            });
46          });
47        </script>
47
48
49
50
51  </head>
52
53  <body>
54  <!-- Google Tag Manager -->
vendor: 66 bytes, lines 54-56
54
55  <noscript><iframe
56  src="//www.googletagmanager.com/ns.html?id=
56GTM-KKPT75
vendor: 89 bytes, lines 56-59
56"
57  height="0" width="0"
58  style="display:none;visibility:hidden"></iframe></noscript>
59  
59<script>
vendor: 331 bytes, lines 59-64
59(function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':
60  new Date().getTime(),event:'gtm.js'});var
61  f=d.getElementsByTagName(s)[0],
62  j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src=
63  '//www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);
64  })(window,document,'script','dataLayer','
64GTM-KKPT75
vendor: 3 bytes, line 64
64');
64</script>
64
65  <!-- End Google Tag Manager -->
66
67    <div class="navbar navbar-inverse navbar-fixed-top">
68      <div class="navbar-inner">
69        <div class="container">
70          <a class="btn btn-navbar" data-toggle="collapse" data-target=".nav-collapse">
71            <span class="icon-bar"></span>
72            <span class="icon-bar"></span>
73            <span class="icon-bar"></span>
74          </a>
75          <a class="brand" href="/">My Musings</a>
76          <div class="nav-collapse collapse">
77            <ul class="nav">
78              <li class="active"><a href="/">Home</a></li>
79              <li class=""><a href="/posts">Posts</a></li>
80              <li class=""><a href="/about">About</a></li>
81              <li class=""><a href="/talks-given">Talks</a></li>
82            </ul>
83          </div><!--/.nav-collapse -->
84        </div>
85      </div>
86    </div>
87
88    <div class="container">
89
90		<div class="row">
91			<div class="span2">
92				<h3>About Me!</h1>
93				<img src="/img/me2.jpg">
94				<p class="mini-about-me">This blog is about my musings and thoughts.  I hope you find it useful, at most, and entertaining, at least.</p>
95                <p class="resume"><a href="/resume">R&eacute;sum&eacute;</a> [<a href="/pdf/James_Keener_Resume.pdf">PDF</a>]</a></p>
96				<h4>Other Pages</h4>
97				<p><a href="/quotes">Quotes</a></p>
98				<p><a href="/links">Links</a></p>
99				<h4>Presence Elsewhere</h4>
100				<p><a href="&#109;&#97;&#105;&#108;&#116;&#111;&#58;&#106;&#105;&#109;&#64;&#106;&#105;&#109;&#107;&#101;&#101;&#110;&#101;&#114;&#46;&#99;&#111;&#109;">&#106;&#105;&#109;&#64;&#106;&#105;&#109;&#107;&#101;&#101;&#110;&#101;&#114;&#46;&#99;&#111;&#109;</a></p>
101				<p><a href="https://www.github.com/jimktrains">GitHub</a></p>
102				<p><a href="https://bitbucket.org/jimktrains">BitBucket</a></p>
103			</div>
104			<div class="span8">
105				<h1><a href="/posts/luddite2">Why I&#x27;m a Luddite, kind-of (Part 2)</a></h1>
106<p>Date: <a href="/posts/by-date/2022">2022</a><del><a href="/posts/by-date/2022/05">05</a></del><a href="/posts/by-date/2022/05/27">27</a><br />
107<br />
108Tags: <a href="/posts/tags/technology">technology</a> <a href="/posts/tags/Luddite">Luddite</a> <a href="/posts/tags/freedom">freedom</a> <a href="/posts/tags/liberty">liberty</a></p>
109
110<p>In 2013 I wrote a ~~ramble~~ post called <a href="/posts/luddite">Why I&#8217;m a Luddite,
111kind-of</a>
111.  Now, some 9 years later, I want to revisit<br />
112some of my opinions then and add some more color to my basic premise. If<br />
113anything, I think my opinion has gone from &#8220;sometimes technology doesn&#8217;t<br />
114actually make things better&#8221; to &#8220;usually technology doesn&#8217;t actually<br />
115make things better.&#8221;</p><br />
116<p>I&#8217;ve come to loathe &#8220;smart&#8221; interfaces.</p><br />
117<p>In 2013 I was still two years away from having a child. While Anne and I<br />
118weren&#8217;t able to get our families to strictly adhere to the &#8220;no batteries<br />
119allowed&#8221; rule, they did follow it pretty well and I&#8217;m pretty happy with<br />
120the result. Our children watch a little more TV than I&#8217;d prefer. A<br />
121significant portion of it is made up of history and science<br />
122documentaries and they are (usually) more than willing or even beg to go<br />
123outside to play, so it&#8217;s a mixed bag.</p><br />
124<p>One sentence, of many, I do think has aged very poorly and doesn&#8217;t really reflect<br />
125my current beliefs, though, is &#8220;I just wish so many people wouldn&#8217;t use<br />
126them as a substitute for real life interaction.&#8221; Before 2020, but<br />
127specially in a post-2020 World, that sentence was just me talking out of<br />
128my butt. While I still think real-life interactions are important,<br />
129interacting with others any way in which we&#8217;re comfortable is a<br />
130positive.</p>
131
132<h1><a href="/posts/visualizing-system-wide-headway">Visualizing System-wide Headway</a></h1>
133<p>Date: <a href="/posts/by-date/2019">2019</a><del><a href="/posts/by-date/2019/11">11</a></del><a href="/posts/by-date/2019/11/13">13</a><br />
134<br />
135Tags: <a href="/posts/tags/postgis">postgis</a> <a href="/posts/tags/postgres">postgres</a> <a href="/posts/tags/gis">gis</a> <a href="/posts/tags/gtfs">gtfs</a> <a href="/posts/tags/transit">transit</a></p>
136
137	<p><p>I was bored this evening and decided to visualized the system-wide headway,<br />
138that is the amount of time between buses at a stop for every stop in the<br />
139system.</p><br />
140<p>Given the <a href="/posts/stop-spacing-part-1">stop-distance table</a> from the previous<br />
141post, I set out to compute the headway for each stop segment (the portion<br />
142of a route between stops).</p><br />
143<pre><code>&#8212; Compute the max headway for each stop per hour.<br />
144create table overall_headway as <br />
145select stop_id,                             
146       hour,                 
147       max(headway) as max_headway<br />
148from (              
149  select stop_id,                                                                     
150         extract(hour from departure_time) as hour,
151         departure_time,
152         lag(departure_time, 1) over w as prev_departure_time,
153         extract(epoch from (departure_time &#8211; lag(departure_time, 1) over w))/60 as headway
154  from gtfs.stops_3365
155  join gtfs.stop_times using (stop_id)
156  window w as (partition by stop_id order by departure_time asc)<br />
157)x                              <br />
158group by stop_id, hour;</p>
159
160
161	<p>&#8212; Assign the segments created for the stop-distance calculation<br />
162&#8212; a headway from overall_headway.<br />
163create table stop_line_headway as<br />
164select coalesce(ohs.hour, ohp.hour) as hour,
165       sdl.stop_id,
166       sdl.dist,
167       prev_stop_id,
168       min(least(ohs.max_headway, ohp.max_headway)) as headway,
169       the_geom<br />
170from stop_dist_line sdl<br />
171left join overall_headway ohs
172  on sdl.stop_id = ohs.stop_id<br />
173left join overall_headway ohp
174  on sdl.prev_stop_id = ohp.stop_id<br />
175where (ohs.hour is null or ohp.hour is null or (ohs.hour = ohp.hour))
176  and st_length(the_geom) &lt; 2500<br />
177group by coalesce(ohs.hour, ohp.hour), sdl.stop_id, sdl.dist, prev_stop_id, the_geom;</p>
178
179	<p>alter table stop_line_headway add id serial primary key;<br />
180</code></pre><br />
181<p>I then pulled this into <span class="caps">QGIS</span>, colored by graduation, and filtered per hour.<br />
182For each our, I exported a <span class="caps">PNG</span> of the map and tried by hand at <span class="caps">KDE</span>nlive for the<br />
183first time. (Although <code>ffmpeg</code> could have gotten the job done as well.)</p><br />
184<iframe width="560" height="315" src="https://www.youtube.com/embed/hZrRQW25Z_0" frameborder="0" allow="accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture" allowfullscreen></iframe></p>
185
186<p>Raw video, images, and more details on building the map to come.</p>
187
188<h1><a href="/posts/stop-spacing-part-1">Stop Spacing, Part 1</a></h1>
189<p>Date: <a href="/posts/by-date/2019">2019</a><del><a href="/posts/by-date/2019/11">11</a></del><a href="/posts/by-date/2019/11/12">12</a><br />
190<br />
191Tags: <a href="/posts/tags/postgis">postgis</a> <a href="/posts/tags/postgres">postgres</a> <a href="/posts/tags/gis">gis</a> <a href="/posts/tags/gtfs">gtfs</a> <a href="/posts/tags/transit">transit</a></p>
192
193	<p><p>For reasons I&#8217;ll get into in another post, I wanted to have a rough idea of the<br />
194distribution of the distance between stops. To do this I used the<br />
195<a href="https://gtfs.org/"><span class="caps">GTFS</span></a> import for the <a href="https://www.portauthority.org/business-center/developer-resources/">Port Authority of Allegheny County<br />
196GTFS</a> that<br />
197I did earlier this year when <a href="/posts/calculating-service-area-part-1-introduction">calculating the service<br />
198area</a>.</p><br />
199<p>A quick note about projections. Since I had created a table with a geometry in<br />
200the Pennsylvania South (ftUS) state plane( <a href="https://epsg.io/3365"><span class="caps">EPGS</span> 3365</a>),<br />
201any distance calculations I do will be done in feet.</p><br />
202<p>A quick, rough and dirty experiment went like this:</p><br />
203<pre><code>create table stop_dist_line as <br />
204select stop.stop_id,
205       stop.stop_name,
206       prev.stop_id as prev_stop_id,
207       prev.stop_name as prev_stop_name,
208       &#8212; Just using euclidean distance right now.
209       &#8212; This could be improved by using the on-the-street distance.
210       stop.the_geom &lt;-&gt; prev.the_geom as dist,
211       st_setsrid(st_makeline(stop.the_geom, prev.the_geom), 3365) as the_geom<br />
212from (
213  &#8212; This gets unique pairs of stops, regardless of ordering.
214  select distinct least(stop_id, prev_stop_id) as prev_stop_id,
215                  greatest(stop_id, prev_stop_id) as stop_id
216  &#8212; This gets all pairs of stops. In theory there could be duplicate pairs
217  &#8212; based on direction.
218  from (select stop_id,
219               lag(stop_id, 1) over stop_order prev_stop_id
220        from gtfs.stop_times
221        join gtfs.trips using (trip_id)
222        window stop_order as (partition by trip_id, direction_id order by stop_sequence)
223  ) unordered_pairs<br />
224) pairs<br />
225join gtfs.stops_3365 stop on stop.stop_id = pairs.stop_id<br />
226join gtfs.stops_3365 prev on prev.stop_id = pairs.prev_stop_id;</p>
227
228	<p>&#8212; Remove segments that don&#8217;t have a line because there is no<br />
229&#8212; previous stop to have a line or have a distance from.<br />
230delete from stop_dist_line where prev_stop_id is null;</p>
231
232	<p>&#8212; Register the new geometry column<br />
233select Populate_Geometry_Columns();<br />
234&#8212; Add a PK to play more nicely with QGIS<br />
235alter table stop_dist_line add id serial primary key;<br />
236&#8212; Add an index on the geometry.<br />
237create index on stop_dist_line using gist(the_geom);<br />
238</code></pre><br />
239<p>We can then get some idea of the distribution by bucketing by 100s.</p><br />
240<pre><code>select dist_group,
241       ttlsf,
242       ttl,
243       round(ttlsf/ttl, 2) as perc,
244       c
245from
246  (select dist_group,
247          sum(c) over (rows between unbounded preceding and current row) as ttlsf,
248          sum(c) over () as ttl,
249          c
250   from (
251      select 100*floor(dist / 100) as dist_group,
252             count(*) as c
253      from stop_dist_line
254      group by dist_group
255      order by dist_group) x)y;
256</code></pre><br />
257<p>Giving</p><br />
258<pre><code> dist_group | ttlsf | ttl  | perc |  c
259------------+-------+------+------+------
260          0 |   149 | 7732 | 0.02 |  149
261        100 |   223 | 7732 | 0.03 |   74
262        200 |   681 | 7732 | 0.09 |  458
263        300 |  1568 | 7732 | 0.20 |  887
264        400 |  2528 | 7732 | 0.33 |  960
265        500 |  3579 | 7732 | 0.46 | 1051
266        600 |  4462 | 7732 | 0.58 |  883
267        700 |  5135 | 7732 | 0.66 |  673
268        800 |  5589 | 7732 | 0.72 |  454
269        900 |  5946 | 7732 | 0.77 |  357
270       1000 |  6234 | 7732 | 0.81 |  288
271       1100 |  6456 | 7732 | 0.83 |  222
272       1200 |  6628 | 7732 | 0.86 |  172
273       1300 |  6766 | 7732 | 0.88 |  138
274       1400 |  6852 | 7732 | 0.89 |   86
275       1500 |  6921 | 7732 | 0.90 |   69
276       ...
277</code></pre><br />
278<p>However, because PostgreSQL is awesome, we can also just compute percentiles<br />
279directly.</p><br />
280<pre><code>select percentile_cont(0.25) within group (order by dist) as p_25,
281       percentile_cont(0.35) within group (order by dist) as p_35,
282       percentile_cont(0.45) within group (order by dist) as p_45,
283       percentile_cont(0.50) within group (order by dist) as p_50
284from stop_dist_line;
285</code></pre><br />
286<p>Giving</p><br />
287<pre><code>        p_25       |       p_35       |       p_45       |       p_50
288------------------+------------------+------------------+------------------
289 438.775035648918 | 518.085656645151 | 590.352075165022 | 631.007847560984
290</code></pre><br />
291<p>Next steps would be to correct many of the little corners cut, e.g. not using<br />
292on-the-street distance, and rerunning the calculations.</p></p>
293
294<h1><a href="/posts/calculating-service-area-part-4-next-steps">Calculating Service Area, Part 4: Next Steps</a></h1>
295<p>Date: <a href="/posts/by-date/2019">2019</a><del><a href="/posts/by-date/2019/01">01</a></del><a href="/posts/by-date/2019/01/27">27</a><br />
296<br />
297Tags: <a href="/posts/tags/postgis">postgis</a> <a href="/posts/tags/postgres">postgres</a> <a href="/posts/tags/transit">transit</a> <a href="/posts/tags/census">census</a></p>
298
299<p>In <a href="/posts/calculating-service-area-part-1-introduction">Part 1</a> we<br />
300talked about the project and got some of the tools and data we need to<br />
301start. A quick reminder of our motivation: to figure out the percentage<br />
302area and population the county transit agency serves by using a<br />
303&#8220;walkshed&#8221;. </p><br />
304<h2>Parts</h2><br />
305<ul><br />
306<li><a href="/posts/calculating-service-area-part-1-introduction">Part 1: Introduction and Setup</a></li><br />
307<li><a href="/posts/calculating-service-area-part-2-defining-the-walkshed">Part 2: Defining the Walkshed</a></li><br />
308<li><a href="/posts/calculating-service-area-part-3-demographics">Part 3: Demographics</a></li><br />
309<li>Part 4: Next Steps</li><br />
310</ul><br />
311<h2>Walkshed</h2><br />
312<p>There are two changes I should make to the walkshed calculation before I dig<br />
313into per-route statistics: (1) excluding interstates and (2) light rail and<br />
314busways stations should have a half mile walkshed.  I don&#8217;t believe they made a<br />
315difference at the county level; once I redo the calculations I&#8217;ll come back to<br />
316see if I have to eat those words.</p><br />
317<h2>More In-Depth Demographic Analysis</h2><br />
318<p>Since we have all of the census data, all of the route data, and we broke our<br />
319walkshed down by stop it&#8217;d be quite easy to look at the service area as a<br />
320whole, and routes in particular against other socio-economic indicators.</p><br />
321<p>Also, more interestingly, who <em>isn&#8217;t</em> being served?</p><br />
322<h2>Filter Walkshed by Reliable Transit</h2><br />
323<p>I&#8217;d love to look at what the walkshed is when only stops that are served more<br />
324frequently than every 30 min mid-day is included. Another alternative would be<br />
325late at night, or weekends. Essentially, what&#8217;s the service area that you could<br />
326do a good deal in without having to get into a car?</p><br />
327<h2>Build a web-based tool to export this data</h2><br />
328<p>I would love to make this data more usable and accessible.  I previously made a<br />
329<a href="https://github.com/jimktrains/pptmap">little tool</a> using python, mapnik, and a<br />
330file from GDAL2Tiles. I&#8217;d like to explore if <span class="caps">QGIS</span>&#8217; &#8220;save project to database&#8220;<br />
331functionality could be used (<em>i.e.</em> mangled, forced, hacked?) to do this more<br />
332cleanly.</p><br />
333<h2>Isochron / Isotrip maps</h2><br />
334<p>I&#8217;ve been wanting to build an application where you could pick a location<br />
335and it&#8217;d generate two (interactive) maps of all other stops: an isochron map<br />
336showing how long it&#8217;d take to get anywhere else (and that has 15-min isochron<br />
337lines showing the boundry of how long it takes to get how far) and an isotrip<br />
338map showing how many trips/transfers it&#8217;d take to get to the other stop, with<br />
339similar isotrip lines.</p><br />
340<h2>Create Subject Tables</h2><br />
341<p>The Access DBs also contain information on what all the fields mean and which<br />
342subject tables they belong to.  We can extract it with the <code>mdb-export</code> tool.</p>
343<pre><code>mdb-export SF1_Access2003.mdb data_field_descriptors &amp;gt; data_field_descriptors.csv
344</code></pre>
345
346<h1><a href="/posts/calculating-service-area-part-3-demographics">Calculating Service Area, Part 3: Demographics</a></h1>
347<p>Date: <a href="/posts/by-date/2019">2019</a><del><a href="/posts/by-date/2019/01">01</a></del><a href="/posts/by-date/2019/01/27">27</a><br />
348<br />
349Tags: <a href="/posts/tags/postgis">postgis</a> <a href="/posts/tags/postgres">postgres</a> <a href="/posts/tags/transit">transit</a> <a href="/posts/tags/census">census</a></p>
350
351	<p><p>In <a href="/posts/calculating-service-area-part-1-introduction">Part 1</a> we<br />
352talked about the project and got some of the tools and data we need to<br />
353start. A quick reminder of our motivation: to figure out the percentage<br />
354area and population the county transit agency serves by using a<br />
355&#8220;walkshed&#8221;. </p><br />
356<p>In <a href="/posts/calculating-service-area-part-2-defining-the-walkshed">Part 2</a>, we<br />
357defined what our walkshed is and defined how we would be linking it to the<br />
358Census 2010. In this post, we&#8217;ll be doing the linking</p><br />
359<p>I want to point out <a href="http://sproke.blogspot.com/2012/01/importing-2010-sf1-census-in-postgresql.html">sproke&#8217;s blogpost from<br />
3602012</a><br />
361where I got many of these ideas many years back.</p><br />
362<h2>Parts</h2><br />
363<ul><br />
364<li><a href="/posts/calculating-service-area-part-1-introduction">Part 1: Introduction and Setup</a></li><br />
365<li><a href="/posts/calculating-service-area-part-2-defining-the-walkshed">Part 2: Defining the Walkshed</a></li><br />
366<li>Part 3: Demographics</li><br />
367<li><a href="/posts/calculating-service-area-part-4-next-steps">Part 4: Next Steps</a></li><br />
368</ul><br />
369<h2>Importing the Census 2010 Data</h2><br />
370<p>First we need to gather the data from <a href="https://www.census.gov/">census.gov</a> and<br />
371initialized our database. I&#8217;m going to grab the documentation (sf1.pdf), along with<br />
372some MS Access database files they use to store the schema.</p><br />
373<pre><code>wget --continue -e robots=off --force-directories --no-parent https://www.census.gov/prod/cen2010/doc/sf1.pdf
374wget --recursive --continue -e robots=off --force-directories --no-parent https://www2.census.gov/census_2010/04-Summary_File_1/SF1_Access2003.mdb
375wget --recursive --continue -e robots=off --force-directories --no-parent https://www2.census.gov/census_2010/04-Summary_File_1/SF1_Access2007.accdb
376wget --recursive --continue -e robots=off --force-directories --no-parent https://www2.census.gov/census_2010/04-Summary_File_1/Pennsylvania/
377wget --recursive --continue -e robots=off --force-directories --no-parent https://www2.census.gov/census_2010/05-Summary_File_2/SF2_MSAccess_2003.mdb
378wget --recursive --continue -e robots=off --force-directories --no-parent https://www2.census.gov/census_2010/05-Summary_File_2/SF2_MSAccess_2007.accdb
379wget --recursive --continue -e robots=off --force-directories --no-parent https://www2.census.gov/census_2010/05-Summary_File_2/Pennsylvania/
380</code></pre><br />
381<p>Some notes, <code>--recursive</code> will grab linked files, <code>--continue</code> will try to<br />
382start mid-file if it&#8217;s been partially downloaded. <code>-e robots=off</code> is naughty<br />
383and tels wget to download things even if robots.txt tells you not to. (I don&#8217;t<br />
384think it&#8217;s an issue here, I just add it all the time.) <code>--force-directories</code><br />
385will create a directory structure begining with the domain so<br />
386&#8220;https://example.com/a/b/c.d&#8221; will be downloaded to a file<br />
387&#8220;example.com/a/b/c.d&#8221; &#8212; I find this makes keeping track of what things are and<br />
388where they came from a lot easier. <code>--no-parent</code> tells wget not to go up the<br />
389directory tree to get assets, which could cause it to begin going down sibling<br />
390directories of the root we gave it.</p><br />
391<p>Now, we are going import that data. The first step <br />
392is to unzip the datafile. Be warned, the unziped data can swell in size. For<br />
393example, the pa dataset is compressed 90%: the zip is 436MiB and uncompressed<br />
394is 4.4GiB.</p><br />
395<pre><code>unzip pa2010.sf1.zip
396</code></pre><br />
397<p>Now that we have the files, let&#8217;s build the schema. Mostly because I love<br />
398using the commandline whenever possible, I decided to stream the schema<br />
399from the Access sample databases through some fix-ups and then into a file.</p><br />
400<pre><code>mdb-schema SF1_Access2003.mdb |
401  # MS Access uses [] as idenifier quotes, let&#8217;s fix that
402  tr &#8216;[]&#8217; &#8216;&#8221;&#8220;&#8217; |</p>
403
404  # Some tables have a name suffix of `mod` because the first few have been
405  # removed to keep the table under 255 fields; let&#8217;s add them back in.
406  sed &#8216;/mod&#8221;$/{$!{N;s/mod&#8221;\n (/&#8221;\n (\n        &#8220;<span class="caps">FILEID</span>&#8221;                        text,\n        &#8220;<span class="caps">STUSAB</span>&#8221;                        text,\n        &#8220;<span class="caps">CHARITER</span>&#8221;                      text,\n        &#8220;<span class="caps">CIFSN</span>&#8221;                  text,/g}}&#8217; |
407
408  # Remove spaces from identifiers.
409  sed -E &#8216;s/(\s+)&#8221;([<sup>&#8221;]+) ([</sup>&#8221;]+)&#8221;/\1\2_\3/g&#8217; |
410
411  # Convert to PG type names
412  sed -E &#8216;s/Text\s*\([0-9]+\)/text/g&#8217; |
413  sed -E &#8216;s/(Long )?Integer/bigint/g&#8217; |
414  sed -E &#8216;s/Double/double precision/g&#8217; |
415
416  # Table 45 is split into 2 parts; let&#8217;s merge them by deleting the trailing
417  # ); of the one and the <span class="caps">CREATE</span> <span class="caps">TABLE</span> along with the first few fields of the
418  # second.
419  sed &#8216;
419s/_PT1//g&#8217; |
420  sed &#8216;/);$/{$!{N;N;N;N;N;N;N;N;s/.*_PT2.*LOG[^,]*//g}}&#8217; |
421
422  # Just remove the quotes &#8212; they&#8217;re not needed anymore
423  tr -d &#8216;&#8221;&#8217; &gt; SF1_pg.schma.sql
424</code></pre>
425<p>Then, we can send that into postgres, but I want to wrap it all in a commit<br />
426so that we don&#8217;t end up with a half-built database if something fails.</p><br />
427<pre><code>cat &amp;lt;(echo "set role census;") \
428    &amp;lt;(echo "begin; create schema census2010; set search_path = census2010;") \
429    SF1_pg.schma.sql \
430    &amp;lt;(echo "commit;") |
431  psql -v ON_ERROR_STOP=1 census
432</code></pre><br />
433<p>Now the fun part! We get to import the datafiles! I&#8217;m going to use the<br />
434<a href="https://www.postgresql.org/docs/10/sql-copy.html">copy</a> statement and the <br />
435psql <a href="https://www.postgresql.org/docs/10/app-psql.html#APP-PSQL-META-COMMANDS-COPY">\copy</a> meta-command. However, like above,<br />
436I&#8217;m going to build up a sql file before running it so that everything can<br />
437be wrapped in a transaction. Yes, this step can take a bit and if it fails<br />
438towards the end it sucks, but in the end it doesn&#8217;t take that long (a few<br />
439minutes) and it&#8217;s significantly easier than trying to scrub duplicates. (We<br />
440haven&#8217;t built indices yet as having them during the input slows it down<br />
441considerably.</p><br />
442<pre><code>echo "begin;" &amp;gt; import_script
443for input_file in pa000*2010.sf1; do
444  table_name=SF1_`echo $input_file | sed 's/pa\(.\{5\}\)2010.sf1/\1/'`
445  # Provides some feedback as `copy` will only output the number of rows imported.
446  echo "\\\echo $input_file  $table_name" &amp;gt;&amp;gt; import_script
447  echo "\\\copy census2010.${table_name} from '${input_file}' with (format csv, header false);\n" &amp;gt;&amp;gt; import_script
448done
449echo "commit;" &amp;gt;&amp;gt; import_script
450cat import_script | psql -v ON_ERROR_STOP=1 census
451</code></pre><br />
452<p>Now comes the more interesting part: the geoheader (links demographic records<br />
453with the geography) is a fixed-width file (<em>i.e.</em> fields are not delimited, but<br />
454are defined as between character positions on a line). What Sophia did was<br />
455to build a staging table that would contain the raw lines, and then do a<br />
456&#8220;insert into select&#8221; with a select statement that teased out the fields using<br />
457the substring function.  We&#8217;ll do a very similar thing. The fields are<br />
458defined in the <a href="https://www.census.gov/prod/cen2010/doc/sf1.pdf">SF1</a><br />
459documentation</p><br />
460<pre><code>echo "set role census; CREATE TABLE census2010.geo_header_staging(data text);" | psql census
461echo "\\\copy census2010.geo_header_staging from 'pageo2010.sf1';" | psql census
462</code></pre><br />
463<p>Then we&#8217;ll tease apart that staging table.</p><br />
464<pre><code>INSERT INTO census2010.geo_header_sf1 (fileid, stusab, sumlev, geocomp, chariter, cifsn, logrecno, region, division, state, county, countycc, countysc, cousub, cousubcc, cousubsc, place, placecc, placesc, tract, blkgrp, block, iuc, concit, concitcc, concitsc, aianhh, aianhhfp, aianhhcc, aihhtli, aitsce, aits, aitscc, ttract, tblkgrp, anrc, anrccc, cbsa, cbsasc, metdiv, csa, necta, nectasc, nectadiv, cnecta, cbsapci, nectapci, ua, uasc, uatype, ur, cd, sldu, sldl, vtd, vtdi, reserve2, zcta5, submcd, submcdcc, sdelem, sdsec, sduni, arealand, areawatr, name, funcstat, gcuni, pop100, hu100, intptlat, intptlon, lsadc, partflag, reserve3, uga, statens, countyns, cousubns, placens, concitns, aianhhns, aitsns, anrcns, submcdns, cd113, cd114, cd115, sldu2, sldu3, sldu4, sldl2, sldl3, sldl4, aianhhsc, csasc, cnectasc, memi, nmemi, puma, reserved)
465SELECT
466nullif(trim(substring(data,1,6)), '') AS fileid,
467nullif(trim(substring(data,7,2)), '') AS stusab,
468nullif(trim(substring(data,9,3)), '') AS sumlev,
469nullif(trim(substring(data,12,2)), '') AS geocomp,
470nullif(trim(substring(data,14,3)), '') AS chariter,
471nullif(trim(substring(data,17,2)), '') AS cifsn,
472nullif(trim(substring(data,19,7)), '')::bigint AS logrecno,
473nullif(trim(substring(data,26,1)), '') AS region,
474nullif(trim(substring(data,27,1)), '') AS division,
475nullif(trim(substring(data,28,2)), '') AS state,
476nullif(trim(substring(data,30,3)), '') AS county,
477nullif(trim(substring(data,33,2)), '') AS countycc,
478nullif(trim(substring(data,35,2)), '') AS countysc,
479nullif(trim(substring(data,37,5)), '') AS cousub,
480nullif(trim(substring(data,42,2)), '') AS cousubcc,
481nullif(trim(substring(data,44,2)), '') AS cousubsc,
482nullif(trim(substring(data,46,5)), '') AS place,
483nullif(trim(substring(data,51,2)), '') AS placecc,
484nullif(trim(substring(data,53,2)), '') AS placesc,
485nullif(trim(substring(data,55,6)), '') AS tract,
486nullif(trim(substring(data,61,1)), '') AS blkgrp,
487nullif(trim(substring(data,62,4)), '') AS block,
488nullif(trim(substring(data,66,2)), '') AS iuc,
489nullif(trim(substring(data,68,5)), '') AS concit,
490nullif(trim(substring(data,73,2)), '') AS concitcc,
491nullif(trim(substring(data,75,2)), '') AS concitsc,
492nullif(trim(substring(data,77,4)), '') AS aianhh,
493nullif(trim(substring(data,81,5)), '') AS aianhhfp,
494nullif(trim(substring(data,86,2)), '') AS aianhhcc,
495nullif(trim(substring(data,88,1)), '') AS aihhtli,
496nullif(trim(substring(data,89,3)), '') AS aitsce,
497nullif(trim(substring(data,92,5)), '') AS aits,
498nullif(trim(substring(data,97,2)), '') AS aitscc,
499nullif(trim(substring(data,99,6)), '') AS ttract,
500nullif(trim(substring(data,105,1)), '') AS tblkgrp,
501nullif(trim(substring(data,106,5)), '') AS anrc,
502nullif(trim(substring(data,111,2)), '') AS anrccc,
503nullif(trim(substring(data,113,5)), '') AS cbsa,
504nullif(trim(substring(data,118,2)), '') AS cbsasc,
505nullif(trim(substring(data,120,5)), '') AS metdiv,
506nullif(trim(substring(data,125,3)), '') AS csa,
507nullif(trim(substring(data,128,5)), '') AS necta,
508nullif(trim(substring(data,133,2)), '') AS nectasc,
509nullif(trim(substring(data,135,5)), '') AS nectadiv,
510nullif(trim(substring(data,140,3)), '') AS cnecta,
511nullif(trim(substring(data,143,1)), '') AS cbsapci,
512nullif(trim(substring(data,144,1)), '') AS nectapci,
513nullif(trim(substring(data,145,5)), '') AS ua,
514nullif(trim(substring(data,150,2)), '') AS uasc,
515nullif(trim(substring(data,152,1)), '') AS uatype,
516nullif(trim(substring(data,153,1)), '') AS ur,
517nullif(trim(substring(data,154,2)), '') AS cd,
518nullif(trim(substring(data,156,3)), '') AS sldu,
519nullif(trim(substring(data,159,3)), '') AS sldl,
520nullif(trim(substring(data,162,6)), '') AS vtd,
521nullif(trim(substring(data,168,1)), '') AS vtdi,
522nullif(trim(substring(data,169,3)), '') AS reserve2,
523nullif(trim(substring(data,172,5)), '') AS zcta5,
524nullif(trim(substring(data,177,5)), '') AS submcd,
525nullif(trim(substring(data,182,2)), '') AS submcdcc,
526nullif(trim(substring(data,184,5)), '') AS sdelem,
527nullif(trim(substring(data,189,5)), '') AS sdsec,
528nullif(trim(substring(data,194,5)), '') AS sduni,
529nullif(trim(substring(data,199,14)), '')::double precision AS arealand,
530nullif(trim(substring(data,213,14)), '')::double precision AS areawatr,
531nullif(trim(substring(data,227,90)), '') AS name,
532nullif(trim(substring(data,317,1)), '') AS funcstat,
533nullif(trim(substring(data,318,1)), '') AS gcuni,
534nullif(trim(substring(data,319,9)), '')::bigint AS pop100,
535nullif(trim(substring(data,328,9)), '')::bigint AS hu100,
536nullif(trim(substring(data,337,11)), '') AS intptlat,
537nullif(trim(substring(data,348,12)), '') AS intptlon,
538nullif(trim(substring(data,360,2)), '') AS lsadc,
539nullif(trim(substring(data,362,1)), '') AS partflag,
540nullif(trim(substring(data,363,6)), '') AS reserve3,
541nullif(trim(substring(data,369,5)), '') AS uga,
542nullif(trim(substring(data,374,8)), '') AS statens,
543nullif(trim(substring(data,382,8)), '') AS countyns,
544nullif(trim(substring(data,390,8)), '') AS cousubns,
545nullif(trim(substring(data,398,8)), '') AS placens,
546nullif(trim(substring(data,406,8)), '') AS concitns,
547nullif(trim(substring(data,414,8)), '') AS aianhhns,
548nullif(trim(substring(data,422,8)), '') AS aitsns,
549nullif(trim(substring(data,430,8)), '') AS anrcns,
550nullif(trim(substring(data,438,8)), '') AS submcdns,
551nullif(trim(substring(data,446,2)), '') AS cd113,
552nullif(trim(substring(data,448,2)), '') AS cd114,
553nullif(trim(substring(data,450,2)), '') AS cd115,
554nullif(trim(substring(data,452,3)), '') AS sldu2,
555nullif(trim(substring(data,455,3)), '') AS sldu3,
556nullif(trim(substring(data,458,3)), '') AS sldu4,
557nullif(trim(substring(data,461,3)), '') AS sldl2,
558nullif(trim(substring(data,464,3)), '') AS sldl3,
559nullif(trim(substring(data,467,3)), '') AS sldl4,
560nullif(trim(substring(data,470,2)), '') AS aianhhsc,
561nullif(trim(substring(data,472,2)), '') AS csasc,
562nullif(trim(substring(data,474,2)), '') AS cnectasc,
563nullif(trim(substring(data,476,1)), '') AS memi,
564nullif(trim(substring(data,477,1)), '') AS nmemi,
565nullif(trim(substring(data,478,5)), '') AS puma,
566nullif(trim(substring(data,483,18)), '') AS reserved
567FROM census2010.geo_header_staging;
568</code></pre><br />
569<p>Now, let&#8217;s delete that staging table.</p><br />
570<pre><code>drop table census2010.geo_header_staging;
571</code></pre><br />
572<p>Let&#8217;s build some indices and foreign keys.</p>
573<pre><code>echo &#8220;begin;&#8221; &gt; create_index
574for table in `echo &#8220;\\\dt census2010.*&#8221; | psql -tA census | grep -v data | cut -f 2 -d \|`; do
575  cat &lt;&lt;<span class="caps">EOS</span> &gt;&gt; create_index
576create index on census2010.${table}(fileid);
577create index on census2010.${table}(stusab);
578create index on census2010.${table}(chariter);
579create index on census2010.${table}(cifsn);
580create unique index on census2010.${table}(logrecno);
581
582	<p>alter table census2010.${table} alter column fileid set not null;<br />
583alter table census2010.${table} alter column stusab set not null;<br />
584alter table census2010.${table} alter column chariter set not null;<br />
585alter table census2010.${table} alter column cifsn set not null;<br />
586alter table census2010.${table} alter column logrecno set not null;</p>
587
588	<p>alter table census2010.${table} add primary key using index ${table}_logrecno_idx;</p>
589
590	<p>EOS<br />
591done<br />
592echo &#8220;commit;&#8221; &gt;&gt; create_index<br />
593cat create_index | psql -v ON_ERROR_STOP=1 census</p>
594
595	<p>create unique index on census2010.geo_header_sf1(fileid, stusab, chariter, cifsn, logrecno);</p>
596
597	<p>echo &#8220;begin;&#8221; &gt; create_fk<br />
598for table in `echo &#8220;\\\dt census2010.*&#8221; | psql -tA census | grep -v data | grep -v geo | cut -f 2 -d \|`; do
599  cat &lt;&lt;<span class="caps">EOS</span> &gt;&gt; create_fk<br />
600alter table census2010.${table} add constraint ${table}_geo_fk foreign key (logrecno) references census2010.geo_header_sf1(logrecno);<br />
601EOS<br />
602done<br />
603echo &#8220;commit;&#8221; &gt;&gt; create_fk<br />
604cat create_fk | psql -v ON_ERROR_STOP=1 census<br />
605</code></pre><br />
606<h2>Service Area Population</h2><br />
607<p>In order to speed up lookups between the <span class="caps">TIGER</span> data, let&#8217;s  create an index to<br />
608use. <code>sumlev</code> (Summary Level) is the aggregation that the row represents,<br />
609<em>e.g.</em> state, county, track, block group, block, &amp;c. They&#8217;re defined in the<br />
610<a href="https://www.census.gov/prod/cen2010/doc/sf1.pdf">SF1</a> documentation. <code>101</code> is<br />
611the block level.</p><br />
612<pre><code>create index on census2010.geo_header_sf1 (state, county, tract, block) where sumlev = '101';
613</code></pre><br />
614<p>Calculating the total population of the county can be done by simply using<br />
615field <code>p0010001</code> (total population) in table <code>sf1_00001</code> (also defined in the<br />
616SF1 documentation), but only for the rows representing the block summary level.<br />
617Without filtering on the summary level the population would be counted for<br />
618every summary level in the county.</p><br />
619<pre><code>select sum(p0010001) 
620from census2010.sf1_00001 
621join census2010.geo_header_sf1 using (logrecno)
622where county = '003' 
623  and sumlev = '101';
624</code></pre><br />
625<p>For the population in our walkshed, we will use the <code>stops_walkshed_blocks</code><br />
626table that we generated in Part 2.  As a recap, that is every block that is<br />
627adjacent to a segment of road that is within a quarter mile of a bus stop. (The<br />
628reason for the <code>group by</code> is that <code>stops_walkshed_blocks</code> is broken down by<br />
629stop, which obviously leads to some blocks being counted multiple times.)</p><br />
630<pre><code>select sum(p0010001)
631from (
632  select max(p0010001) as p0010001 
633  from census2010.sf1_00001
634  join census2010.geo_header_sf1 using (logrecno) 
635  join tiger2017.tabblock 
636    on  statefp10 = state 
637    and countyfp10 = county
638    and tractce10 = tract 
639    and blockce10 = block
640  join stops_walkshed_blocks 
641    on stops_walkshed_blocks.geoid10 = tabblock.geoid10
642  where county = '003' and sumlev = '101'
643  group by tabblock.geoid10) x;
644</code></pre><br />
645<h2>Results</h2><br />
646<p>This yields a total population of 1,223,348 and a service area population of<br />
647775,835 which is just over 63% of the county&#8217;s population.</p></p>
648
649			</div>
650		</div>
651
652	<hr>
653	<footer class="footer" id="footer">
654	<div
655	<p>This site is made possible by <a href="https://github.com/jimktrains/gus">gus</a> and <a href="http://twitter.github.com/bootstrap/">Bootstrap</a>. <a href="http://jimkeener.com/rss.xml">RSS</a></p>
656<p>My words are copyright Jim Keener 2023. The technology is all <a href="http://www.gnu.org/philosophy/free-software-for-freedom.html">F/OSS</a>.</p>
657    </footer>
658    </div> <!-- /container -->
659  </body>
660</html>

Line numbers count LF bytes from the start of the resource, as the search results do. Vendor segments are library code the classifier recognised; they are stored but not indexed. Bytes are shown as Latin1 characters, one per byte.