BUGSrad: Radiative Transfer with Gases & Clouds



Atmospheric Scenario


Parameters


Clouds

Vertical Profiles: Control & Perturbation Runs

Control Atmosphere:

Vertical Profile

Control Global Properties


Modified Atmosphere:

Vertical Profile

Modified Global Properties


Download

Control Atmosphere:

Net (Downward) Radiative Fluxes

Control Heating Rates


Modified Atmosphere:

Net (Downward) Radiative Fluxes

Modified Heating Rates

MODEL DESCRIPTION

Basic Overview of the Project

In this project, you will use a simple “off-the-shelf” broadband radiation code developed at CSU to calculate shortwave (SW) and longwave (LW) fluxes and heating rates for model atmospheric profiles. In doing so, you will:

  • Gain understanding of what flux & heating rate profiles look like.

  • Explore how these flux & heating rates change due to changes in atmospheric consituents, such as changes in the gas amounts (water vapor, co2, ozone) and in clouds.

For more information about the underlying physics and numerics of the radiative transfer model, please see:

  • Fu, Q., and K-N. Liou, 1992: On the correlated k-distribution method for radiative transfer in nonhomogeneous atmospheres. J. Atmos. Sci, 49, 2139–2156.

  • Stephens, G  L., P  M. Gabriel, and P  T. Partain, 2001: Parameterization of atmospheric radiative transfer. Part I: Validity of simple models. J. Atmos. Sci., 58, 3391–3409.

A doubling of CO2 is expected to increase downwelling forcing at 200 mbar by about 5.3 W/m2. See e.g.

  • Collins, W. D., et al. “Radiative forcing by well‐mixed greenhouse gases: Estimates from climate models in the Intergovernmental Panel on Climate Change (IPCC) Fourth Assessment Report (AR4).” Journal of Geophysical Research: Atmospheres 111.D14 (2006).

There is currently some argument above what the global mean downwelling longwave flux is at the surface. Part of this argument stems from uncertainties in cloud measurements. See:

  • Stephens, Graeme L., et al. “The global character of the flux of downward longwave radiation.” Journal of Climate 25.7 (2012): 2329-2340.

The program that controls this website is surprisingly simple! It is all written in the programming language R. You can read all about it on the “Website Code” tab to the right.

How this website works (including all the code!)

This website is controlled using the R package “shiny.” There are three important components:

  • The user interface (ui.R) that controls the sliders and passes inputs to the model
  • The “server” (server.R) that responds to user interface events in the web browser
  • A fortran program (BUGSrad) that performs the broadband radiative transfer calculations
  • A module that takes input and runs BUGSrad (run.BUGSrad.R)
  • A module that reads input files and returns content for display on the website (modelInput.R)

Scroll down or click links in the list above to read all about it!


This is the user interface that actually controls this website

ui.R:

library(shiny)
library(markdown)

# Define user interface for atmospheric radiative transfer model
shinyUI(pageWithSidebar(

  #  Application title
  headerPanel("BUGSrad: Radiative Transfer with Gases & Clouds"),

  # Sidebar on the left with controls for the user to specify atmospheric properties
  sidebarPanel(

    br(),
    # Decide whether to plot with log of pressure for vertical coordinate
    checkboxInput('logp', 'Plot log(pressure)', value = FALSE),

    br(),
    # Choose a baseline atmospheric profile for the control run
    h4('Atmospheric Scenario'),
    selectInput("atmosphere", "McClatchey Profile:",
        list("Tropical" = "tropical",
             "Subarctic Winter" = "subarctic_winter")
             ),

    # Have the mu0 slider initial conditions be reactive to whichever atmosphere you have selected
    br(),
    uiOutput("mu0_slider"),

    # Vary the albedo
    br(),
    sliderInput("alpha", "Albedo",
        min = 0, max = 1.0, value = 0.3, step= 0.01),
    br(),
    # Vary concentrations of absorbing gases for perturbation run
    h4('Gases'),
    sliderInput("CO2.ppm", "CO2 (ppm)",
        min = 0, max = 2000, value = 380, step= 5),
    sliderInput("CH4.ppm", "CH4 (ppm)",
        min = 0, max = 5, value = 1.8, step= 0.1),
    sliderInput("vapor.mult", "Water Vapor (multiple of baseline)",
        min = 0, max = 3, value = 1, step= 0.1),
    sliderInput("O3.mult", "Ozone (multiple of baseline)",
        min = 0, max = 3, value = 1, step= 0.1),

    # Add liquid water clouds: base, top, mixing ratio, and cloud fraction
    br(),
    h4('Clouds'),

    checkboxInput("waterbutton", "Water Cloud"),
    conditionalPanel(
      condition = "input.waterbutton == true",
      sliderInput("water.cloud", "Water Cloud Base and Top (mb)",
        min = 200, max = 900, value = c(700,900), step= 50),
      sliderInput("cloud.liq.mix", "Water Mixing Ratio (g/kg)",
        min = 0, max = 2, value = 0, step= .01),
      sliderInput("cloud.liq.fraction", "Water cloud fraction (0 to 100%)",
        min = 0, max = 100, value = 0, step= 1)
    ),

    checkboxInput("icebutton", "Ice Cloud"),
    conditionalPanel(
      condition = "input.icebutton == true",
      sliderInput("ice.cloud", "Ice Cloud Base and Top (mb)",
        min = 100, max = 700, value = c(200,300), step= 50),
      sliderInput("cloud.ice.mix", "Ice Mixing Ratio (g/kg)",
        min = 0, max = 1, value = 0, step= 0.01),
      sliderInput("cloud.ice.fraction", "Ice cloud fraction (0 to 100%)",
        min = 0, max = 100, value = 0, step=.01)
    )
  ),

  # Show a tabset in the main panel of the browser that displays model output
  mainPanel(
    tabsetPanel(

      # First tab shows plots of model output for both control and modified atmosphere
      tabPanel("Profile Plots",
           h4('Vertical Profiles: Control & Perturbation Runs'), plotOutput('four.panel'),
           tags$style(type="text/css", ".tab-content { overflow: visible; }")),

      # Second tab displays tables of model input
      tabPanel("Model Input",
           h4('Control Atmosphere:'),
           h4('Vertical Profile'),
           tableOutput("control.profile"),
           h4('Control Global Properties'),
           tableOutput("control.global"),
           br(),
           h4('Modified Atmosphere:'),
           h4('Vertical Profile'),
           tableOutput("modified.profile"),
           h4('Modified Global Properties'),
           tableOutput("modified.global")
      ),

      # Third tab displays tables of model ouput
      tabPanel("Model Output",
           downloadButton('downloadData', 'Download'),
           h4('Control Atmosphere:'),
           h4('Net (Downward) Radiative Fluxes'),
           tableOutput("control.flux"),
           h4('Control Heating Rates'),
           tableOutput("control.rate"),
           br(),
           h4('Modified Atmosphere:'),
           h4('Net (Downward) Radiative Fluxes'),
           tableOutput("modified.flux"),
           h4('Modified Heating Rates'),
           tableOutput("modified.rate")
      ),

      # Fourth tab displays a brief model description
      tabPanel("Model Description",
           includeMarkdown('doc/model.description.md')),

      # Fifth tab displays website code
      tabPanel("Website code",
           includeMarkdown('doc/website.code.md'))
    )
  )
))

server.R:

library(shiny)

# Source required R scripts
source('model/run.BUGSrad.R')
source('model/plot.output.R')
source('model/modelInput.R')

shinyServer(function(input, output) {

  # Capture model input for later display
  model.input <- reactive({
    modelInput(input$atmosphere)
  })

  # Capture model output for later display
  model.output <- reactive({
    run.BUGSrad(alpha=input$alpha, mu0=input$mu0, CO2.ppm=input$CO2.ppm,
        CH4.ppm=input$CH4.ppm, O3.mult=input$O3.mult,
        vapor.mult=input$vapor.mult, logp=input$logp,
        cloud.base.liq=input$water.cloud[2],
        cloud.top.liq=input$water.cloud[1],
        cloud.liq.mix=input$cloud.liq.mix,
        cloud.liq.fraction=input$cloud.liq.fraction,
        cloud.base.ice=input$ice.cloud[2],
        cloud.top.ice=input$ice.cloud[1],
        cloud.ice.mix=input$cloud.ice.mix,
        cloud.ice.fraction=input$cloud.ice.fraction,
        atmosphere=input$atmosphere)
  })

  # 4-panel plot of model output
  output$four.panel <- renderPlot(
      plot.output(model.output()$control.flux,
          model.output()$control.rate,
          model.output()$modified.flux,
          model.output()$modified.rate,
          control.global=model.input()$control.global,
          alpha=input$alpha, mu0=input$mu0,
          CO2.ppm=input$CO2.ppm, CH4.ppm=input$CH4.ppm, O3.mult=input$O3.mult,
          vapor.mult=input$vapor.mult, logp=input$logp,
          cloud.base.liq=input$water.cloud[2],
          cloud.top.liq=input$water.cloud[1],
          cloud.liq.mix=input$cloud.liq.mix,
          cloud.liq.fraction=input$cloud.liq.fraction,
          cloud.base.ice=input$ice.cloud[2],
          cloud.top.ice=input$ice.cloud[1],
          cloud.ice.mix=input$cloud.ice.mix,
          cloud.ice.fraction=input$cloud.ice.fraction,
          atmosphere=input$atmosphere), height=800, width=600)

  # Display model input:
  output$control.profile <- renderTable(model.input()$control.profile)
  output$control.global <- renderTable(model.input()$control.global)  output$modified.profile <- renderTable(model.input()$modified.profile)
  output$modified.global <- renderTable(model.input()$modified.global)
  # Display model output:
  output$control.flux <- renderTable(model.output()$control.flux)
  output$control.rate <- renderTable(model.output()$control.rate)  output$modified.flux <- renderTable(model.output()$modified.flux)
  output$modified.rate <- renderTable(model.output()$modified.rate)

  # Have the mu0 slider initial conditions be reactive to whichever atmosphere you have selected
  output$mu0_slider <- renderUI({
    switch(input$atmosphere,
      "tropical" = sliderInput("mu0", "Solar Zenith Angle (deg)",
        min = 0, max = 90.0, value = 50., step= 1),
      "subarctic_winter" = sliderInput("mu0", "Solar Zenith Angle (deg)",
        min = 0, max = 90.0, value = 84., step= 1)
    )
  })

  # A button to download the output
  output$downloadData <- downloadHandler(
    filename = function() {
          paste('output.csv')
    },

    content = function(file) {
      write("Control Net (Downward) Radiative Fluxes", file, append = TRUE)
      write.table(model.output()$control.flux, file, sep = ",", row.names = FALSE, quote = FALSE, append = TRUE)
      write("Control Heating Rates", file, append = TRUE)
      write.table(model.output()$control.rate, file, sep = ",", row.names = FALSE, quote = FALSE, append = TRUE)
      write("Modified Net (Downward) Radiative Fluxes", file, append = TRUE)
      write.table(model.output()$modified.flux, file, sep = ",", row.names = FALSE, quote = FALSE, append = TRUE)
      write("Modified Heating Rates", file, append = TRUE)
      write.table(model.output()$modified.rate, file, sep = ",", row.names = FALSE, quote = FALSE, append = TRUE)
    }
  )


})

run.BUGSrad.R:

run.BUGSrad = function(logp=logp, alpha=alpha, mu0=mu0, CO2.ppm=CO2.ppm, CH4.ppm=CH4.ppm, O3.mult=O3.mult,
               vapor.mult=vapor.mult,
               cloud.base.liq=cloud.base.liq, cloud.top.liq=cloud.top.liq,
               cloud.liq.mix=cloud.liq.mix, cloud.liq.fraction=cloud.liq.fraction,
               cloud.base.ice=cloud.base.ice, cloud.top.ice=cloud.top.ice,
               cloud.ice.mix=cloud.ice.mix, cloud.ice.fraction=cloud.ice.fraction,
               atmosphere='tropical'){

  # read, modify, and write control file


  # Set filenames
  model.file <- 'model/bugsrad'

  control.file <- paste('data/mc_',atmosphere,'_bugsrad.txt',sep="")

  data.file <- paste('data/',atmosphere,'.data.txt',sep="")
  bottom.file <- paste('data/',atmosphere,'.bottom.txt',sep="")
  header.file <- paste('data/',atmosphere,'.header.txt',sep="")
  trailer.file <- 'data/trailer.txt'
  modified.data.file <- '/tmp/BUGSrad.modified.data.txt'
  modified.bottom.file <- '/tmp/BUGSrad.modified.bottom.txt'
  modified.file <- '/tmp/BUGSrad.modified.txt'

  control.output.file <- '/tmp/BUGSrad.control.out'
  modified.output.file <- '/tmp/BUGSrad.modified.out'

  #################

  # Read profile of control atmosphere
  profile <- read.table(data.file,
            col.names=c('layer','Pcenter','Pupper','temp','SpecHum',
                    'O3MixRatio','CloudWaterMix','CloudIceMix',
                    'RainMix','SnowMix','CloudFraction'),
            row.names=1, nrows=32)

  # Read global quantities for control atmosphere
  bottom <- read.table(bottom.file,
               col.names=c('Psurf','Tsurf','mu0','alvdr','alvdf','alndr',
                   'alndf','co2','ch4','n2o'), nrows=1)

  #################

  # Modify parameters as requested
  bottom$mu0 <- cos(mu0*pi/180)

  # Set a universal albedo
  bottom$alvdr <- alpha
  bottom$alvdf <- alpha
  bottom$alndr <- alpha
  bottom$alndf <- alpha

  bottom$co2 <- CO2.ppm
  bottom$ch4 <- CH4.ppm
  profile$SpecHum <- profile$SpecHum * vapor.mult
  profile$O3MixRatio <- profile$O3MixRatio * O3.mult

  # Add ice clouds
  profile <- within(profile,
            CloudIceMix[Pcenter <= 100*cloud.base.ice
                & Pupper >= 100*cloud.top.ice]
            <- cloud.ice.mix * .001)
  profile <- within(profile,
            CloudFraction[Pcenter <= 100*cloud.base.ice
                  & Pupper >= 100*cloud.top.ice]
            <- cloud.ice.fraction * .01)

  # Add liquid clouds
  profile <- within(profile,
            CloudWaterMix[Pcenter <= 100*cloud.base.liq
                  & Pupper >= 100*cloud.top.liq]
            <- cloud.liq.mix * .001)
  profile <- within(profile,
            CloudFraction[Pcenter <= 100*cloud.base.liq
                  & Pupper >= 100*cloud.top.liq]
            <- cloud.liq.fraction * .01)

  #################

  # Write modified data
  write.table(profile, file=modified.data.file, col.names=FALSE, quote=FALSE)
  write.table(bottom, file=modified.bottom.file, col.names=FALSE, row.names=FALSE)
  system(paste('cat',header.file,modified.data.file,trailer.file,
          modified.bottom.file,'>',modified.file))

  #################

  # Run radiative transfer for control atmosphere
  run.command <- paste(model.file, '<', control.file, '>', control.output.file)
  system(run.command)

  #################

  # Run radiative transfer for modified atmosphere
  run.command <- paste(model.file, '<', modified.file, '>', modified.output.file)
  system(run.command)

  #################

  # Read output files

  # read control fluxes from first half of BUGSrad output file
  control.flux <- read.table(control.output.file, skip=4,
                 col.names=c('level','p','SW_DN','SW_UP','LW_DN','LW_UP'),
                 row.names=1, nrows=33)

  # read control heating rates from second half of BUGSrad output file
  control.rate <- read.table(control.output.file, skip=39,
                 col.names=c('level','p','SW', 'LW'),
                 row.names=1, nrows=32)

  # read modified fluxes from first half of BUGSrad output file
  modified.flux <- read.table(modified.output.file, skip=4,
                  col.names=c('level','p','SW_DN','SW_UP','LW_DN','LW_UP'),
                  row.names=1, nrows=33)

  # read modified heating rates from second half of BUGSrad output file
  modified.rate <- read.table(modified.output.file, skip=39,
                  col.names=c('level','p','SW', 'LW'),
                  row.names=1, nrows=32)

  # Return model output to calling program
  return(list(control.flux=control.flux, control.rate=control.rate,
          modified.flux=modified.flux, modified.rate=modified.rate))

}


modelInput.R:

modelInput <- function(atmosphere){

  # Read the model input files and return contents for display on website

  control.profile.file <- paste('data/',atmosphere,'.data.txt',sep="")
  control.global.file <- paste('data/',atmosphere,'.bottom.txt',sep="")
  modified.profile.file <- '/tmp/BUGSrad.modified.data.txt'
  modified.global.file <- '/tmp/BUGSrad.modified.bottom.txt'

  # Read profile of control atmosphere
  control.profile <- read.table(control.profile.file,
            col.names=c('layer','Pcenter','Pupper','temp','SpecHum',
                    'O3MixRatio','CloudWaterMix','CloudIceMix',
                    'RainMix','SnowMix','CloudFraction'),
            row.names=1, nrows=32)

  # Read global quantities for control atmosphere
  control.global <- read.table(control.global.file,
               col.names=c('Psurf','Tsurf','mu0','alvdr','alvdf','alndr',
                   'alndf','co2','ch4','n2o'), nrows=1)

  # Read profile of modified atmosphere
  modified.profile <- read.table(modified.profile.file,
                col.names=c('layer','Pcenter','Pupper','temp','SpecHum',
                        'O3MixRatio','CloudWaterMix','CloudIceMix',
                        'RainMix','SnowMix','CloudFraction'),
                row.names=1, nrows=32)

  # Read global quantities for modified atmosphere
  modified.global <- read.table(modified.global.file,
                   col.names=c('Psurf','Tsurf','mu0','alvdr','alvdf','alndr',
                       'alndf','co2','ch4','n2o'), nrows=1)

  return(list(control.profile=control.profile, control.global=control.global,
     modified.profile=modified.profile, modified.global=modified.global))

}