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.