Estimated Subsidence in the San Joaquin Valley between 1949 – 2005

San Joaquin Valley Subsidence Analysis README.  
Written: Joel Dudas, 3/12/2017. Amended: Ben Brezing, 4/2/2019. DWR’s Division of Engineering Geodetic Branch received a request in 1/2017 from Jeanine Jones to produce a graphic of historic subsidence in the entirety of the San Joaquin Valley. The task was assigned to the Mapping & Photogrammetry Office and the Geospatial Data Support Section to complete by early February. After reviewing the alternatives, the decision was made to produce contours from the oldest available set of quad maps for which there was reasonable certainty about quality and datum, and to compare that to the most current Valley-wide DEM. For the first requirement, research indicated that the 1950’s vintage quad maps for the Valley were the best alternative. Prior quad map editions are uneven in quality and vintage, and the actual control used for the contour lines was extremely suspect. The 1950’s quads, by contrast, were produced primarily on the basis of 1948-1949 aerial photography, along with control corresponding to that period, and referenced to the National Geodetic Vertical Datum of 1929. For the current set, the most recent Valley-wide dataset that was freely available, in the public domain, and of reasonable accuracy was the 2005 NextMap SAR acquisition (referenced to NAVD88). The primary bulk of the work focused on digitizing the 1950’s contours. First, all of the necessary quads were downloaded from the online USGS quad source https://ngmdb.usgs.gov/maps/Topoview/viewer/#4/41.13/-107.51. Then the entire staff of the Mapping & Photogrammetry Lab (including both the Mapping Office and GDDS staff) proceeded to digitize the contours. Given the short turnaround time constraint and limited budget, certain shortcuts occurred in contour development. While efforts were made to digitize accurately, speed really was important. Contours were primarily focused only on agricultural and other lowland areas, and so highlands were by and large skipped. The tight details of contours along rivers, levees, and hillsides was skipped and/or simplified. In some cases, only major contours were digitized. The mapping on the source quads itself varied….in a few cases on spot elevations on benchmarks were available in quads. The contour interval sometimes varied, even within the quad sheet itself. In addition, because 8 different people were creating the contours, variability exists in the style and attention to detail. It should be understood that given the purpose of the project (display regional subsidence patterns), that literal and precise development of the historic contour sets leaves some things to be desired. These caveats being said, the linework is reasonably accurate for what it is (particularly given that the contours of that era themselves were mapped at an unknown and varying actual quality). The digitizers tagged the lines with Z values manually entered after linework that corresponded to the mapped elevation contours. Joel Dudas then did what could be called a “rough” QA/QC of the contours. The individual lines were stitched together into a single contour set, and exported to an elevation raster (using TopoToRaster in ArcGIS 10.4). Gross blunders in Z values were corrected. Gaps in the coverage were filled. The elevation grid was then adjusted to NAVD88 using a single adjustment for the entire coverage area (2.5’, which is a pretty close average of values in this region). The NextMap data was extracted for the area, and converted into feet. The two raster sets were fixed to the same origin point. The subsidence grid was then created by subtracting the old contour-derived grid from the NextMAP DEM. The subsidence grid that includes all of the values has the suffix “ALL”. Then, to improve the display fidelity, some of the extreme values (above +5’ and below -20’*) were filtered out of the dataset, and the subsidence grid was regenerated for these areas and suffixed with “cut.” The purpose of this cut was to extract some of the riverine and hilly areas that produced more extreme values and other artifacts purely due to the analysis approach (i.e. not actual real elevation change). * - some of the areas with more than 20 feet of subsidence were omitted from this clipping, because they were in heavily subsided areas and may be “real subsidence.”The resulting subsidence product should be perceived in light of the above. Some of the collar of the San Joaquin Valley shows large changes, but that is simply due to the analysis method. Also, individual grid cells may or may not be comparing the same real features. Errors are baked into both comparison datasets. However, it is important to note that the large areas of subsidence in the primary agriculture area agree fairly well with a cruder USGS subsidence map of the Valley based on extensometer data. We have confidence that the big picture story these results show us is largely correct, and that the magnitudes of subsidence are somewhat reasonable. The contour set can serve as the baseline to support future comparisons using more recent or future data as it becomes available. It should be noted there are two key versions of the data. The “Final Deliverables” from 2/2017 were delivered to support the initial Public Affairs press release. Subsequent improvements were made in coverage and blunder correction as time permitted (it should be noted this occurred in the midst of the Oroville Dam emergency) to produce the final as of 3/12/2017. Further improvements in overall quality and filtering could occur in the future if time and needs demand it.  
Update (4/3/2019, Ben Brezing): The raster was further smoothed to remove artifacts that result from comparing the high resolution NextMAP DEM to the lower resolution DEM that was derived from the 1950’s quad map contours. The smoothing was accomplished by removing raster cells with values that are more than 0.5 feet different than adjacent cells (25 meter cell size), as well as the adjacent cells. The resulting raster was then resampled to a raster with 100 meter cell size using cubic resampling technique and was then converted to a point feature class. The point feature class was then interpolated to a raster with 250 meter cell size using the IDW technique, a fixed search radius of 1250 meters and power=2. The resulting raster was clipped to a smaller extent to remove noisier areas around the edges of the Central Valley while retaining coverage for the main area of interest.