[ncl-talk] Fwd: problem with plotting "t-test" pattern with precipitation change
Moustapha Tall
tall.moustapha89 at gmail.com
Mon Oct 26 12:22:40 MDT 2015
---------- Forwarded message ----------
From: Moustapha Tall <tall.moustapha89 at gmail.com>
Date: 2015-10-26 16:17 GMT+00:00
Subject: problem with plotting "t-test" pattern with precipitation change
To: ncl-talk-request at ucar.edu
Hi NCL webmasters and all,
I am tryning to plot the mid-century precipitation change for
a given model, thus plotting the rainfall change patterns with "stippling"
where
they are statistically significant at the 95 percent confidence level
(based on a t-test).
I get this message:
(0) Error: scalar_field: If the input data is 1-dimensional, you must
set sfXArray and sfYArray to 1-dimensional arrays of the same length.
warning:create: Bad HLU id passed to create, ignoring it
My script is below,
;**************************************************************************************
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
begin
;*****************************************************************************************
; Creating Function to Read at once Precip and Convert Calendar to Human
Readable YYYYMM
;*****************************************************************************************
undef("getPrecip")
function getPrecip(fname[1]:string, vname[1]:string)
begin
f = addfile(fname, "r")
x = f->$vname$(0,:,:) ; read specified time period
return(x)
end
;*****************************************************************************************
; Reading Precip for each File for the Defined Period ymStrt - ymLast
;*****************************************************************************************
gcmpre1 = getPrecip("remap_reg_jjas_mm_mmean_EC-EARTH_1976-2005.nc","pr")
gcmfut1 =
getPrecip("remap_pr_wafr_jjas_EC-EARTH_rcp45_2021-2050.nc","pr")
;*******************************************************************************
;Taking JJAS difference
;*******************************************************************************
diff1 = gcmfut1 - gcmpre1
diff1 = 100*diff1/mask(gcmpre1,(gcmpre1.ge.0.9),True)
diff1 at long_name = " RCP4.5 - Present JJAS Precipitation (%)"
diff1 at units = " "
printVarSummary(diff1) ; print a summary of the variable
copy_VarCoords(gcmpre1,diff1)
;************************************************************************************************************
; Create graphic array for panel plot and Open a Work Station
;************************************************************************************************************
plots = new(1,graphic) ; Open a variale
plots to store 4 Figures
wks = gsn_open_wks("ps","Precip_JJAS_Change") ; open a
workstation
gsn_define_colormap(wks, "CBR_drywet")
;*************************************************************************************************************
; Create plots options
;*************************************************************************************************************
res = True ; plot options desired
res at gsnSpreadColors = True ; span full colormap
res at gsnDraw = False ; don't draw
res at gsnFrame = False ; don't advance frame
res at cnFillOn = True ; turn on color fill
res at cnLevelSelectionMode = "ExplicitLevels" ;
set explicit contour levels
res at cnLevels = (/-75.,-50.,-25.,-10.,-5.,5.,10.,25.,50.,75./)
res at cnLinesOn = False ; turn off contour lines
res at cnLineLabelsOn = False ; turn off the values that
would be on the contours
res at mpFillOn = False ; turn off gray continents
res at gsnMaximize = True ; Use the full page
res at gsnPaperOrientation = "portrait" ; Change orientation
res at mpProjection = "Mercator"
res at mpLimitMode = "LatLon"
res at gsnAddCyclic = False ; Regional data
; res at mpMaxLatF = max(latitude) ; Taking the whole
domain
; res at mpMinLatF = min(latitude)
; res at mpMaxLonF = max(longitude)
; res at mpMinLonF = min(longitude)
res at mpMaxLatF = 25.0 ; Zoom in West Africa
res at mpMinLatF = -5.0
res at mpMaxLonF = 25.0
res at mpMinLonF = -25.0
res at mpOutlineBoundarySets = "National" ; Display the
conutries borders
res at mpOutlineOn = True ; Turn on the outline of
the map
res at mpGeophysicalLineThicknessF = 3.0 ; Increase the thickness
of map border
res at pmTickMarkDisplayMode = "Always" ; Turn on map
tickmarks, i.e. writing lons and lats along the axes
res at lbLabelBarOn = False ; turn off
individual colorbars
;*********************************************************************************************************************
; Changing the font seize
;*********************************************************************************************************************
res at tmXBLabelFontHeightF = 0.02 ; resize tick labels
for X axis
res at tmYLLabelFontHeightF = 0.02 ; resize tick labels
for Y axis
; change label spacing
to avoid overlap
;res at lbLabelStride = 1 ; and write every
other label
res at lbLabelFontHeightF = 0.018 ; change the size of
the label bar
;res at tiMainString = "1st Plot Test" ; Provide a main Tile
res at gsnStringFontHeightF = 0.025 ; Increase the title
font size (long_name)
;*******************************************************************************************************
; setting and plotting the t-test significance
;*******************************************************************************************************
aveX = dim_avg_n_Wrap(gcmpre1,0) ; average over the 0th dim
varX = dim_variance_n_Wrap(gcmpre1,0) ; compute variance
printVarSummary(aveX)
printVarSummary(varX)
aveY = dim_avg_n_Wrap(gcmfut1,0) ; average over the 0th dim
varY = dim_variance_n_Wrap(gcmfut1,0) ; compute variance
printVarSummary(aveY)
printVarSummary(varY)
sX = dimsizes (aveX)
sY = dimsizes (aveY) ; different sizes
alphat = 100.*(1. -ttest(aveX,varX,sX, aveY,varY,sY, False, False))
;;;;;;;;;;;;;creating ressources for T-TEST;;;;;;;;;;;;;;;;;;;;;;;;;;;;
rest_test = True ; plot options desired
rest_test at gsnSpreadColors = True ; span full colormap
rest_test at gsnDraw = False ; don't draw
rest_test at gsnFrame = False ; don't advance frame
rest_test at cnFillOn = True ; turn on color fill
rest_test at lbLabelBarOn = False ; turn off individual
colorbars
rest_test at cnLinesOn = False ; turn off individual
colorbars
rest_test at cnLineLabelsOn = False ; turn off individual
colorbars
rest_test at cnInfoLabelOn = False ; turn off individual
colorbars
rest_test at cnLevelSpacingF = 0.05
rest_test at cnFillDotSizeF = 0.004
t_test_plot = gsn_csm_contour(wks,alphat,rest_test)
opt_t_test = True
opt_t_test at gsnShadeFillType = "Pattern"
opt_t_test at gsnShadeLow = 17
plot_t_test_plot =
gsn_contour_shade(t_test_plot,0.05,999.,opt_t_test)
;*******************************************************************************************************
; Displaying the Panel Plot
;*******************************************************************************************************
plots(0) = gsn_csm_contour_map(wks,diff1,res)
overlay(plots(0), plot_t_test_plot)
pnlres = True ; Add panel
options
pnlres at gsnFrame = False ; Don't frame
yet
;pnlres at txString = "Default X and Y positions of all plots" ; add main
title for the panel
pnlres at gsnPanelLabelBar = True ; Setting a
common colorbar
pnlres at lbLabelFontHeightF = 0.02 ; Increase
font size of the numbers in the label bar
pnlres at pmLabelBarWidthF = 0.9 ; label bar
width
pnlres at pmLabelBarHeightF = 0.08 ; label bar
height
pnlres at pmLabelBarOrthogonalPosF = 0.002 ; move the
labelbar up and down
pnlres at pmLabelBarParallelPosF = 0.03 ;move the
labelbar left or right
pnlres at lbLabelStride = 1 ; and write
every other label
; pnlres at lbLabelAutoStride = True ; Optimal
labels, i.e. avoid overlap
pnlres at gsnPanelDebug = True ; To get
information from panel
pnlres at gsnPanelFigureStrings= (/"a)"/) ; add strings to panel
pnlres at gsnPanelFigureStringsFontHeightF = 0.016 ; Increase
the size of these strings
pnlres at gsnPanelFigureStringsPerimOn = False ;
Remove the box around these strings
pnlres at amJust = "TopLeft" ; Where to put
these strings
pnlres at gsnPanelYWhiteSpacePercent = 5 ; Adds the
white space to the panel plot
pnlres at gsnPanelXWhiteSpacePercent = 5 ; Adds the
white space to the panel plot
pnlres at gsnPanelTop = 1.
pnlres at gsnPanelBottom = 0.
gsn_panel(wks,plots,(/1,1/),pnlres)
end
;**************************************************************************************
and
Variable: diff1
Type: float
Total Size: 33600 bytes
8400 values
Number of Dimensions: 2
Dimensions and sizes: [70] x [120]
Coordinates:
Number Of Attributes: 3
units :
long_name : RCP4.5 - Present JJAS Precipitation (%)
_FillValue : 1e+20
Variable: aveX
Type: float
Total Size: 480 bytes
120 values
Number of Dimensions: 1
Dimensions and sizes: [lon | 120]
Coordinates:
lon: [-29.75..29.75]
Number Of Attributes: 12
associated_files : baseURL:
http://cmip-pcmdi.llnl.gov/CMIP5/dataLocation gridspecFile:
gridspec_atmos_fx_EC-EARTH_historical_r0i0p0.nc areacella:
areacella_fx_EC-EARTH_historical_r0i0p0.nc
history : 2012-05-10T11:28:56Z altered by CMOR: replaced missing value
flag (1e+28) with standard missing value (1e+20). 2012-05-10T11:28:56Z
altered by CMOR: Converted type from 'd' to 'f'. 2012-05-10T11:28:56Z
altered by CMOR: Inverted axis: lat.
cell_methods : time: mean (interval: 3 hours)
original_name : TP=LSP+CP
comment : at surface; includes both liquid and solid phases from all
types of clouds (both large-scale and convective)
missing_value : 1e+20
_FillValue : 1e+20
units : kg m-2 s-1
long_name : Precipitation
standard_name : precipitation_flux
time : 7699.5
average_op_ncl : dim_avg_n over dimension(s): lat
Variable: varX
Type: float
Total Size: 480 bytes
120 values
Number of Dimensions: 1
Dimensions and sizes: [lon | 120]
Coordinates:
lon: [-29.75..29.75]
Number Of Attributes: 12
associated_files : baseURL:
http://cmip-pcmdi.llnl.gov/CMIP5/dataLocation gridspecFile:
gridspec_atmos_fx_EC-EARTH_historical_r0i0p0.nc areacella:
areacella_fx_EC-EARTH_historical_r0i0p0.nc
history : 2012-05-10T11:28:56Z altered by CMOR: replaced missing value
flag (1e+28) with standard missing value (1e+20). 2012-05-10T11:28:56Z
altered by CMOR: Converted type from 'd' to 'f'. 2012-05-10T11:28:56Z
altered by CMOR: Inverted axis: lat.
cell_methods : time: mean (interval: 3 hours)
original_name : TP=LSP+CP
comment : at surface; includes both liquid and solid phases from all
types of clouds (both large-scale and convective)
missing_value : 1e+20
_FillValue : 1e+20
units : kg m-2 s-1
long_name : Precipitation
standard_name : precipitation_flux
time : 7699.5
variance_op_ncl : dim_variance_n over dimension(s): lat
Variable: aveY
Type: float
Total Size: 480 bytes
120 values
Number of Dimensions: 1
Dimensions and sizes: [lon | 120]
Coordinates:
lon: [-29.75..29.75]
Number Of Attributes: 12
associated_files : baseURL:
http://cmip-pcmdi.llnl.gov/CMIP5/dataLocation gridspecFile:
gridspec_atmos_fx_EC-EARTH_rcp45_r0i0p0.nc areacella:
areacella_fx_EC-EARTH_rcp45_r0i0p0.nc
history : 2012-03-19T15:13:27Z altered by CMOR: replaced missing value
flag (1e+28) with standard missing value (1e+20). 2012-03-19T15:13:27Z
altered by CMOR: Inverted axis: lat.
cell_methods : time: mean (interval: 3 hours)
original_name : TP=LSP+CP
comment : at surface; includes both liquid and solid phases from all
types of clouds (both large-scale and convective)
missing_value : 1e+20
_FillValue : 1e+20
units : kg m-2 s-1
long_name : Precipitation
standard_name : precipitation_flux
time : 10987
average_op_ncl : dim_avg_n over dimension(s): lat
Variable: varY
Type: float
Total Size: 480 bytes
120 values
Number of Dimensions: 1
Dimensions and sizes: [lon | 120]
Coordinates:
lon: [-29.75..29.75]
Number Of Attributes: 12
associated_files : baseURL:
http://cmip-pcmdi.llnl.gov/CMIP5/dataLocation gridspecFile:
gridspec_atmos_fx_EC-EARTH_rcp45_r0i0p0.nc areacella:
areacella_fx_EC-EARTH_rcp45_r0i0p0.nc
history : 2012-03-19T15:13:27Z altered by CMOR: replaced missing value
flag (1e+28) with standard missing value (1e+20). 2012-03-19T15:13:27Z
altered by CMOR: Inverted axis: lat.
cell_methods : time: mean (interval: 3 hours)
original_name : TP=LSP+CP
comment : at surface; includes both liquid and solid phases from all
types of clouds (both large-scale and convective)
missing_value : 1e+20
_FillValue : 1e+20
units : kg m-2 s-1
long_name : Precipitation
standard_name : precipitation_flux
time : 10987
variance_op_ncl : dim_variance_n over dimension(s): lat
(0) Error: scalar_field: If the input data is 1-dimensional, you must
set sfXArray and sfYArray to 1-dimensional arrays of the same length.
warning:create: Bad HLU id passed to create, ignoring it
Can anyone help me for that!
cheers
--
--
Tall Moustapha,
Laboratoire de Physique de l'Atmosphère et de l'Océan Siméon-Fongang
(LPAO-SF)
École Supérieure Polytechnique (ESP) / Université Cheikh Anta DIOP (UCAD)
BP: 5085 Dakar fann ; Fax : (+221) 33 825 93 64
Tel perso : (+221) 77 730 80 11
Dakar / Sénégal
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20151026/4ff84488/attachment.html
More information about the ncl-talk
mailing list