From 3508368298d77531d90a556ad70d3dea0461341d Mon Sep 17 00:00:00 2001 From: Alexander Mahura Date: Mon, 9 Oct 2023 17:20:18 +0300 Subject: [PATCH] add updated bash-scripts and fortran-files for synop STATS --- .../plots_for_one_station.f95 | 13 +- .../rank_histogram_summary_statistics.f95 | 4 +- .../rank_histogram_summary_statistics.f95~ | 380 ------------------ .../rank_histograms_one_station.exe | Bin 29688 -> 0 bytes .../rank_histograms_one_station.f95 | 14 +- .../produce_rank_histograms_all_stations.sh | 11 +- .../produce_standard_plots_all_stations.sh | 55 ++- .../summary_rank_histograms_all_stations.sh | 36 +- SYNOP/main_synop.sh | 2 + 9 files changed, 82 insertions(+), 433 deletions(-) delete mode 100644 SYNOP/STATS/fortran-programs/rank_histogram_summary_statistics.f95~ delete mode 100755 SYNOP/STATS/fortran-programs/rank_histograms_one_station.exe diff --git a/SYNOP/STATS/fortran-programs/plots_for_one_station.f95 b/SYNOP/STATS/fortran-programs/plots_for_one_station.f95 index 0f5cfe093..0b3626418 100644 --- a/SYNOP/STATS/fortran-programs/plots_for_one_station.f95 +++ b/SYNOP/STATS/fortran-programs/plots_for_one_station.f95 @@ -210,6 +210,7 @@ program plots_for_one_station integer,parameter :: n_colors=13 ! number of colours for GrADS output integer :: color(n_colors) ! colours for GrADS output + integer :: num_int_arr ! fix for Warning : Legacy Extension: REAL array index at (1) integer :: irec ! record number for GrADS output ! @@ -250,8 +251,8 @@ program plots_for_one_station ! f=miss open(unit=1,form='formatted',file=sim_file,status='old') - do while(.true.) -11 read(1,*,err=11,end=12)yyear,mmonth,dday,hhour,f1 + do while(.true.) +11 read(1,*,err=11,end=12) yyear,mmonth,dday,hhour,f1 if(l_round_6h)then call round_6h(yyear,mmonth,dday,hhour,year,month,day,hour) ! Rounding to nearest 6-h? else @@ -363,9 +364,13 @@ program plots_for_one_station write(1,'(A11)')'set cmark 3' write(1,'(A15)')'set digsiz 0.07' if(nyear.eq.1)then - write(1,'(A11,I2)')'set ccolor ',color(nint(1.+n_colors)/2.) + num_int_arr=nint(1.+n_colors)/2. + write(1,'(A11,I2)')'set ccolor ',color(num_int_arr) + !write(1,'(A11,I2)')'set ccolor ',color(nint(1.+n_colors)/2.) else - write(1,'(A11,I2)')'set ccolor ',color(nint(1+(i-1.)/(nyear-1.)*(n_colors-1.))) + num_int_arr=nint(1+(i-1.)/(nyear-1.)*(n_colors-1.)) + write(1,'(A11,I2)')'set ccolor ',color(num_int_arr) + !write(1,'(A11,I2)')'set ccolor ',color(nint(1+(i-1.)/(nyear-1.)*(n_colors-1.))) endif write(1,'(A7,I2,A1)')'d &0(z=',i,')' enddo diff --git a/SYNOP/STATS/fortran-programs/rank_histogram_summary_statistics.f95 b/SYNOP/STATS/fortran-programs/rank_histogram_summary_statistics.f95 index 5f129d6c6..46a8a608a 100644 --- a/SYNOP/STATS/fortran-programs/rank_histogram_summary_statistics.f95 +++ b/SYNOP/STATS/fortran-programs/rank_histogram_summary_statistics.f95 @@ -175,7 +175,7 @@ ! open(unit=1,form='formatted',status='old',file=infile) read(1,*) ! the first line in infile is a header - write(*,*)'Opened infile' + write(*,*)' Opened infile: ', infile do while(.true.) ! Loop continued until end of file reached if(l_code_in_char)then @@ -224,7 +224,7 @@ enddo 3 continue close(1) - write(*,*)'Read infile. n_line,n_station=',n_line,n_station + write(*,*)' Read infile: # of lines, # of stations: ',n_line,n_station ! if(n_line.ne.n_station*nhour_used)then ! write(*,*)'Something wrong !!!' diff --git a/SYNOP/STATS/fortran-programs/rank_histogram_summary_statistics.f95~ b/SYNOP/STATS/fortran-programs/rank_histogram_summary_statistics.f95~ deleted file mode 100644 index 743dc2712..000000000 --- a/SYNOP/STATS/fortran-programs/rank_histogram_summary_statistics.f95~ +++ /dev/null @@ -1,380 +0,0 @@ - program rank_histogram_summary_statistics -! -! Calculation of all-station rank histogram summary statistics: -! -! 1) distribution of p-values (in nbins_p bins) -! 2) Average frequencies in NQUANT+1 = 100 quantile bins -! (the quantile bin number is hard-coded!) -! -! The calculation is made for the UTC hours defined by -! the namelist parameters HOUR1, HOUR2, HOUR_STEP. It is assumed -! that the input file (INFILE) only includes data for these -! hours. -! -! Jouni Räisänen, University of Helsinki, August 2023 -! -!------------------------------------------------------------------------- -! -! INPUT FILE: -! -! infile: a text file including the rank histogram statistics for all stations -! -! The first line is a header. The other lines include: -! -! 1) station code (as a 3-character string, if l_code_in_char =.true. Else as an integer) -! 2) longitude -! 3) latitude -! 4) UTC hour (24 for 'all-UTC' statistics -! 5) total number of observations for this UTC hour -! 6) Mean squared deviation (MSD) from a flat quantile frequency histogram -! 7) p-value for MSD (= probability for getting at least as large MSD by chance) -! 8-107) frequencies in quantile space: 0-1%, 1-2% ... 98-99%, 99-100% -! -! Order of data assumed in the input file: -! -! line 1: header -! line 2: station 1, HOUR1 -! line 3: station 1, HOUR1 + HOUR_STEP -! line n = 1 + (HOUR2-HOUR1)/HOUR_STEP + 1 : station 2, HOUR2 -! line n+1: station2, HOUR1 -! etc. -!-------------------------------------------------------------------------- -! -! OUTPUT FILES: -! -! outfile_p_values: p-values at individual stations (This far: only the all-UTC p-values!) -! This file can be used for plotting the p-values on a map in python -! -! The first line is a header for eventual conversion to ODB. The other lines include -! -! 1) Station code as integer. If the original code was a 3-char string (GRUAN), -! the ordinal number of the station is written (for plotting in python) -! 2) longitude -! 3) latitude -! 4) p-value for MSD -!---------------------------------------------------------------------------- -! outfile_p_freq: frequency histogram of p-values (for different UTC hours separately) -! This file can be used for plotting a histogram of the p-value frequencies in python! -! -! The first line is a header. The others include -! -! 1) UTC hour -! 2-nbins_p+1) frequencies of p-values in nbins_p bins -!----------------------------------------------------------------------------- -! outfile_q = histogram of all-station-mean quantile bin frequencies (for different UTC hours separately) -! -! The first line is a header. The others include -! -! 1) UTC hour -! 2-101) quantile bin frequencies for 0-1%, 1-2% ... 99-100%. -!------------------------------------------------------------------------------ -! outfile_grads (only if grads_out=.true.): p-value frequencies and quantile bin frequencies in GrADS binary format -!------------------------------------------------------------------------------ -! text_output: frequencies of selected small p-values and some other data in free-text format, for each UTC hour separately. -!----------------------------------------------------------------------------------------------------------------------------- -! -! NAMELIST PARAMETERS: -! -! infile: see above -! outfile_p_values: see above -! outfile_p_freq: see above -! outfile_q: see above -! outfile_grads: see above -! text_output: see above -! nbins_p: number of p-value bins in outfile_p_freq -! grads_out: if .true. (.false.), outfile_grads is written (not written) -! -! hour1,hour2,hour_step: UTC hours in output = hour1, hour1+hour_step ... hour2 -! -! l_code_in_char: if .true., the station codes in infile are assumed to be 3-char string (else: integer) -! miss: missing value code. Default = -9.99e6. All values f with |f| >= |miss| are treated as missing. - - IMPLICIT NONE - INTEGER :: I,J ! loop variables - INTEGER,PARAMETER :: NHOUR=24 ! number of stations - INTEGER,PARAMETER :: NQUANT=99 ! number of quantiles from 1% to 99% - INTEGER,PARAMETER :: NPMAX=100 ! maximum bin number for frequency diagram of p-values - INTEGER,PARAMETER :: NSTATMAX=10000 ! maximum number of stations - REAL :: frequency_q1(nquant+1) ! quantile frequencies at a single station and UTC hour - REAL :: frequency_q(nquant+1,nhour+1) !all-station mean quantile frequencies for each UTC hour - REAL :: frequency_p(npmax,nhour+1) !frequencies of p-values for each UTC hours - REAL :: max_freq_p,max_freq_q ! maximum of p-value and quantile frequencies (for GrADS only) - REAL :: frequency_p01(nhour+1) ! frequency of p_values <= 0.1 % - REAL :: frequency_p1(nhour+1) ! same, <= 1 % - REAL :: frequency_p5(nhour+1) ! same, <= 5 % - character*2 :: number_of_bin ! two-digit code for quantiles for outfile_q and outfile_p_values headers - character*16 :: frequency_value ! frequency written in F16.6 format - character*160 :: infile,outfile_p_freq,outfile_p_values,outfile_q ! input and output files (see above) - character*160 :: outfile_grads,text_output ! output files (see above) - INTEGER,PARAMETER :: L_dataline=1700 ! maximum length of lines in text files - character*1700 :: dataline,headerline,emptyline ! string variables for input and output - - INTEGER :: nbins_p ! number of p-value bins in outfile_p_freq - INTEGER :: hour1,hour2,hour_step ! input / output UTC hours (see above) - INTEGER :: hour ! UTC hour - INTEGER :: n_line ! line count for lines read from infile - INTEGER :: n_station ! station count - INTEGER :: pbin ! index for p-value bin - - LOGICAL :: HOUR_USED(nhour+1) ! .true. for hours included in (hour1, hour1+hour_step,...,hour2) - INTEGER :: nhour_used ! number of UTC hours used - INTEGER :: station1_int ! the station code read from infile (RHARM, FMI) - CHARACTER*3 :: station1_char ! the station code read from infile (GRUAN) - LOGICAL :: l_code_in_char ! if .true., station codes in 3-char strings assumed - REAL :: longitude1,latitude1 !longitude and latitude - INTEGER :: station(nstatmax) !station codes of all stations - REAL :: longitude(nstatmax),latitude(nstatmax) ! station longitudes and latitudes - REAL :: p_value1, msd ! p-value and MSD read from infile - REAL :: p_value(nhour+1,nstatmax) ! p-values at individual stations for map - REAL :: total_number ! number of observations read from infile - - LOGICAL :: GRADS_OUT ! Output also as GrADS binaries, just for testing - INTEGER :: IREC ! record number for GrADS output -! - REAL :: MISS ! Missing value code - INTEGER :: valid_stations(nhour+1) ! stations with valid data - - NAMELIST/param/infile,outfile_p_values,outfile_p_freq,outfile_q,& - outfile_grads,text_output,& - nbins_p,grads_out,& - hour1,hour2,hour_step,& - l_code_in_char,miss - - miss=-9.99e6 ! default for missing value - READ(*,NML=PARAM) - valid_stations=0 -! -! Which UTC hours are used and what is their total number? -! - hour_used=.false. - nhour_used=0. - do i=hour1,hour2,hour_step - hour_used(i+1)=.true. - nhour_used=nhour_used+1 - enddo -! -! Empty data line -! - do i=1,L_dataline - emptyline(i:i)=' ' - enddo -!---------------------------------------------------------------------- -! -! Initialize counters etc. -! - frequency_p=0. - frequency_q=0. - frequency_p01=0. - frequency_p1=0. - frequency_p5=0. - n_station=0 - n_line=0 - p_value=miss -! -! Read the contents of the input file -! - open(unit=1,form='formatted',status='old',file=infile) - read(1,*) ! the first line in infile is a header - write(*,*)'Opened infile' - - do while(.true.) ! Loop continued until end of file reached - if(l_code_in_char)then - 1 read(1,*,err=1,end=3)& - station1_char,longitude1,latitude1,hour,total_number,msd,p_value1,(frequency_q1(i),i=1,nquant+1) - else - 2 read(1,*,err=2,end=3)& - station1_int,longitude1,latitude1,hour,total_number,msd,p_value1,(frequency_q1(i),i=1,nquant+1) - endif -! write(*,*)'n_station,Hour step:',n_station,hour_step - if(hour_used(hour+1))then ! Only use the statistics for the UTC hours define by hour1, hour2, hour_step - n_line=n_line+1 - if(mod(n_line,nhour_used).eq.1)then ! This assumes that all the required UTC hours are always included in infile - n_station=n_station+1 - if(l_code_in_char)then - station(n_station)=n_station - else - station(n_station)=station1_int - endif - longitude(n_station)=longitude1 - latitude(n_station)=latitude1 - endif - if(total_number.gt.0.and.abs(frequency_q1(1)).lt.abs(miss))then ! Only include stations with some valid data - p_value(hour+1,n_station)=p_value1 ! for map of p-values - valid_stations(hour+1)=valid_stations(hour+1)+1 - do i=1,nquant+1 - frequency_q(i,hour+1)=frequency_q(i,hour+1)+frequency_q1(i) ! update the quantile bin frequencies - enddo - pbin=min(1+(p_value1*nbins_p),real(nbins_p)) - frequency_p(pbin,hour+1)=frequency_p(pbin,hour+1)+1 ! update the p-value bin frequencies -! -! Frequencies of small p-values: -! - if(p_value1.le.0.001.and.abs(p_value1).lt.abs(miss))then - frequency_p01(hour+1)=frequency_p01(hour+1)+1. - endif - if(p_value1.le.0.01.and.abs(p_value1).lt.abs(miss))then - frequency_p1(hour+1)=frequency_p1(hour+1)+1. - endif - if(p_value1.le.0.05.and.abs(p_value1).lt.abs(miss))then - frequency_p5(hour+1)=frequency_p5(hour+1)+1. - endif - endif - - endif - enddo -3 continue - close(1) - write(*,*)'Read infile. n_line,n_station=',n_line,n_station - -! if(n_line.ne.n_station*nhour_used)then -! write(*,*)'Something wrong !!!' -! stop -! endif -!--------------------------------------------------------------------- -! Divide all the frequencies by the number of stations -! - do hour=hour1,hour2,hour_step - do i=1,nquant+1 - frequency_q(i,hour+1)=frequency_q(i,hour+1)/valid_stations(hour+1) - enddo - do i=1,nbins_p - frequency_p(i,hour+1)=frequency_p(i,hour+1)/valid_stations(hour+1) - enddo - frequency_p01(hour+1)=frequency_p01(hour+1)/valid_stations(hour+1) - frequency_p1(hour+1)=frequency_p1(hour+1)/valid_stations(hour+1) - frequency_p5(hour+1)=frequency_p5(hour+1)/valid_stations(hour+1) - enddo -! -! Find the normalized maximum frequencies for GrADS axis limits -! - max_freq_p=0 - max_freq_q=0 - do i=1,nbins_p - if(frequency_p(i,nhour+1).gt.max_freq_p/nbins_p)then - max_freq_p=frequency_p(i,nhour+1)*nbins_p - endif - enddo - max_freq_p=1.02*max_freq_p - - do i=1,nquant+1 - if(frequency_q(i,nhour+1).gt.max_freq_q/(nquant+1))then - max_freq_q=frequency_q(i,nhour+1)*(nquant+1) - endif - enddo - max_freq_q=1.02*max_freq_q -! -!--------------------------------------------------------------------- -! -! Write the text output file: -! 1) number of stations -! 2) frequencies of selected small p-values - - open(unit=1,form='formatted',file=text_output) - write(1,'(A20,I3)')'Number of stations: ',n_station - write(1,'(A24,I3)')'Number of p-value bins: ',nbins_p - if(grads_out)then - write(1,'(A23,F6.3)')'Axis_limit_p_for_GrADS ',max_freq_p - write(1,'(A23,F6.3)')'Axis_limit_q_for_GrADS ',max_freq_q - endif - do hour=hour1,hour2,hour_step - write(1,'(A22,I2.2,A6,F5.3)')'Frequency, p <= 0.001 ',hour,' UTC: ',& - frequency_p01(hour+1) - write(1,'(A22,I2.2,A6,F5.3)')'Frequency, p <= 0.01 ',hour,' UTC: ',& - frequency_p1(hour+1) - write(1,'(A22,I2.2,A6,F5.3)')'Frequency, p <= 0.05 ',hour,' UTC: ',& - frequency_p5(hour+1) - enddo - close(1) -! -! Write the average quantile bin frequencies and p-value bin frequencies -! to ODB compatible text files - -!-------------------------------------------------------------------------- -! -! Open the quantile bin frequency output file and write its first line. -! - open(unit=2,form='formatted',status='unknown',file=outfile_q) - headerline=emptyline - headerline='hour@hdr:integer ' - do j=0,nquant - write(number_of_bin,'(i2.2)')j - headerline=trim(headerline)//' f'//number_of_bin//'@body:real' - enddo - write(2,*)trim(headerline) -! -! Write the data lines for each UTC hour used -! - do hour=hour1,hour2,hour_step - j=hour+1 - write(dataline,'(I6)')& - hour - do i=1,nquant+1 - write(frequency_value,'(F16.6)')frequency_q(i,j) - dataline=trim(dataline)//frequency_value - enddo - write(2,*)trim(dataline) - enddo - close(2) -! -! Write the file giving the p-values for individual stations -! - open(unit=2,form='formatted',status='unknown',file=outfile_p_values) - headerline=' station@hdr:integer longitude@hdr:real latitude@hdr:real p_value@body:real' - write(2,*)trim(headerline) - do i=1,n_station - write(2,'(7X,I6,3F16.6)')& ! this far: just fwrite the station ordinal number, in case - ! the station code was a character string - station(i),longitude(i),latitude(i),p_value(nhour+1,i) - enddo - close(2) -! -! Open the p_value frequency output file and write its first line. -! - open(unit=2,form='formatted',status='unknown',file=outfile_p_freq) - headerline=emptyline - headerline='hour@hdr:integer ' - do j=1,nbins_p - write(number_of_bin,'(i2.2)')j-1 - headerline=trim(headerline)//' p'//number_of_bin//'@body:real' - enddo - write(2,*)trim(headerline) -! -! Write the data lines for each UTC hour used -! - do hour=hour1,hour2,hour_step - j=hour+1 - write(dataline,'(I6)')& - hour - do i=1,nbins_p - write(frequency_value,'(F16.6)')frequency_p(i,j) - dataline=trim(dataline)//frequency_value - enddo - write(2,*)trim(dataline) - enddo - close(2) -! -!------------------------------------------------------------------------------- -! -! Open the file for GrADS output and write its contents, for easier visualisation? -! - if(grads_out)then - open(unit=11,form='unformatted',file=outfile_grads,access='DIRECT',& - recl=4,status='unknown') -! - do hour=hour1,hour2,hour_step - j=(hour-hour1)/hour_step+1 - do i=1,nbins_p - irec=(j-1)*(nbins_p+nquant+1)+i - write(11,rec=irec)frequency_p(i,hour+1) - enddo - do i=1,nquant+1 - irec=(j-1)*(nbins_p+nquant+1)+nbins_p+i - write(11,rec=irec)frequency_q(i,hour+1) - write(*,*)hour,i,frequency_q(i,hour+1) - enddo - enddo - close(11) - - endif ! if (grads_out) - - END program rank_histogram_summary_statistics diff --git a/SYNOP/STATS/fortran-programs/rank_histograms_one_station.exe b/SYNOP/STATS/fortran-programs/rank_histograms_one_station.exe deleted file mode 100755 index 68af780d96994881fe44a9f93119733f5a2562a7..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 29688 zcmeHwdw5jU+3(&n2}A>#aFbB!=d?uyYC=H4prA8Aa7P0~NHlH1AtaMaNoLZ?3=%(u z*f7@JVKnUnEp2V}X|3|LUt7MOVu3211SGMw8m(1qX{8ozB37wbh*Zt_z3Z}P&&=Kd zJUxBRA7{brwchJ`*X>=GJ+rd8w0x;URTUW=%5@4s?y(+*;!j1{_<05jT#0gl;!!SB zE>%VWlZB6mQ7Gmq84n57GR_fvPQ)Zv;^2r9hleF(Je8vr8A}T(&K5Ktc?85`&oY4( zv69bWk|SGp3jEGB9!8LH2^oZN_#bn5LdcPMB)3P%?GbV^b_@GuEXzk_qfe!vS2{d^ z2=-8hhX_>c~1FOugJWZxQX;?6~y>ehwN7=f4dJ|6G|_)uA<;6wSz50d8`o`MkRN)ON2 z;CI;I15jY4=Q}odvek`GdicDJ{@>Z)kJ;c)+299}p_QFQHh8a%o@qAtN*nwx8~xAO z;E&koxyuGmuGFer-8OjYH&*=BHu$AB_#fKn|B?;A;`24M1OvXda7`%e^C`ZD`d}y= zstNe$(N^yd`9l7hX5a0hrm){ctqlfhYr+r__i?d*QH>Bt(JUk{;)65((Kz% z6G~xc3Ve0JfH^u8~-7Nt29j0EcHmDb?x^-7Ds z1+}WQ)HF8-YYA!zZh$>LIO$r^#4R;V0i{0V_bWik#wLG#6(L$EHII>MAwke=g1!)` zscWeTHxiJ@8 z)BG#LS`(h;OBotXcxne3T1|MGuVvV1!khhLhY3$}ybSl6@U%XVq05AqYY3w4FyU$4 zAj3l@JgqHc*lEJkxx zJomF84PUw%PuC?LqzG;#`ry%FgiBWvOkHoVpTRd0>>xP7;6(&eR~y{L;A;q`E;jfO zgXa)TU2Cw5!7~Y_u0Pnp;7<}vU2L$G!5<@-y4K)Y29GD0hQwe6gD)hQy4awH!486{ zYYi4N_{=K6)TIX941R}T>Pmz841R-P>OzAGgI^_>y3XK<4}rM!_XJay89d71=Lx2+ zGT6`HUlB}QWH7Na>LKm955toT4Y^4T zij?nII0g19k;zq>{@TDMRM?^H(ZHx{4lw;E57gHexh~m4R-i*6EB&=K$90OlmD=5J zYUnwd{>0%AhKDtMcR!1^_pvCRKvb?r_Yf3|5)@Ok9hzQxSc{b&-bWhIhc%bhQ>GvA zKF(mt*O^cmYd=`~Q(3#c9F9dPpOdcq9DK6>vH#{Pz3f8{Fb4h6N| z=m1E?sQ=#R38eZl3VNeQk?KI+z0qSxQAgI|45{h817m>_;JyhU?XRHL1L0>$0)dkN z?gw4@B|QcX5j7U2>LB}&6s9V$WW$)=Qd;gRgHPj;WTT`}N%0C&Af(20Q2{0zbX1OB zdTcNb@!rxmNtHp;O9x`5Z-QZd>A)s`tm@5J>44g=_jJE;iQ50%kY4&UScP3m4A2hE zBCCrm1m8UVhy`}wcvlzC5EjKsj~V)?N1*IuKov+p*~_j~l#{M%Y9`JwE~~3cyZaz@ zBk$@p2dH1c#kFesU3)Y=vP;vec53>vhsyQR9p(D6u05=gKK&o8N?zd6Q63HO=pY>hpf~M_DV;~1(M04j2ITC!BzpX2z3=oJ-Dk4(zV4w* zX7uR&>l1Y|y7j(jV)Kj)uF(cpQ=$v@aK^~!-qM2~ms#XPYKKnGT`_ugvGh)s&S&X| zSlS&^D7%xcnTdMW+fPx|xTg0(pg4wlnwdE1x~H8qKnHJ>w8B^jFDv&(ad~OXJhn)jzY)LdRzh~1zqW7}(2vS{Agc}VI z^}T5dQ4AhbB3H``iOi#k_5fI>vkXjmqZ=D}gDFxbb0C6|0|+7m^7Npv42Xqi03*<| z$l^+P45qQG>nYY05=mQXw{+IRq?us1riPY!zJ7sJEfCH=k8 zwMf=7Pp}#(n$F-8z^L2}a(s6OiN+bmySErr@1BS7+)LeWU!(8={X@jGyFW1kZ&2_) z%z3+>7=ib8!Ml<39_=21H!66`Iq$ZECSICsz0B+HrFnbb4#DddHjKpk9l`tVX;Qd* z>qxv03*J7?yU|>>bQ^vwcpu=r;gNWMCU{#p?`m^fr0d-$c$aYA>XCT61#bc8E&sKt z&S-&PuCqSD`{pTD55*(3@$-T=!Fe~DZAjO9RPb)&yy1~}-xj=$oR@WH*(VFq?0sMG zE&y+D^j-91>eHeN(+stk6`J*E0xKD~-HrC;J%x=s{iT!I4|_{rfT$9SQ>ThO(M31XslMmlMr@0=N{I?~yrYPHX6 zoo}eUTUeenjyg+UP}Rs0xN$6Gv-JMsFLw2o(s_w;_gzOG#zgFTFyVSE0oQO*tn>)R zzP97gmh9GdhKJpq;5+_u7X~{U#gAGu#q%H#J=U5_c5v}v6mUb6XC zt?AIi7!Awf6V-hjduGh$Yc;)JQxlrnk2=x1|C+6-L(Ae>4mr)#Zk}*sd;vL=tS6Kj z51(iLs-ZZ?f2nmIQjdR6>+DvKZ^`g7Kc{@-Q8wK5s2eV(H#!y1hX~&rosOh>qX#~x z7%5t81L;9Biy<|nXocJxJ&Y6*(mCZSDsNm?0^S#>HE($rZ4{-B7frafL_p23XuqLl zD@#!qNz>x~2le&1o*bL|m}&Yml*c;c0cvxy$L=HVqH$`NQ(LX6WFKVLT2rezb?SS> z#CoV%^EI5hBk76E)Ty9~UK3-YnA zOuj{qJ5LS|vmW*%)A;!mGD_w!q`J0z1Izo)A%|tQ*A2!FPv&q=}&U1O{G5y>c6_uuY}HwO8+sCHkIx{eG3~> z8v1t8&>95Ee7=!r>L;%di?ZVG(3%qD3$*xEj#!k{Pc0mP;R$-sH{c2SwfMwQTpQAb z$F{N6#g;hr9`XoL@(5(|;D2J<)WtPXhlqaReKMQ!>KvNtimn4_W$GPgHB}h@ueEj~ zSv&Zff1>ZCFOHry`YRG}Iu=girW4fXqH8gE?#sU%J0jvF|M1~MC^a>VWn258-(ndMr*RcM z<+^tjIldqsu_!HMFpbxhb1z>WkKV~_j`M(r9o!@h(A;)1?bh}B8 zmZgxp|8C&0Xo*vM(oUu53OcmtLvlNO@cPFu&RAwL;@0a^ zri3~<_q(DOfM7P7cLJ(Mbb|Z7ol?uSeVX3K%hQMTfhZ5h=-pyE=)+djkJj4C{GhdJ z+P)_#=cATw`7y(6VcL(Eel#t7X>6!r6V+>2vS6a3v8A5igCo9ZIQo2ewB9q=O+@N~ zyJoQcb7{XA`@7i|;M+7(_4u9Wr_n#Oe-e#!f-UubufKHq6~5d*?bQ#?f5vtH25P0L z*dOrCKq_@>HERvh!>6~Ax)99awSYvNh6u*R9W(KWTpEieRiHt^!1e;Hz$bdpmHN@H z$!XWr9Ufy1E~WX@Z`Xm1`JXqqH0HKaoY~}$-87wQYIS#}Ee}M=1{m}0H{qR$_crhj zv<#pK(`mqGuVk{Utuw@{qg)DIfYh+IE8&X_GFv^2(T&&}UovgDGmo4--Xam1~-Dsg*j8}|Cqe-G9 z##-d-Y&3nedObK&y*`@KL2*K>!O7Y~*{tVT2Z|PxJi#6Q##4j^fei~f01N-HyjO9l$`XuH7cF6T_SomG_8k^(M zbTscjMwzEC1Lp#HmSy0F(1Hz@%?{basWzLnm7Hp`StI+Pc`K$x)C#TW%)8uqppM{- z|HIxHZ@%+*fwd#PbDi1kkQ2tN3}-XK#`ykYREcDqvCICr=p5-DDH|T?_;2BCRz4$EwxlpS{&ImW|yO`}cogHE+{u4{)kYt8L;` zn^x-p_5XoZd->IWZL8_fXKJ+zsgy>iNNF`U+NQJguxd5;@I6HPC6;p!E%pAJJmi#& zr5;FI*}SuzG$aoe$q#O){u|3AY~4j@J7J*}MY_qDPSr;yBF_l;2bW+&DZ&oL!Mr7bV=_dhEeRs2kgp)`0a%>#Y`L z*vjmtIrVTR)mUiS7LJZ(|BEIBQ3$9JQ!`BtREP{WT1BHY9%XMuY#jSUEc*G;Cw)%< z8qVkF6P|$d%T4d+UZLK&{ULM)VCi&tFCAZN`b*eToYcNgFT?dSu^j!v~fA>!Le>&rvs~TY#GNE16u;@A`W_F zMx6<5?`oWCYx~_IXnL{!D-G}PI(u?7J#+$ZG;_4hSBJD2{c3mTpNH_Calol}v!a~9 z!TjJ8DRvLyg(mENgCfK3V;uP!LmGDfiDO^kSZViRj(wS9rQOeS(Bpb+(z8h3TcIUO zP`3AT|3@V-8fteVpXll+(x2N$$Ne+~{OC{QObFYVn5^yKC-w9lfh9dezVDB;YdomM zwv#+(jKzsP-ags8_kYl@ynE=o7tpEqlq%&n)vJ(80@0 z@-ZeJixYdieY|(?i>%D_ogJgh4?}|Ui88Z{QRejqAznTHY-Qe0wx36tfBbu1=24LA zWqT3Wgrf$w{TVXeKDeTEQkLN8Qq&c{uJLNgY3$!&;ksZ3*ZD z%zcv2L~j&JR;2d8EG~HXl*?WQ?q&uNSu`GFLAAX95>qE< z=MMl{ErHqU&;AtYNjgPxp?6!JA`x}#-|>FiI2hwqlAPebN+TY3bRGS?}GpAC0Sq$X2kT35S#=A{A z)f{09_jslRt&K;Jcq3DS*2^PE>|wQGxR>{Shox@Zsztu=5nE)^O|*U_2WZW*dy_!y_!?D>W`)uu&$*ovKrL$Ymu3m5@%_V>P(5Vw8&-5W}`(KUNFCq zG@XndwrVBzb%J=bg1#$i95+Aytx;R&cnlk__!h9^25XMIWY(173aA{UU%q8@l`nq8iaiHU{ex$z0Ct{w z3@fGG+_{`nZ7Q}{)R1)kM#cV|XK7Qh-(Xbs;>D(7U*a4iRP1h{16fnWo_fZ{;r`UB zvZK7R|8F?lg)oG_0f=|eE1R?X(&LaY`O+&GD1GVBC;d175FsB zF;}|4(v{weCOeBO?SZDGD;=OJi${-9WU9FzQCYS#Y%DH#Os@1gR#YO3I$$iQgD6vO zOz`&*n7;n~vSA~8@*m;B-*@tL~g8kqVK_806e_22O2M$7#LQMWz| zDnBueMd>>!^sbwK4fY$SKjiz3j2hm78s>|?^@%3_a$Ap5nZ)7;waqv4X{C7dSJWBo zoLmzvc9Pq-9}shT+kPW8XKc{XIOmXaAxY-WVj;G@Dsu^JqJ zfaOSh&}t<9AkDJ#9gH_3xAR-_4e=&>L6*^UKZ1Q?eq%E?oj9h=Z3gc$HWGpiDniK7)Q>n4@)n z;2ih~^%?9h-LwB7f8tuMpTVYkiT<{Cc;XSQ^MI<&|6}A0`i=GH)^Ov!KJ!mv)3Wb6 zc>vCo^ow6*tHzBqC<_8E$5`)cY-$S!8$vZLZN6Z@PZuoERS!k=^XDozM`{A$re=Se zyQQhE4YxIX;`&zHD$&vus0UUFM#3z{y$X!Y$Do$BI^VirFpOI@YFas}F&GKWR_Hm0 zJ&}pORcZIv087s~3U2QRH_qnK90k_Y!J0@|X5nKPrP)^-tn>St0zSIB0wB9K!bge~ zvQP;&*Wu<5-0=|h*P(p1_(THh1Hsz^R9r=r&C^&Hx`to+;cgBF8k)k9IzL0wB^&N$ zR6r_*GGKIZv4vTWKz%et?r<=Sdr1P3mUaG+XI-$a{TgNh>d*wW;&zedh(C#$-h>+} zX1ME381I!=%vm~nUeUZOW>;TPHs^}6s|7q`df7Z%5|fz0?jBlFwz70_Wh7v(za^+M z_JbhA7+CCsNV$X_2pEE!c}hMzu5S$%{QAO|KD7WCkD~M~QZII1s{z3BU3nizwfOd~ zvK8){wGX=d8}N;k$@efnPDgwOaXI2UK7btJiQu82-my|qH|DDsj~SiQ zrH;-gJl+l~Zs`3qF^FKzuPD<%C;qjtrwe+Oy!@qk1vk2KZ_nvauK(DzGq0TXNmAj) zhw9)#=w)|CdGqq`aV#D?s%)#%d7}flb3KjV`7J!dEZV_H9_LKPKzZGZ&#Q2|d*_v$ke;j&eRS7nfER%VGD8Ohvi>1U*@# z=ieMojP?X)ePV_$N`t4R!E~1jg)ux*T5~-tCGYx}SjJ&_w}ref zjSlW9P#>l6F&-(q)_7Qo`U(X}mvy$o;UTp0j}M1~6jwMs6d8BA4P|qfHklM!1)bjA zQK%92(hNnRfJestU;@OAFgy(SJS-u^6pu*Aeoj*q1$oy;@dX|ROS~drS&m{x<#&Ej zQs$3S1y-J}91r8bIA^f58q+%mAcF6qPiLGtBdC1DkB;2YYj2O^A?u(9)4Ofq{Y;)9ZC_NvPYyc zEMW8&gD%IDyu(J0CpjKcr4@zSbH5=k%k$NAMgjBuZa}Z*_Iy#a6a9LaLV>6U33F`n zITJ}fByyMwk>$V1RKGN~OYu=V$jfH2IVDHE8R;>TAWP)c_&De&L7hAw$qSrrk*_l= z8`Hg#H3wLQE^=@vH{hA0b`k<%4}_8#{0vZ*@)_WB)Gg3GW}L_%y&eX@;>gK<5&X_q z@NgVn0!-Fb3_W@(q_TO=V0w<|`B}#?7DufI1ARs-IXPcIe8D|1+?YTMULw9wZP<_h9Y$)OTG@(>J==RizZczC^@Ph5Bg9hFCre3`yufaKZ$hi zb);xECC3&~(nlA{P>@|heIPGPgvE5l421^v1m=PNh8oxB8Ieq;;Ch{`Sxk196ZjtG z#K0W@^Xb=J>MhqJc?_R%-9qK104UHg!T2Ucy_G6%G=K?z6dV)PD#u;QM77ef3D27y zv}&0o#mNrkqD|0EAN6X{vh+SUN4*D659m`-pk&Lc>TIxIM3yX=63FA$Sk_cEuug@1pQ`>dZ!Dg4TOlnS)S{wU zZS8G!{?@ixxXIHu58Z*?3|O>^uPiQbNXOv3PIcxNE&9}T%cf0>Y;fj>8y&clP5l(p zjjv2*>Q0ulVw7s&B$a{q%9<_E^Vl-6=f}Ki?sRqJ*Km+p4AknZN}zCQT)yIg-O%f!gY22 z`l3czA98Ry4=Ad0e4jHAc8+(mXV@87;2eMB4XYJLmh&}w^nTSbeTwtKB4+WL^R0qF zvI?YkcV4dE#Ht;|mmc2vX@SdV%UQB1XEg>}{IhD-MVgxHW`$Z?cw$jiS^1J#4RpC* z;hdr?XV0o_i_o2{Ztp5LVBGMB^y1mG`ON7qYzapK{)G+xfIrk!3$P|s+qjV5|LQJm zVEkDg}YuXwMTSLLHzn1pjkZP?7W19;Vq>v#C{ef_( zec^0u@~&bzhmbF>Ni4)gmQWK2Q)X;WjJ-O_5^e;{^$X3x5Q+jrP_lZW7y3h?KoCG6 z2m%>mFajBFZYylTj=Omwb)BY~=0g8QKWeGIqO4T$m94-&z7CslbhtXcx2`P=`Eh9@ z#z$tmiT+O;Znjx=(^T5lbcf#;#^sKb0d9y4)wJ^~9u6o z>~q#VNZ_>evH`*;c+@62WK zAMsQ%I_Nx`!c_KDl)oXi(r;I9L*Q3&)ANF59-$!_S?L+YDrB5-|GVU#0=$aT3VE&} z@v|8Jgp7Lg2z+||tYUg5De30|_hCKn#;Ha^w2RbFKinTL0aM&;gQstDjYmDFZ*K0l z;itC{R{Cl0W`+MM>?HkiS^8;V(_13n)1tr0d+V|EQrsrHWs>%nR{C|)<8Uc*yh?sL zue9P{Brh ze`ML_0B>dgQX72G2H$0azZU(9^6gl)L1xL~AKCCfV}m~dd_Ml-fb^eI*=fd~r_jB$ z6xg{q;_XxqjxomhE+f(d>yB_!auJFpXInUeJ7a5=I)BLDfHz&Z4?W!C!xiTNKQ1H3 z?diUT=HNQKgQ*J#Lv6mA$VR0$*wWhU#~Y-gdGqJbw`MZ0OWC3vw~^Ph_yi zm>~1gmb0Giwc9dp*>c=sFJG&*6-Uz38c1UrQ*Axe8nPfxo^eNrK?ISeU;ut<)y@ES_EgTuuzvFtjw2uX1RA+ z*5HY$?M9|Dl~<(n1o8|G7rkbI#?x@Z^4wU(Uytm(MOp*ngx#V9Zux3GN*e(6m)XcvDkizGyT(;QgDs^o^Ut>+6j@H~M z%_X{9ilo6{)k2YB9NwVIKM!GT7a9yY>gq_7EFD2DPbD#Kw#gTiQyIEYY+>e>U)$K!B_ zQ#SZRZ8(}p0eqlg&tF3VBGuX)R*KlhxhU-4h&1itK?^cti~NnkJK`<8VnAGi@ASEZ zOhYlgaexzrh9Gp|QN%J}%CJq8ha^6ENa7GL3qaN&1r%LRM%pCZCURAH#6y7OHdHEi-0`oitS^iNRG~@-SEJq?2A1N>2v$hKP zdj)~YNcAl9FF=|C?VO2K%FFk*LsTdP`JR?Tne~4g2-xcy_RII8pA{KV{Zf$j%Y0Bp$+t=y+Lheud)8wg*D5T4eN9Kf+<(K=YqeA{6asI=E@o>w8h$lmy>@QBkOL;lYC7F0B z$k;99Io9-~*Xk4`BZrJ}D?o+Ld(v1qX(f-E(#B6ZFYQCqaU|JP$>lrKh-siwW+VS5 z9pEASHBBH5ROpr&!ll>WZ*OJd%G|U8E9GUpR>*ViuHB7JF9%<%m}zF0W+X~0rnko`~g3-XV$?=eC7sN;$EZOSj%#vx%P|9&Iq LXiqC}8mRnlF)P#< diff --git a/SYNOP/STATS/fortran-programs/rank_histograms_one_station.f95 b/SYNOP/STATS/fortran-programs/rank_histograms_one_station.f95 index c73168c0d..2540ae9b0 100644 --- a/SYNOP/STATS/fortran-programs/rank_histograms_one_station.f95 +++ b/SYNOP/STATS/fortran-programs/rank_histograms_one_station.f95 @@ -186,12 +186,10 @@ if(l_code_in_char)then ! The err=1 specifier ensures that header lines are skipped 1 read(1,*,err=1,end=3)& - station_char,longitude,latitude,month,day,hour,& - (q1(i),i=1,nquant) + station_char,longitude,latitude,month,day,hour,(q1(i),i=1,nquant) else 2 read(1,*,err=2,end=3)& - station_int,longitude,latitude,month,day,hour,& - (q1(i),i=1,nquant) + station_int,longitude,latitude,month,day,hour,(q1(i),i=1,nquant) endif do i=1,nquant quantiles(i,month,day,hour+1)=q1(i) @@ -217,7 +215,7 @@ open(unit=1,form='formatted',status='old',file=infile) do while(.true.) -11 read(1,*,err=11,end=12)yyear,mmonth,dday,hhour,value1 +11 read(1,*,err=11,end=12) yyear,mmonth,dday,hhour,value1 ! ! Rounding of hours? ! @@ -305,12 +303,10 @@ do while(.true.) if(l_code_in_char)then 21 read(1,*,err=21,end=23)station_char,longitude,latitude,realization,& - (msd_bootstrap(j+1),j=hour1,hour2,hour_step),& - msd_bootstrap(nhour+1) + (msd_bootstrap(j+1),j=hour1,hour2,hour_step),msd_bootstrap(nhour+1) else 22 read(1,*,err=22,end=23)station_int,longitude,latitude,realization,& - (msd_bootstrap(j+1),j=hour1,hour2,hour_step),& - msd_bootstrap(nhour+1) + (msd_bootstrap(j+1),j=hour1,hour2,hour_step),msd_bootstrap(nhour+1) endif ! ! Update the p-value counters diff --git a/SYNOP/STATS/produce_rank_histograms_all_stations.sh b/SYNOP/STATS/produce_rank_histograms_all_stations.sh index b8b064ed2..29ea15fff 100755 --- a/SYNOP/STATS/produce_rank_histograms_all_stations.sh +++ b/SYNOP/STATS/produce_rank_histograms_all_stations.sh @@ -40,6 +40,7 @@ echo "Bash version ${BASH_VERSION}..." # 1. Variable code # variable=$1 +echo " Meteorological variable code = $1" # # The following variables are included: # @@ -168,6 +169,8 @@ echo " List of synop stations (with geographical coordinates): $station_list" # Compile the Fortran program that calculates the rank histograms # echo " Compiling fortran-program to calculate rank histograms ..." +echo " fortran-programs/rank_histograms_one_station.f95" + gfortran fortran-programs/rank_histograms_one_station.f95 -o rank_histograms_one_station # ################################################################# @@ -186,6 +189,7 @@ while [ ${line_number} -le `expr $n_lines` ] do head -`expr ${line_number}` ${station_list} | tail -1 > input.txt read station longitude latitude < input.txt +echo " **********" echo " Synop station ID, longitude, latitude: ${station} ${longitude} ${latitude}" #echo ${station} @@ -193,17 +197,16 @@ echo " Synop station ID, longitude, latitude: ${station} ${longitude} ${latitu #echo ${latitude} # ################################################################ -# # Select the simulation data for the station (mimicked this far by observations!). # When simulation data are available, a means to retrieve the correct data points from it must be designed +echo " Selecting the simulation data (mimicked his far by observations) for synop station: ${station}" # odb_command="odb sql 'select year,month,day,hour,value where station=${station} and (year>=${year1} and year<=${year2}) and (hour=0 or hour=6 or hour=12 or hour=18) and variable=${variable}' -i ${sim_file} -o sim_data" eval ${odb_command} # ############################################################### -# # Select the quantiles for the station: -# +echo " Selecting the quantiles for synop station: ${station}" # 00, 06, 12 and 18 UTC are picked separately, since # select \* does not allow for parentheses (very strange) # @@ -224,9 +227,9 @@ cat quantile_selection_* > quantile_selection odb sql select \* where station=${station} -i ${msd_bootstrap_file} -o msd_bootstrap_selection # ############################################################### -# # Produce the rank histogram for one station. # Include data for 00, 06, 12 and 18 UTC. +echo " Producing the rank histogram for one synop station: ${station}" # ############################################################## # diff --git a/SYNOP/STATS/produce_standard_plots_all_stations.sh b/SYNOP/STATS/produce_standard_plots_all_stations.sh index 4398a3658..b45052afc 100755 --- a/SYNOP/STATS/produce_standard_plots_all_stations.sh +++ b/SYNOP/STATS/produce_standard_plots_all_stations.sh @@ -1,5 +1,7 @@ #!/bin/bash echo "Bash version ${BASH_VERSION}..." +echo "Python version" +which python ########################################################################## # # Produce the following plots: @@ -33,8 +35,10 @@ echo "Bash version ${BASH_VERSION}..." # their relevant parts are converted to text format for processing. # # Execution (e.g): ./produce_standard_plot_all_stations 39 2020 2020 1 12 -# +# ORIGINAL: # Jouni Räisänen, August 2023 +# MODIFIED: +# Alexander Mahura, Sep-Oct 2023 # ################################################################## # @@ -43,7 +47,8 @@ echo "Bash version ${BASH_VERSION}..." # 1. Variable code # variable=$1 -# +echo " Meteorological variable code = $1" + # The following variables are included: # # 91 'total amount of clouds' @@ -168,11 +173,13 @@ fi # and is available as a .txt file. # station_list=list_of_stations_${variable}.txt +echo " List of synop stations (with geographical coordinates): $station_list" # # ################################################################ # echo " Compiling the Fortran program needed for creating the figures ..." +echo " fortran-programs/plots_for_one_station.f95" # gfortran fortran-programs/plots_for_one_station.f95 -o plots_for_one_station # @@ -192,24 +199,26 @@ while [ ${line_number} -le `expr $n_lines` ] do head -`expr ${line_number}` ${station_list} | tail -1 > input.txt read station longitude latitude < input.txt +echo " **********" echo " Synop station ID, longitude, latitude: ${station} ${longitude} ${latitude}" # ################################################################ -# # Select the simulation data for the station (mimicked this far by observations!) # When simulation data are available, a means to retrieve the correct data points from it must be designed. +echo " Selecting the simulation data (mimicked his far by observations) for synop station: ${station}" # odb_command="odb sql 'select year,month,day,hour,value where station=${station} and (year>=${year1} and year<=${year2}) and (hour=0 or hour=6 or hour=12 or hour=18) and variable=${variable}' -i ${sim_file} -o sim_data" eval ${odb_command} + ################################################################ -# # Select the rank histogram data for the station +echo " Selecting rank histogram data for synop station: ${station}" # odb sql select \* where station=${station} -i ${rh_file} -o rank_histograms_${variable}_${station}_${year1}${month1}-${year2}${month2} # ############################################################### -# -# Select the quantiles for the station. +# Select the quantiles for the station +echo " Selecting the quantiles for synop station: ${station}" # # 00, 06, 12 and 18 UTC are picked separately, since # select \* does not allow for parentheses (very strange) @@ -245,8 +254,8 @@ cat quantile_selection_* > quantile_selection l_code_in_char=.false.,l_round_6h=.false. &END EOF -########################################################### -# + +################################################################ # Produce the figure in Python # cat < standard_plot.cfg @@ -268,26 +277,24 @@ rank_hists=rank_histogram_${variable}_${station}_${year1}${month1}-${year2}${mon [output] fig_name=${figure_dir}/standard_plot_${variable}_${station}_${year1}${month1}-${year2}${month2}_python.png EOF -######################## + +################################################################ # -#which python -echo " Calling python to plot quantiles rank histogram ..." -#echo Calling Python -#python3 python/plot_quantiles_rankhist.py standard_plot.cfg +echo " Calling python to plot quantiles rank histogram for synop station: ${station}" python3 python/plot_quantiles_rankhist.py standard_plot.cfg -#echo Called Python -####################### -# + +################################################################ # Remove unnecessary files -# +echo " Removing unnecessary temporary files ..." + rm input.txt rm vrange_* rm msd_and_p-value -rm msd_bootstrap_selection rm quantile_selection rm rank_histograms_${variable}_${station}_${year1}${month1}-${year2}${month2} rm sim_data -rm standard_plot_one_station +###rm msd_bootstrap_selection +###rm standard_plot_one_station rm time_series_commands rm time_series_${variable}_${station}_${year1}${month1}-${year2}${month2}.grads rm rank_histogram_${variable}_${station}_${year1}${month1}-${year2}${month2}.grads @@ -298,7 +305,11 @@ rm rank_histogram_${variable}_${station}_${year1}${month1}-${year2}${month2}.txt # ((line_number++)) done -################################################ + +################################################################ # The Fortran executable is not needed any more: -############################################## -rm plots_for_one_station +# +rm plots_for_one_station +rm quantile_selection_* +rm msd_and_p-value_* +rm coordinates diff --git a/SYNOP/STATS/summary_rank_histograms_all_stations.sh b/SYNOP/STATS/summary_rank_histograms_all_stations.sh index 7f038235c..ce8e6e465 100755 --- a/SYNOP/STATS/summary_rank_histograms_all_stations.sh +++ b/SYNOP/STATS/summary_rank_histograms_all_stations.sh @@ -1,5 +1,7 @@ #!/bin/bash echo "Bash version ${BASH_VERSION}..." +echo "Python version" +which python ################################################################## # # Calculation and plotting of summary statistics for quantile @@ -23,8 +25,10 @@ echo "Bash version ${BASH_VERSION}..." # (for the time range specified by the script arguments) # # Execution (e.g): ./summary_rank_histograms_all_stations 39 2020 2020 1 12 -# +# ORIGINAL: # Jouni Räisänen, August 2023 +# MODIFIED: +# Alexander Mahura, Sep-Oct 2023 # ################################################################## # @@ -33,6 +37,7 @@ echo "Bash version ${BASH_VERSION}..." # 1. Meteorological variable code # variable=$1 +echo " Meteorological variable code = $1" # # The following variables are included: # @@ -95,12 +100,11 @@ fi # ################################################################## # -# Compile the Fortran program that produces the rank histogram summary statistics +echo " Compiling the Fortran program that produces the rank histogram summary statistics ..." +echo " fortran-programs/rank_histogram_summary_statistics.f95" # gfortran fortran-programs/rank_histogram_summary_statistics.f95 -o rank_histogram_summary_statistics -echo "COMPILED - program that produces the rank histogram summary statistics" - # ################################################################## # @@ -108,7 +112,9 @@ echo "COMPILED - program that produces the rank histogram summary statistics" # # 1) Rank histogram directory and input name file name without extension # +echo " Directory for rank histogram ..." rh_dir=rank_histograms/${variable}_${year1}${month1}-${year2}${month2} +echo " $rh_dir" rh_file=${rh_dir}/rank_histograms_${variable}_${year1}${month1}-${year2}${month2}_all_stations rh_summary_file=${rh_dir}/rank_histograms_${variable}_${year1}${month1}-${year2}${month2}_summary # @@ -118,7 +124,10 @@ out_file=${rh_dir}/rh_summary_${variable}_${year1}${month1}-${year2}${month2} # # 3) Directory for figures # +echo " Directory for figures ..." figure_dir=figures +echo " $figure_dir" + # ################################################################## # @@ -129,7 +138,7 @@ odb sql select \* -i ${rh_file}.odb -o ${rh_file}.txt # ################################################################## # -# Calculate the rank histogram summary statistics +echo " Calculating the rank histogram summary statistics ..." # ./rank_histogram_summary_statistics <