correct_annual_mean_chemistry Subroutine

public subroutine correct_annual_mean_chemistry()

Arguments

None

Calls

proc~~correct_annual_mean_chemistry~~CallsGraph proc~correct_annual_mean_chemistry correct_annual_mean_chemistry proc~uemep_annual_mean_pdf_correction_no2_o3 uEMEP_annual_mean_pdf_correction_NO2_O3 proc~correct_annual_mean_chemistry->proc~uemep_annual_mean_pdf_correction_no2_o3 cosd cosd proc~uemep_annual_mean_pdf_correction_no2_o3->cosd proc~get_sun_angles get_sun_angles proc~uemep_annual_mean_pdf_correction_no2_o3->proc~get_sun_angles proc~date_to_julian date_to_julian proc~get_sun_angles->proc~date_to_julian proc~date_to_number date_to_number proc~date_to_julian->proc~date_to_number sngl sngl proc~date_to_number->sngl

Called by

proc~~correct_annual_mean_chemistry~~CalledByGraph proc~correct_annual_mean_chemistry correct_annual_mean_chemistry program~uemep uEMEP program~uemep->proc~correct_annual_mean_chemistry

Source Code

    subroutine correct_annual_mean_chemistry()
        integer :: i, j, t
        integer :: t_start, t_end
        integer :: i_integral, j_integral
        real :: o3_in, nox_in, no2_in, J_photo_in, temperature_in, lon_in, lat_in
        real :: ox_sigma_ratio_in, nox_sigma_ratio_in
        real :: o3_out, no2_out
        real :: sum_no2_in, sum_no2_out
        integer :: no2_count
        logical :: run_all_flag

        write(unit_logfile,'(A)') ''
        write(unit_logfile,'(A)') '================================================================'
        write(unit_logfile,'(A)') 'Correcting annual mean NO2 and O3 (correct_annual_mean_chemistry)'
        write(unit_logfile,'(A)') '================================================================'

        t_start = 1
        t_end = subgrid_dim(t_dim_index)

        sum_no2_in = 0.0
        sum_no2_out = 0.0
        no2_count = 0
        do t = t_start, t_end

            !Run the conversion routine once to get the Jd distribution which is saved. This is to save time as this is slow. Done in the centre of the grid
            if (quick_annual_mean_pdf_chemistry_correction) then
                run_all_flag = .false.
                i = subgrid_dim(x_dim_index)/2
                j = subgrid_dim(y_dim_index)/2
                o3_in = comp_subgrid(i,j,t,o3_index)
                no2_in = comp_subgrid(i,j,t,no2_index)
                nox_in = comp_subgrid(i,j,t,nox_index)
                i_integral = crossreference_target_to_integral_subgrid(i,j,x_dim_index)
                j_integral = crossreference_target_to_integral_subgrid(i,j,y_dim_index)
                J_photo_in = meteo_subgrid(i_integral,j_integral,t,J_subgrid_index)
                temperature_in = meteo_subgrid(i_integral,j_integral,t,t2m_subgrid_index)
                lon_in = lon_subgrid(i,j)
                lat_in = lat_subgrid(i,j)
                ox_sigma_ratio_in = ox_sigma_ratio_pdf
                nox_sigma_ratio_in = nox_sigma_ratio_pdf
                call uEMEP_annual_mean_pdf_correction_NO2_O3(min_bin_pdf, max_bin_pdf, log10_step_bin_pdf, .true., no2_in, nox_in, o3_in, J_photo_in, &
                    temperature_in, ox_sigma_ratio_in, nox_sigma_ratio_in, lon_in, lat_in, no2_out, o3_out)
            else
                run_all_flag = .true.
            end if

            do j = 1, subgrid_dim(y_dim_index)
                do i = 1, subgrid_dim(x_dim_index)

                    o3_in = comp_subgrid(i,j,t,o3_index)
                    no2_in = comp_subgrid(i,j,t,no2_index)
                    nox_in = comp_subgrid(i,j,t,nox_index)
                    i_integral = crossreference_target_to_integral_subgrid(i,j,x_dim_index)
                    j_integral = crossreference_target_to_integral_subgrid(i,j,y_dim_index)
                    J_photo_in = meteo_subgrid(i_integral,j_integral,t,J_subgrid_index)
                    temperature_in = meteo_subgrid(i_integral,j_integral,t,t2m_subgrid_index)
                    lon_in = lon_subgrid(i,j)
                    lat_in = lat_subgrid(i,j)

                    ox_sigma_ratio_in = ox_sigma_ratio_pdf
                    nox_sigma_ratio_in = nox_sigma_ratio_pdf

                    if (o3_in .le. 0 .or. nox_in .le. 0 .or. no2_in .le. 0) then
                        o3_out = o3_in
                        no2_out = no2_in
                    else
                        call uEMEP_annual_mean_pdf_correction_NO2_O3(min_bin_pdf, max_bin_pdf, log10_step_bin_pdf, run_all_flag, no2_in, nox_in, o3_in, J_photo_in, &
                            temperature_in, ox_sigma_ratio_in, nox_sigma_ratio_in, lon_in, lat_in, no2_out, o3_out)
                    end if

                    comp_subgrid(i,j,t,o3_index) = o3_out
                    comp_subgrid(i,j,t,no2_index) = no2_out
                    sum_no2_in = sum_no2_in + no2_in
                    sum_no2_out = sum_no2_out + no2_out
                    no2_count = no2_count + 1

                    if (isnan(no2_out)) then
                        write(unit_logfile,'(a,2i5,4f12.4)') 'NaN in pdf output. Stopping ',i,j,no2_in,no2_out,o3_in,o3_out
                        stop 1
                    end if
                end do
            end do
        end do

        write(unit_logfile,'(a,f12.4)') 'Average NO2 scaling with pdf correction = ',sum_no2_out/sum_no2_in
    end subroutine correct_annual_mean_chemistry