grid_ID formula giving bad results

[nbest937]

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates grid_ID in this commit but I must not be using the floor() function correctly because when I go into the debugger I see that the first term in the grid_ID expression returns zero:

(gdb) print grid_ID
$7 = -651839
(gdb) print 90. -lat( ilat) *12.
$8 = -150.50001525878906
(gdb) print floor( 90. -lat( ilat) *12.)
$9 = 0
(gdb) print nint( 90. -lat( ilat) *12.)
No symbol "nint" in current context.
(gdb) print floor( -150.5)
$10 = 0 
(gdb) print ceiling(( 180. +lon( jlon)) *12.)
No symbol "ceiling" in current context.

As a result the data is all written to a file called literally ***/*******/GENERIC1.WTH because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try nint() or ceiling()? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

[dmcinern]

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90 but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

[joshuaelliott]

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 **

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 **

+! when longitiude runs from 0 to 360:

279 **

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

  280 **

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 **

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90 and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.com wrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit < https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

[joshuaelliott]

hey david. sorry for the barrage of emails, but i wanted to add one more piece of info. i realized that there was a typo in the new grid_ID formula regarding the order of operations in the argument for floor(). it should be (90-lat)_12 rather than 90-lat_12, which means that the argument should never actually run negative anyway. since the same is true for the argument of ceiling(), we don't have to worry about that weird behaviour (though i'm still curious why its happening).

On Fri, Dec 28, 2012 at 10:29 AM, joshua elliott jelliott@ci.uchicago.eduwrote:

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 **

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 **

+! when longitiude runs from 0 to 360:

279 **

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

  280 **

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 **

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.comwrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit < https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

[joshuaelliott]

yet another update on the status of the f90 debugging:

currently, the issue we're trying to figure out has to do with garbage data that appears in the first 24 days of the .WTH files (see the attached example). i can't figure out what is happening but i figure it must have something to do with an array being defined as an inappropriate size or accessed inappropriately. i thought this might have something to do with max_points_all, since there was a compiler error associated with this early on and neil had to bump up the value somewhat to avoid it, but i think it probably has more to do with how min_land_point_proc is getting set, and how that effects ilat and jlon.

! copy data from netcdf file to big array of data ! do counter = 1, n_land_points_proc

  do counter = chunk_start,chunk_end
    ilat = print_point_lat(counter+min_land_point_proc-1)
    jlon = print_point_lon(counter+min_land_point_proc-1)

! all_data(v,counter,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr) all_data(v,counter-chunk_start+1,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr)

i'm assuming that either all_data or data are being allocated or accessed incorectly, but i really can't see why 24 days would be the magic number that get's screwed up. one thing that's interesting is that day_start = 1, meaning that it looks like the last line reduces to

    all_data(v,counter-chunk_start+1,1:1+nday_yr-1) =

data(jlon,ilat,1:nday_yr)

and that looks to me like the time dimension is off by one (the LHS goes from 1 to 1+nday-1 and the RHS to 1+nday). perhaps day_start is supposed to be zero? this has always been like this though, so i'm probably missing something.

hopefully we haven't completely botched up the code. let me know what you think.

joshua

On Fri, Dec 28, 2012 at 12:27 PM, joshua elliott jelliott@ci.uchicago.eduwrote:

hey david. sorry for the barrage of emails, but i wanted to add one more piece of info. i realized that there was a typo in the new grid_ID formula regarding the order of operations in the argument for floor(). it should be (90-lat)_12 rather than 90-lat_12, which means that the argument should never actually run negative anyway. since the same is true for the argument of ceiling(), we don't have to worry about that weird behaviour (though i'm still curious why its happening).

On Fri, Dec 28, 2012 at 10:29 AM, joshua elliott <jelliott@ci.uchicago.edu

wrote:

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 **

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 **

+! when longitiude runs from 0 to 360:

279 **

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

  280 **

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 **

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.comwrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit < https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

[joshuaelliott]

i also wanted to add some of the stdout in case it illuminates anything: nbest@midway-login2 narr]$ nc_wth_gen.sh 1979 1981 data/nc data/wth2 GENERIC1.WTH 100 50 start yr: 1979 end yr: 1981 fileroot:data/nc outdir:data/wth2 fname:GENERIC1.WTH n_procs: 100 proc_num: 50

data/nc/precip/precip_1979.nc precip data/nc/solar/solar_1979.nc solar data/nc/tmax/tmax_1979.nc tmax data/nc/tmin/tmin_1979.nc tmin calculating land points n_land_points_all 460800 min_land_point_proc 225793 max_land_point_proc 230400 n_land_points_proc 4608 reading data chunk start= 1 chunk end= 4608

On Fri, Dec 28, 2012 at 6:08 PM, joshua elliott jelliott@ci.uchicago.eduwrote:

yet another update on the status of the f90 debugging:

currently, the issue we're trying to figure out has to do with garbage data that appears in the first 24 days of the .WTH files (see the attached example). i can't figure out what is happening but i figure it must have something to do with an array being defined as an inappropriate size or accessed inappropriately. i thought this might have something to do with max_points_all, since there was a compiler error associated with this early on and neil had to bump up the value somewhat to avoid it, but i think it probably has more to do with how min_land_point_proc is getting set, and how that effects ilat and jlon.

! copy data from netcdf file to big array of data ! do counter = 1, n_land_points_proc

  do counter = chunk_start,chunk_end

    ilat = print_point_lat(counter+min_land_point_proc-1)

    jlon = print_point_lon(counter+min_land_point_proc-1)

! all_data(v,counter,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr)

    all_data(v,counter-chunk_start+1,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr)

i'm assuming that either all_data or data are being allocated or accessed incorectly, but i really can't see why 24 days would be the magic number that get's screwed up. one thing that's interesting is that day_start = 1, meaning that it looks like the last line reduces to

    all_data(v,counter-chunk_start+1,1:1+nday_yr-1) = data(jlon,ilat,1:nday_yr)

and that looks to me like the time dimension is off by one (the LHS goes from 1 to 1+nday-1 and the RHS to 1+nday). perhaps day_start is supposed to be zero? this has always been like this though, so i'm probably missing something.

hopefully we haven't completely botched up the code. let me know what you think.

joshua

On Fri, Dec 28, 2012 at 12:27 PM, joshua elliott <jelliott@ci.uchicago.edu

wrote:

hey david. sorry for the barrage of emails, but i wanted to add one more piece of info. i realized that there was a typo in the new grid_ID formula regarding the order of operations in the argument for floor(). it should be (90-lat)_12 rather than 90-lat_12, which means that the argument should never actually run negative anyway. since the same is true for the argument of ceiling(), we don't have to worry about that weird behaviour (though i'm still curious why its happening).

On Fri, Dec 28, 2012 at 10:29 AM, joshua elliott < jelliott@ci.uchicago.edu> wrote:

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 **

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 **

+! when longitiude runs from 0 to 360:

279 **

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

  280 **

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 **

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.comwrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit < https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

[dmcinern]

Hi Josh,

Sorry for not replying sooner. I've been away for a couple of days. From your emails it looks like you've sorted out the grid_id issues (let me know if i'm wrong) and you're trying to figure out what's going on with the first 24 days. In a previous email you sent you mentioned the first 8 days were causing problems, so i'm trying to work out why either the first 8 or 24 days would not be printing correctly, and what would cause this number to change.

My guess is that the problem is a time indexing issue, rather than anything to do with the number or indexing of land points. In a previous email I asked:

Does this problem occur at all grid points? And does the problem start occurring after line 249 in the code? i.e. are the arrays solar, tmax, tmin and precip OK? If it occurs before this point we have a problem with the program reading in the data. if it occurs later then it's an issue with processing the data and writing it out.

I couldn't find a response to this. The way to check this would be to add this code to line 263.

!v=1 (precip), position = 1 print*, all_data(1,1,1:nday_yr) stop

and then run the program from the command line

nc_wth_gen.sh 1979 1981 data/nc data/wth2 GENERIC1.WTH 100 1

This should give you the first year of precip data (directly from the netcdf file) for the first land point. If you see problems with the first 24 days then there is an issue reading in the data from the netcdf file. If it looks ok then the problem occurs later in the code, in which case you could check on line 333 to see if the precip is ok

print*,precip stop

If that looks ok then it's a problem writing the data.

I hope this helps. Btw, I just applied for an account on midway, so hopefully that will be available soon and i can help out more.

David

On Dec 29, 2012, at 9:39 AM, Joshua Elliott wrote:

yet another update on the status of the f90 debugging:

currently, the issue we're trying to figure out has to do with garbage data that appears in the first 24 days of the .WTH files (see the attached example). i can't figure out what is happening but i figure it must have something to do with an array being defined as an inappropriate size or accessed inappropriately. i thought this might have something to do with max_points_all, since there was a compiler error associated with this early on and neil had to bump up the value somewhat to avoid it, but i think it probably has more to do with how min_land_point_proc is getting set, and how that effects ilat and jlon.

max_points_all is the total number of land points. this is going to be <= nlat_all * nlon_all, so setting it to nlat_all*nlon_all is fine (but maybe this could be reduced to make things faster once we sort out the current problem).

! copy data from netcdf file to big array of data ! do counter = 1, n_land_points_proc

do counter = chunk_start,chunk_end ilat = print_point_lat(counter+min_land_point_proc-1) jlon = print_point_lon(counter+min_land_point_proc-1) ! all_data(v,counter,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr) all_data(v,counter-chunk_start+1,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr)

i'm assuming that either all_data or data are being allocated or accessed incorectly, but i really can't see why 24 days would be the magic number that get's screwed up. one thing that's interesting is that day_start = 1, meaning that it looks like the last line reduces to

all_data(v,counter-chunk_start+1,1:1+nday_yr-1) = data(jlon,ilat,1:nday_yr)

and that looks to me like the time dimension is off by one (the LHS goes from 1 to 1+nday-1 and the RHS to 1+nday). perhaps day_start is supposed to be zero? this has always been like this though, so i'm probably missing something.

the LHS goes from 1:nday_yr and so does the RHS, right?

hopefully we haven't completely botched up the code. let me know what you think.

joshua

On Fri, Dec 28, 2012 at 12:27 PM, joshua elliott jelliott@ci.uchicago.eduwrote:

hey david. sorry for the barrage of emails, but i wanted to add one more piece of info. i realized that there was a typo in the new grid_ID formula regarding the order of operations in the argument for floor(). it should be (90-lat)_12 rather than 90-lat_12, which means that the argument should never actually run negative anyway. since the same is true for the argument of ceiling(), we don't have to worry about that weird behaviour (though i'm still curious why its happening).

On Fri, Dec 28, 2012 at 10:29 AM, joshua elliott <jelliott@ci.uchicago.edu

wrote:

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 \

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 \

+! when longitiude runs from 0 to 360:

279 \

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

280 \

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 \

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.comwrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit < https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630 — Reply to this email directly or view it on GitHub.

[joshuaelliott]

Awesome david, thanks so much. please don't worry about reply speed. I know its the holidays and I just figured it was better I keep a record of what I was doing on the email thread like a support ticket. Given my impotence with Fortran, all I could really do is bounce around with the logic and hope something funny shook loose. I will try yr suggested fixes as soon as possible. Hopefully during this long wknd though thats definitely wishful thinking.

Would definitely be great for you to have midway access, but I can probably muddle through with regular desperate support requests if you prefer, so either way is fine. If you have half an hour next week where we could Skype audio and you can walk me through the logic let me know. Might save us both a bit of time.

Thanks so much! Joshua

--- Sent from mobile --- Joshua W. Elliott Research Scientist and Fellow Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Tel: 773-834-6812; Fax: 773-834-6818 On Dec 28, 2012 10:48 PM, "dmcinern" notifications@github.com wrote:

Hi Josh,

Sorry for not replying sooner. I've been away for a couple of days. From your emails it looks like you've sorted out the grid_id issues (let me know if i'm wrong) and you're trying to figure out what's going on with the first 24 days. In a previous email you sent you mentioned the first 8 days were causing problems, so i'm trying to work out why either the first 8 or 24 days would not be printing correctly, and what would cause this number to change.

My guess is that the problem is a time indexing issue, rather than anything to do with the number or indexing of land points. In a previous email I asked:

Does this problem occur at all grid points? And does the problem start occurring after line 249 in the code? i.e. are the arrays solar, tmax, tmin and precip OK? If it occurs before this point we have a problem with the program reading in the data. if it occurs later then it's an issue with processing the data and writing it out.

I couldn't find a response to this. The way to check this would be to add this code to line 263.

!v=1 (precip), position = 1 print*, all_data(1,1,1:nday_yr) stop

and then run the program from the command line

nc_wth_gen.sh 1979 1981 data/nc data/wth2 GENERIC1.WTH 100 1

This should give you the first year of precip data (directly from the netcdf file) for the first land point. If you see problems with the first 24 days then there is an issue reading in the data from the netcdf file. If it looks ok then the problem occurs later in the code, in which case you could check on line 333 to see if the precip is ok

print*,precip stop

If that looks ok then it's a problem writing the data.

I hope this helps. Btw, I just applied for an account on midway, so hopefully that will be available soon and i can help out more.

David

On Dec 29, 2012, at 9:39 AM, Joshua Elliott wrote:

yet another update on the status of the f90 debugging:

currently, the issue we're trying to figure out has to do with garbage data that appears in the first 24 days of the .WTH files (see the attached example). i can't figure out what is happening but i figure it must have something to do with an array being defined as an inappropriate size or accessed inappropriately. i thought this might have something to do with max_points_all, since there was a compiler error associated with this early on and neil had to bump up the value somewhat to avoid it, but i think it probably has more to do with how min_land_point_proc is getting set, and how that effects ilat and jlon.

max_points_all is the total number of land points. this is going to be <= nlat_all * nlon_all, so setting it to nlat_all*nlon_all is fine (but maybe this could be reduced to make things faster once we sort out the current problem).

! copy data from netcdf file to big array of data ! do counter = 1, n_land_points_proc

do counter = chunk_start,chunk_end ilat = print_point_lat(counter+min_land_point_proc-1) jlon = print_point_lon(counter+min_land_point_proc-1) ! all_data(v,counter,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr) all_data(v,counter-chunk_start+1,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr)

i'm assuming that either all_data or data are being allocated or accessed incorectly, but i really can't see why 24 days would be the magic number that get's screwed up. one thing that's interesting is that day_start = 1, meaning that it looks like the last line reduces to

all_data(v,counter-chunk_start+1,1:1+nday_yr-1) = data(jlon,ilat,1:nday_yr)

and that looks to me like the time dimension is off by one (the LHS goes from 1 to 1+nday-1 and the RHS to 1+nday). perhaps day_start is supposed to be zero? this has always been like this though, so i'm probably missing something.

the LHS goes from 1:nday_yr and so does the RHS, right?

hopefully we haven't completely botched up the code. let me know what you think.

joshua

On Fri, Dec 28, 2012 at 12:27 PM, joshua elliott jelliott@ci.uchicago.eduwrote:

hey david. sorry for the barrage of emails, but i wanted to add one more piece of info. i realized that there was a typo in the new grid_ID formula regarding the order of operations in the argument for floor(). it should be (90-lat)_12 rather than 90-lat_12, which means that the argument should never actually run negative anyway. since the same is true for the argument of ceiling(), we don't have to worry about that weird behaviour (though i'm still curious why its happening).

On Fri, Dec 28, 2012 at 10:29 AM, joshua elliott < jelliott@ci.uchicago.edu

wrote:

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 **

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 **

+! when longitiude runs from 0 to 360:

279 **

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

280 **

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 **

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.comwrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit <

https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHub< https://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269>.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630 — Reply to this email directly or view it on GitHub.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11747087.

[joshuaelliott]

see below.

Does this problem occur at all grid points? And does the problem start

occurring after line 249 in the code? i.e. are the arrays solar, tmax, tmin and precip OK? If it occurs before this point we have a problem with the program reading in the data. if it occurs later then it's an issue with processing the data and writing it out.

I couldn't find a response to this. The way to check this would be to add this code to line 263.

!v=1 (precip), position = 1 print*, all_data(1,1,1:nday_yr) stop

i do believe it occurs on all grid points, yes. there is definitely problems with the first 24 days here, so i guess the problem is with reading the data:

[joshuaelliott@midway-login1 wth_gen]$ ./nc_wth_gen 1979 1981 /project/joshuaelliott/narr/data/nc /project/joshuaelliott/users/jelliott/wth GENERIC1.WTH 100 1 start yr: 1979 end yr: 1981 fileroot:/project/joshuaelliott/narr/data/nc outdir:/project/joshuaelliott/users/jelliott/wth fname:GENERIC1.WTH n_procs: 100 proc_num: 1 ..... solar /project/joshuaelliott/narr/data/nc/tmax/tmax_1981.nc

tmax /project/joshuaelliott/narr/data/nc/tmin/tmin_1981.nc

tmin 428.87500 553.37500 354.50000 137.56250 19.062500 216.93750 454.87500 266.68750 88.375000 8.0000000 283.06250 375.25000 177.62500 28.812500 1.4375000 94.750000 206.93750 78.437500 12.812500 0.0000000 157.75000 100.81250 25.062500 0.12500000 0.0000000 0.43273774 3.8882337 1.2508581 0.26484609 8.29597712E-02 1.5200049 0.76495135 0.90222436 9.91761684E-03 8.83305073E-03 1.2190942 0.80468750 2.2871804 4.8631272 0.54709667 1.3758547 1.2366278 3.1601200 10.810010 4.99007702E-02 5.24642318E-03 2.4730859 1.6830666 0.70921081 0.33593750 3.06755560E-03 3.4556451 0.57812500 7.81250000E-03 0.22048584 0.74002111 1.56250000E-02 1.76470280E-02 6.42454624E-03 3.8135412 3.8281250 0.11775509 0.47011262 3.8993871 1.3174862 0.20050412 1.70095749E-02 ..............

On Dec 29, 2012, at 9:39 AM, Joshua Elliott wrote:

yet another update on the status of the f90 debugging:

currently, the issue we're trying to figure out has to do with garbage data that appears in the first 24 days of the .WTH files (see the attached example). i can't figure out what is happening but i figure it must have something to do with an array being defined as an inappropriate size or accessed inappropriately. i thought this might have something to do with max_points_all, since there was a compiler error associated with this early on and neil had to bump up the value somewhat to avoid it, but i think it probably has more to do with how min_land_point_proc is getting set, and how that effects ilat and jlon.

max_points_all is the total number of land points. this is going to be <= nlat_all * nlon_all, so setting it to nlat_all*nlon_all is fine (but maybe this could be reduced to make things faster once we sort out the current problem).

! copy data from netcdf file to big array of data ! do counter = 1, n_land_points_proc

do counter = chunk_start,chunk_end ilat = print_point_lat(counter+min_land_point_proc-1) jlon = print_point_lon(counter+min_land_point_proc-1) ! all_data(v,counter,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr) all_data(v,counter-chunk_start+1,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr)

i'm assuming that either all_data or data are being allocated or accessed incorectly, but i really can't see why 24 days would be the magic number that get's screwed up. one thing that's interesting is that day_start = 1, meaning that it looks like the last line reduces to

all_data(v,counter-chunk_start+1,1:1+nday_yr-1) = data(jlon,ilat,1:nday_yr)

and that looks to me like the time dimension is off by one (the LHS goes from 1 to 1+nday-1 and the RHS to 1+nday). perhaps day_start is supposed to be zero? this has always been like this though, so i'm probably missing something.

the LHS goes from 1:nday_yr and so does the RHS, right?

hopefully we haven't completely botched up the code. let me know what you think.

joshua

On Fri, Dec 28, 2012 at 12:27 PM, joshua elliott jelliott@ci.uchicago.eduwrote:

hey david. sorry for the barrage of emails, but i wanted to add one more piece of info. i realized that there was a typo in the new grid_ID formula regarding the order of operations in the argument for floor(). it should be (90-lat)_12 rather than 90-lat_12, which means that the argument should never actually run negative anyway. since the same is true for the argument of ceiling(), we don't have to worry about that weird behaviour (though i'm still curious why its happening).

On Fri, Dec 28, 2012 at 10:29 AM, joshua elliott < jelliott@ci.uchicago.edu

wrote:

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 **

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 **

+! when longitiude runs from 0 to 360:

279 **

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

280 **

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 **

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.comwrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit <

https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHub< https://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269>.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630 — Reply to this email directly or view it on GitHub.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11747087.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

[joshuaelliott]

one thing i just discovered.

apparently the nc input files were printed with 1 more day then you would expect, so that regular years have 366 rather than 365 days. i'm not sure if its an extra day at the beginning or end, and i really don't think it matters, but it seems like its possible that this is screwing up the read data step, so i figured i'd point it out. one obvious potential culprit is that the data array is allocated with only nday steps in the time dimensions, but the read data function pulls in an nc file that is bigger. i will push up the nday_yr intergers and see what happens.

On Sat, Dec 29, 2012 at 9:04 AM, joshua elliott jelliott@ci.uchicago.eduwrote:

see below.

Does this problem occur at all grid points? And does the problem start

occurring after line 249 in the code? i.e. are the arrays solar, tmax, tmin and precip OK? If it occurs before this point we have a problem with the program reading in the data. if it occurs later then it's an issue with processing the data and writing it out.

I couldn't find a response to this. The way to check this would be to add this code to line 263.

!v=1 (precip), position = 1 print*, all_data(1,1,1:nday_yr) stop

i do believe it occurs on all grid points, yes. there is definitely problems with the first 24 days here, so i guess the problem is with reading the data:

[joshuaelliott@midway-login1 wth_gen]$ ./nc_wth_gen 1979 1981 /project/joshuaelliott/narr/data/nc /project/joshuaelliott/users/jelliott/wth GENERIC1.WTH 100 1 start yr: 1979 end yr: 1981 fileroot:/project/joshuaelliott/narr/data/nc outdir:/project/joshuaelliott/users/jelliott/wth fname:GENERIC1.WTH n_procs: 100 proc_num: 1 ..... solar /project/joshuaelliott/narr/data/nc/tmax/tmax_1981.nc

tmax /project/joshuaelliott/narr/data/nc/tmin/tmin_1981.nc

tmin 428.87500 553.37500 354.50000 137.56250 19.062500 216.93750 454.87500 266.68750 88.375000 8.0000000 283.06250 375.25000 177.62500 28.812500 1.4375000 94.750000 206.93750 78.437500 12.812500 0.0000000 157.75000 100.81250 25.062500 0.12500000 0.0000000 0.43273774 3.8882337 1.2508581 0.26484609 8.29597712E-02 1.5200049 0.76495135 0.90222436 9.91761684E-03 8.83305073E-03 1.2190942 0.80468750 2.2871804 4.8631272 0.54709667 1.3758547 1.2366278 3.1601200 10.810010 4.99007702E-02 5.24642318E-03 2.4730859 1.6830666 0.70921081 0.33593750 3.06755560E-03 3.4556451 0.57812500 7.81250000E-03 0.22048584 0.74002111 1.56250000E-02 1.76470280E-02 6.42454624E-03 3.8135412 3.8281250 0.11775509 0.47011262 3.8993871 1.3174862 0.20050412 1.70095749E-02 ..............

On Dec 29, 2012, at 9:39 AM, Joshua Elliott wrote:

yet another update on the status of the f90 debugging:

currently, the issue we're trying to figure out has to do with garbage data that appears in the first 24 days of the .WTH files (see the attached example). i can't figure out what is happening but i figure it must have something to do with an array being defined as an inappropriate size or accessed inappropriately. i thought this might have something to do with max_points_all, since there was a compiler error associated with this early on and neil had to bump up the value somewhat to avoid it, but i think it probably has more to do with how min_land_point_proc is getting set, and how that effects ilat and jlon.

max_points_all is the total number of land points. this is going to be <= nlat_all * nlon_all, so setting it to nlat_all*nlon_all is fine (but maybe this could be reduced to make things faster once we sort out the current problem).

! copy data from netcdf file to big array of data ! do counter = 1, n_land_points_proc

do counter = chunk_start,chunk_end ilat = print_point_lat(counter+min_land_point_proc-1) jlon = print_point_lon(counter+min_land_point_proc-1) ! all_data(v,counter,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr) all_data(v,counter-chunk_start+1,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr)

i'm assuming that either all_data or data are being allocated or accessed incorectly, but i really can't see why 24 days would be the magic number that get's screwed up. one thing that's interesting is that day_start = 1, meaning that it looks like the last line reduces to

all_data(v,counter-chunk_start+1,1:1+nday_yr-1) = data(jlon,ilat,1:nday_yr)

and that looks to me like the time dimension is off by one (the LHS goes from 1 to 1+nday-1 and the RHS to 1+nday). perhaps day_start is supposed to be zero? this has always been like this though, so i'm probably missing something.

the LHS goes from 1:nday_yr and so does the RHS, right?

hopefully we haven't completely botched up the code. let me know what you think.

joshua

On Fri, Dec 28, 2012 at 12:27 PM, joshua elliott jelliott@ci.uchicago.eduwrote:

hey david. sorry for the barrage of emails, but i wanted to add one more piece of info. i realized that there was a typo in the new grid_ID formula regarding the order of operations in the argument for floor(). it should be (90-lat)_12 rather than 90-lat_12, which means that the argument should never actually run negative anyway. since the same is true for the argument of ceiling(), we don't have to worry about that weird behaviour (though i'm still curious why its happening).

On Fri, Dec 28, 2012 at 10:29 AM, joshua elliott < jelliott@ci.uchicago.edu

wrote:

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 **

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 **

+! when longitiude runs from 0 to 360:

279 **

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

280 **

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 **

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.comwrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit <

https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHub< https://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269>.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630 — Reply to this email directly or view it on GitHub.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11747087.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

[joshuaelliott]

ok, looks like that's the problem. if i just do

! allocate(data(nlon_all,nlat_all,nday)) allocate(data(nlon_all,nlat_all,nday+1))

then it seems to work ok. long term its probably best that we take the fake day off the input files, but for now this is hopefully doing the correct thing. let me know if you think there are any obvious shortfalls to this hack.

thanks! joshua

On Sat, Dec 29, 2012 at 9:56 AM, joshua elliott jelliott@ci.uchicago.eduwrote:

one thing i just discovered.

apparently the nc input files were printed with 1 more day then you would expect, so that regular years have 366 rather than 365 days. i'm not sure if its an extra day at the beginning or end, and i really don't think it matters, but it seems like its possible that this is screwing up the read data step, so i figured i'd point it out. one obvious potential culprit is that the data array is allocated with only nday steps in the time dimensions, but the read data function pulls in an nc file that is bigger. i will push up the nday_yr intergers and see what happens.

On Sat, Dec 29, 2012 at 9:04 AM, joshua elliott jelliott@ci.uchicago.eduwrote:

see below.

Does this problem occur at all grid points? And does the problem start

occurring after line 249 in the code? i.e. are the arrays solar, tmax, tmin and precip OK? If it occurs before this point we have a problem with the program reading in the data. if it occurs later then it's an issue with processing the data and writing it out.

I couldn't find a response to this. The way to check this would be to add this code to line 263.

!v=1 (precip), position = 1 print*, all_data(1,1,1:nday_yr) stop

i do believe it occurs on all grid points, yes. there is definitely problems with the first 24 days here, so i guess the problem is with reading the data:

[joshuaelliott@midway-login1 wth_gen]$ ./nc_wth_gen 1979 1981 /project/joshuaelliott/narr/data/nc /project/joshuaelliott/users/jelliott/wth GENERIC1.WTH 100 1 start yr: 1979 end yr: 1981 fileroot:/project/joshuaelliott/narr/data/nc outdir:/project/joshuaelliott/users/jelliott/wth fname:GENERIC1.WTH n_procs: 100 proc_num: 1 ..... solar /project/joshuaelliott/narr/data/nc/tmax/tmax_1981.nc

tmax /project/joshuaelliott/narr/data/nc/tmin/tmin_1981.nc

tmin 428.87500 553.37500 354.50000 137.56250 19.062500 216.93750 454.87500 266.68750 88.375000 8.0000000 283.06250 375.25000 177.62500 28.812500 1.4375000 94.750000 206.93750 78.437500 12.812500 0.0000000 157.75000 100.81250 25.062500 0.12500000 0.0000000 0.43273774 3.8882337 1.2508581 0.26484609 8.29597712E-02 1.5200049 0.76495135 0.90222436 9.91761684E-03 8.83305073E-03 1.2190942 0.80468750 2.2871804 4.8631272 0.54709667 1.3758547 1.2366278 3.1601200 10.810010 4.99007702E-02 5.24642318E-03 2.4730859 1.6830666 0.70921081 0.33593750 3.06755560E-03 3.4556451 0.57812500 7.81250000E-03 0.22048584 0.74002111 1.56250000E-02 1.76470280E-02 6.42454624E-03 3.8135412 3.8281250 0.11775509 0.47011262 3.8993871 1.3174862 0.20050412 1.70095749E-02 ..............

On Dec 29, 2012, at 9:39 AM, Joshua Elliott wrote:

yet another update on the status of the f90 debugging:

currently, the issue we're trying to figure out has to do with garbage data that appears in the first 24 days of the .WTH files (see the attached example). i can't figure out what is happening but i figure it must have something to do with an array being defined as an inappropriate size or accessed inappropriately. i thought this might have something to do with max_points_all, since there was a compiler error associated with this early on and neil had to bump up the value somewhat to avoid it, but i think it probably has more to do with how min_land_point_proc is getting set, and how that effects ilat and jlon.

max_points_all is the total number of land points. this is going to be <= nlat_all * nlon_all, so setting it to nlat_all*nlon_all is fine (but maybe this could be reduced to make things faster once we sort out the current problem).

! copy data from netcdf file to big array of data ! do counter = 1, n_land_points_proc

do counter = chunk_start,chunk_end ilat = print_point_lat(counter+min_land_point_proc-1) jlon = print_point_lon(counter+min_land_point_proc-1) ! all_data(v,counter,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr) all_data(v,counter-chunk_start+1,day_start:day_start+nday_yr-1) = data(jlon,ilat,1:nday_yr)

i'm assuming that either all_data or data are being allocated or accessed incorectly, but i really can't see why 24 days would be the magic number that get's screwed up. one thing that's interesting is that day_start = 1, meaning that it looks like the last line reduces to

all_data(v,counter-chunk_start+1,1:1+nday_yr-1) = data(jlon,ilat,1:nday_yr)

and that looks to me like the time dimension is off by one (the LHS goes from 1 to 1+nday-1 and the RHS to 1+nday). perhaps day_start is supposed to be zero? this has always been like this though, so i'm probably missing something.

the LHS goes from 1:nday_yr and so does the RHS, right?

hopefully we haven't completely botched up the code. let me know what you think.

joshua

On Fri, Dec 28, 2012 at 12:27 PM, joshua elliott jelliott@ci.uchicago.eduwrote:

hey david. sorry for the barrage of emails, but i wanted to add one more piece of info. i realized that there was a typo in the new grid_ID formula regarding the order of operations in the argument for floor(). it should be (90-lat)_12 rather than 90-lat_12, which means that the argument should never actually run negative anyway. since the same is true for the argument of ceiling(), we don't have to worry about that weird behaviour (though i'm still curious why its happening).

On Fri, Dec 28, 2012 at 10:29 AM, joshua elliott < jelliott@ci.uchicago.edu

wrote:

hey david. couple of quick questions. i'm sorry that we are so fortran illiterate. i'm trying, but right at the moment i'm absolutely in a desperate rush to get this working so i'm reaching out and begging for help rather than learning new things. please forgive me.

neil had to adjust the grid_ID calculation because the longitudes in the NARR data run from 0 to 360 rather than -180 to 180. he printed a bunch of data before he did this and the code was running fine, but the grid_ID folder names were wrong of course.

he modified the code very simply as follows, simply shifting the longitudes above 180 to their appropriate values:

277 **

+! grid_ID = nint(( 90. - lat(ilat)) 51840 +( 180 +lon(jlon)) 12)

278 **

+! when longitiude runs from 0 to 360:

279 **

• if( lon( jlon) > 180.) lon( jlon) = lon( jlon) - 360.

280 **

+! when longitude runs from -180 to 180 simply comment out the above statement.

281 **

• grid_ID = floor( 90. -lat( ilat) 12.) 360 12 +ceiling(( 180. +lon( jlon)) 12.)

my first question is whether the syntax on the added if statement is correct. does it need a then or anything? is it rewriting the value in the lon array appropriately?

secondly, it appears than in the latest working version of the code the grid_ID was actually being calculated twice. once above, and then again in the calc_land_points subroutine. neil has modified this function substantially, with the goal i believe being to use no mask at all but instead for the function to simply include all grid cells. can you look at this version: https://github.com/RDCEP/wth_gen/blob/narr/nc_wth_gen.f90and make sure that everything with that function appears copacetic?

my assumption is that the above global version of grid_ID is the one that gets used to determine the folder structure for printing the files. is that correct?

thanks so much for everything you've done on this. cheers, joshua

On Tue, Dec 18, 2012 at 5:18 PM, dmcinern notifications@github.comwrote:

Hi Neil,

I'm not sure why the floor(-150.5) is returning 0. i tried using this test program,

!!!!!!!!!!! program test implicit none print*,floor(-150.5) !!!!!!!!!!!

with both gfortran and ifort, and got -151.

David

On 19/12/2012 6:51 AM, Neil Best wrote:

Joshua reported that grid IDs that nc_wth_gen writes are not valid so I tried to improve the formula that calculates |grid_ID| in this commit <

https://github.com/RDCEP/wth_gen/commit/d9e2ef80062ba85f06a0dfdd17aa7d18b663fb5d#nc_wth_gen.f90>

but I must not be using the |floor()| function correctly because when I go into the debugger I see that the first term in the |grid_ID| expression returns zero:

(gdb) print grid_ID $7 = -651839 (gdb) print 90. -lat( ilat) 12. $8 = -150.50001525878906 (gdb) print floor( 90. -lat( ilat) 12.) $9 = 0 (gdb) print nint( 90. -lat( ilat) 12.) No symbol "nint" in current context. (gdb) print floor( -150.5) $10 = 0 (gdb) print ceiling(( 180. +lon( jlon)) 12.) No symbol "ceiling" in current context.

As a result the data is all written to a file called literally |_/_****/GENERIC1.WTH| because of how the negative grid_ID gets formatted.

Any idea what is going on here? Is this a type mismatch problem? Any guess why I couldn't try |nint()| or |ceiling()|? I must need some other inlcude file or compiler directive to get these functions to work properly. From my research I suspect that this has something to do with the various flavors/standards of the language.

I am learning Fortran and GDB on the fly here so assume nothing about my proficiency. However, it is clear that the logic in this program is very dependent on the particulars of the input data.

— Reply to this email directly or view it on GitHub https://github.com/RDCEP/wth_gen/issues/5.

— Reply to this email directly or view it on GitHub< https://github.com/RDCEP/wth_gen/issues/5#issuecomment-11510269>.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal website< https://sites.google.com/site/joshuawrightelliott/> ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630 — Reply to this email directly or view it on GitHub.

— Reply to this email directly or view it on GitHubhttps://github.com/RDCEP/wth_gen/issues/5#issuecomment-11747087.

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630

Joshua W. Elliott Research Scientist and Fellow University of Chicago and ANL Computation Institute 5735 S. Ellis Ave. Chicago, IL 60637 Links: Personal websitehttps://sites.google.com/site/joshuawrightelliott/ ; SSRN papers http://ssrn.com/author=1655092; E-mail: jelliott@ci.uchicago.edu

Adjunct Research Scientist Columbia University Center for Climate Systems Research NASA GISS; 2880 Broadway. NY, NY 10025; Rm 341 Tel: 212-678-5630