Extracting water quality data from DataRetrieval R Package
In this blog post, I will be showing some capabilities of dataRetrieval
package by USGS to extract water quality data.
To install the dataRetrieval
package use the following command:
devtools::install_github(repo = "USGS-R/dataRetrieval")
The main function I will be discussing is readWQPdata
. The users can search the data using various options as shown below:
- bBox = Bounding box that uses the coordinates of lower left corner and upper right corner
- lat / long = lat / long will be specified by the user if they are interested to see if any data is available within radial distance.
- countycode
- startDate:
- endDate:
- characteristicType
- characteristicName : CharacteristicName can be one of following
read.table(file = "../../data/0219_waterquality/characteristic_name.txt",sep = "\n", header = F ) %>%
knitr::kable()
V1 |
---|
Conductivity |
Turbidity |
Alkalinity, total |
Biochemical oxygen demand, standard conditions |
Chloride |
Chemical oxygen demand |
Hardness, Ca, Mg |
Inorganic nitrogen (nitrate and nitrite) |
Total dissolved solids |
Kjeldahl nitrogen |
Phosphorus |
Total suspended solids |
Temperature, air |
Temperature, water |
pH |
Dissolved oxygen (DO) |
Weather condition (WMO code 4501) (choice list) |
Cadmium |
Cyanide |
Chromium |
Copper |
Fecal Coliform |
Iron |
Manganese |
Lead |
Phenol |
Organic Nitrogen |
Aluminum |
Arsenic |
Mercury |
Ammonia-nitrogen |
Silver |
Nickel |
Zinc |
Depth, bottom |
RBP Water Odors (choice list) |
RBP Water Surface Oils (choice list) |
Chlorophyll a |
Orthophosphate |
Specific conductance |
Antimony |
Selenium |
Thallium |
Calcium |
Magnesium |
Escherichia coli |
Dissolved oxygen saturation |
Depth, data-logger (non-ported) |
Light attenuation, depth at 99% |
Depth, Secchi disk depth |
Flow, stream stage (choice list) |
Extracting walker county Temperature data
Let us try to see how many stations are there is walker county Alabama that has data for 2011 period.
To find out the countycode we can use the function:
dataRetrieval::countyCdLookup(state = "AL", county = "walker")
## [1] "127"
We can see that the county code of Walker County
is 127
. Now, let us use readWQPdata
function with countycode as follows:
tempwalkercounty <- readWQPdata(countycode="US:01:127",startDate="2011-01-01",
endDate = "2011-12-31",
characteristicName="Temperature, water")
Use of bounding box to extract the data
temp_bbox <- readWQPdata(bBox=c(-87.26, 33.59, -87.07, 33.825),
startDate="2011-01-01",
endDate = "2011-12-31",
characteristicName="Temperature, water")
head(temp_bbox)
## # A tibble: 6 x 65
## OrganizationIden~ OrganizationFormal~ ActivityIdentif~ ActivityTypeCode
## <chr> <chr> <chr> <chr>
## 1 21AWIC ALABAMA DEPT. OF E~ 21AWIC-154896_8~ Field Msr/Obs-Po~
## 2 21AWIC ALABAMA DEPT. OF E~ 21AWIC-154889_8~ Field Msr/Obs-Po~
## 3 21AWIC ALABAMA DEPT. OF E~ 21AWIC-154889_8~ Field Msr/Obs-Po~
## 4 21AWIC ALABAMA DEPT. OF E~ 21AWIC-155417_9~ Field Msr/Obs-Po~
## 5 21AWIC ALABAMA DEPT. OF E~ 21AWIC-154896_8~ Field Msr/Obs-Po~
## 6 21AWIC ALABAMA DEPT. OF E~ 21AWIC-154911_8~ Field Msr/Obs-Po~
## # ... with 61 more variables: ActivityMediaName <chr>,
## # ActivityMediaSubdivisionName <chr>, ActivityStartDate <date>,
## # ActivityStartTime.Time <chr>, ActivityStartTime.TimeZoneCode <chr>,
## # ActivityEndDate <date>, ActivityEndTime.Time <chr>,
## # ActivityEndTime.TimeZoneCode <chr>,
## # ActivityDepthHeightMeasure.MeasureValue <dbl>,
## # ActivityDepthHeightMeasure.MeasureUnitCode <chr>,
## # ActivityDepthAltitudeReferencePointText <chr>,
## # ActivityTopDepthHeightMeasure.MeasureValue <chr>,
## # ActivityTopDepthHeightMeasure.MeasureUnitCode <chr>,
## # ActivityBottomDepthHeightMeasure.MeasureValue <chr>,
## # ActivityBottomDepthHeightMeasure.MeasureUnitCode <chr>,
## # ProjectIdentifier <int>, ActivityConductingOrganizationText <chr>,
## # MonitoringLocationIdentifier <chr>, ActivityCommentText <chr>,
## # SampleAquifer <chr>, HydrologicCondition <chr>, HydrologicEvent <chr>,
## # SampleCollectionMethod.MethodIdentifier <chr>,
## # SampleCollectionMethod.MethodIdentifierContext <chr>,
## # SampleCollectionMethod.MethodName <chr>,
## # SampleCollectionEquipmentName <chr>,
## # ResultDetectionConditionText <chr>, CharacteristicName <chr>,
## # ResultSampleFractionText <chr>, ResultMeasureValue <dbl>,
## # ResultMeasure.MeasureUnitCode <chr>, MeasureQualifierCode <chr>,
## # ResultStatusIdentifier <chr>, StatisticalBaseCode <chr>,
## # ResultValueTypeName <chr>, ResultWeightBasisText <chr>,
## # ResultTimeBasisText <chr>, ResultTemperatureBasisText <chr>,
## # ResultParticleSizeBasisText <chr>, PrecisionValue <chr>,
## # ResultCommentText <chr>, USGSPCode <chr>,
## # ResultDepthHeightMeasure.MeasureValue <chr>,
## # ResultDepthHeightMeasure.MeasureUnitCode <chr>,
## # ResultDepthAltitudeReferencePointText <chr>,
## # SubjectTaxonomicName <chr>, SampleTissueAnatomyName <chr>,
## # ResultAnalyticalMethod.MethodIdentifier <chr>,
## # ResultAnalyticalMethod.MethodIdentifierContext <chr>,
## # ResultAnalyticalMethod.MethodName <chr>, MethodDescriptionText <chr>,
## # LaboratoryName <chr>, AnalysisStartDate <date>,
## # ResultLaboratoryCommentText <chr>,
## # DetectionQuantitationLimitTypeName <chr>,
## # DetectionQuantitationLimitMeasure.MeasureValue <dbl>,
## # DetectionQuantitationLimitMeasure.MeasureUnitCode <chr>,
## # PreparationStartDate <date>, ProviderName <chr>,
## # ActivityStartDateTime <dttm>, ActivityEndDateTime <dttm>
Extracting the site information from temp_bbox
siteInfo <- attr(temp_bbox, "siteInfo")
head(as.tibble(siteInfo))
## # A tibble: 5 x 42
## station_nm agency_cd site_no dec_lat_va dec_lon_va hucCd
## <chr> <chr> <chr> <dbl> <dbl> <chr>
## 1 CLCJ-1 21AWIC 21AWIC-2228 33.6 -87.1 03160111
## 2 CLCJ-3 21AWIC 21AWIC-7015 33.6 -87.1 03160111
## 3 CLCJ-4 21AWIC 21AWIC-7016 33.6 -87.1 03160111
## 4 CHMJ-47 21AWIC 21AWIC-7033 33.6 -87.1 03160111
## 5 WIMJ-1 21AWIC 21AWIC-7034 33.6 -87.1 03160111
## # ... with 36 more variables: OrganizationIdentifier <chr>,
## # OrganizationFormalName <chr>, MonitoringLocationIdentifier <chr>,
## # MonitoringLocationName <chr>, MonitoringLocationTypeName <chr>,
## # MonitoringLocationDescriptionText <chr>, HUCEightDigitCode <chr>,
## # DrainageAreaMeasure.MeasureValue <chr>,
## # DrainageAreaMeasure.MeasureUnitCode <chr>,
## # ContributingDrainageAreaMeasure.MeasureValue <chr>,
## # ContributingDrainageAreaMeasure.MeasureUnitCode <chr>,
## # LatitudeMeasure <dbl>, LongitudeMeasure <dbl>,
## # SourceMapScaleNumeric <int>,
## # HorizontalAccuracyMeasure.MeasureValue <chr>,
## # HorizontalAccuracyMeasure.MeasureUnitCode <chr>,
## # HorizontalCollectionMethodName <chr>,
## # HorizontalCoordinateReferenceSystemDatumName <chr>,
## # VerticalMeasure.MeasureValue <chr>,
## # VerticalMeasure.MeasureUnitCode <chr>,
## # VerticalAccuracyMeasure.MeasureValue <chr>,
## # VerticalAccuracyMeasure.MeasureUnitCode <chr>,
## # VerticalCollectionMethodName <chr>,
## # VerticalCoordinateReferenceSystemDatumName <chr>, CountryCode <chr>,
## # StateCode <chr>, CountyCode <chr>, AquiferName <chr>,
## # FormationTypeText <chr>, AquiferTypeName <chr>,
## # ConstructionDateText <chr>, WellDepthMeasure.MeasureValue <dbl>,
## # WellDepthMeasure.MeasureUnitCode <chr>,
## # WellHoleDepthMeasure.MeasureValue <dbl>,
## # WellHoleDepthMeasure.MeasureUnitCode <chr>, ProviderName <chr>
Plot stations using leaflet
library.
library(leaflet)
m <- leaflet(siteInfo) %>%
addProviderTiles(providers$OpenStreetMap) %>%
addCircleMarkers(~dec_lon_va, ~dec_lat_va, color = "red", popup = ~station_nm)
m
I have just scratched the surface here and would like to encourage you guys to visit the full tutorial at https://owi.usgs.gov/R/dataRetrieval.html#1
.