MESAHub / mesa

Modules for Experiments in Stellar Astrophysics
http://mesastar.org
GNU Lesser General Public License v2.1
145 stars 39 forks source link

Color files break when using [Fe/H] instead of M_div_H #379

Closed rjfarmer closed 2 years ago

rjfarmer commented 2 years ago

Well this was a fun bug to find. This is with 21.12.1

If you have a custom colors file where the header is:

#Teff logg M_div_H bessell_u Bessell_B Bessell_V Bessell_R Bessell_I

Everything works fine

But if you instead have:

#Teff logg [Fe/H] bessell_u Bessell_B Bessell_V Bessell_R Bessell_I

(Which is the mist way of writing the metallicity) then instead of correctly reading in the filter names (bessell_u Bessell_B Bessell_V Bessell_R Bessell_I) we end up with the partial contents of data/chem_data/lodders03.data file as the filter names.

If you add

write(*,*) trim(fname),tmp,tmp,tmp,col_names(1:n_colors)

at line 249 in colors/private/mod_colors.f90 then in the working case we have:

colors.dataM_div_H                                                                                                                                                                                                                                                         M_div_H                                                                                                                                                                                                                                                         M_div_H                                                                                                                                                                                                                                                         bessell_u                                                                                                                                                                                                                                                       Bessell_B                                                                                                                                                                                                                                                       Bessell_V                                                                                                                                                                                                                                                       Bessell_R                                                                                                                                                                                                                                                       Bessell_I 

and in the broken case we have:

colors.data[Fe                                                                                                                                                                                                                                                             [Fe                                                                                                                                                                                                                                                             [Fe                                                                                                                                                                                                                                                             �1���
                                                                                                                                                     �7�
                                                                                                                                                        �7                                                                                 
 60 Nd   150      5.62        0.04695           
 62 Sm   144      3.0734      0.00781   
 62 Sm   147     14.9934      0.03811   
 62 Sm   148     11.2406      0.0286    
 62 Sm   149     13.8189      0.0351    
 62 Sm   150      7.3796      0.0188    
 62 Sm   152     26.7421      0.0680    
 62 Sm   154     22.752       0.0578    
 63 Eu   151     47.81        0.04548   
 63 Eu   153     52.19        0.04965   
 64 Gd   152      0.2029      0.00067   
 64 Gd   154      2.1809      0.00724   
 64 Gd   155     14.7998      0.04915   
 64 Gd   156     20.4664      0.06797   
 64 Gd   157     15.6518      0.05198   
 64 Gd   158     24.8347      0.08248   
 64 Gd   160     21.8635      0.07261   
 65 Tb   159    100           0.05907   
 66 Dy   156      0.056       0.000216  
 66 Dy   158      0.096       0.000371  
 66 Dy   160      2.34        0.00904   
 66 Dy   161     18.91        0.0730    
 66 Dy   162     25.51        0.0985    
 66 Dy   163     24.9         0.0962    
 66 Dy   164     28.19        0.1089    
 67 Ho   165    100           0.08986   
 68 Er   162      0.137       0.000350  
 68 Er   164      1.609       0.004109  
 68 Er   166     3

Some how [Fe/H] is breaking gfortran's read routines and i guess we are getting a memory corruption and then start reading a different file in?

Maybe [Fe is being read as a ANSI escape code perhaps? Then /h] means something?

https://en.wikipedia.org/wiki/ANSI_escape_code#Fe_Escape_sequences

Here is a minimal breaking test case file:

#Teff logg [Fe/H] bessell_u Bessell_B Bessell_V Bessell_R Bessell_I
2500 -4 -0.25 -8.003592 -8.111115 -6.819048 -4.321354 -1.940325
2500 -3 -0.25 -8.003592 -8.111115 -6.819048 -4.321354 -1.940325
2500 -2 -0.25 -8.003592 -8.111115 -6.819048 -4.321354 -1.940325
2500 -1.5 -0.25 -8.003592 -8.111115 -6.819048 -4.321354 -1.940325
2500 -1 -0.25 -8.003592 -8.111115 -6.819048 -4.321354 -1.940325
2500 -0.5 -0.25 -8.003592 -8.111115 -6.819048 -4.321354 -1.940325
2500 0 -0.25 -8.003592 -8.111115 -6.819048 -4.321354 -1.940325
2500 0.5 -0.25 -6.360487 -6.760914 -6.465778 -3.81477 -1.366576
2500 1 -0.25 -5.950181 -6.418117 -6.413375 -3.694483 -1.235039
2500 1.5 -0.25 -5.700899 -6.176356 -6.30592 -3.569254 -1.111708
2500 2 -0.25 -5.710777 -6.036413 -5.909205 -3.319799 -0.89653
2500 2.5 -0.25 -5.961649 -6.053194 -5.548792 -3.125635 -0.743097
2500 3 -0.25 -6.358657 -6.20274 -5.434846 -3.061187 -0.692906
2500 3.5 -0.25 -6.358657 -6.20274 -5.434846 -3.061187 -0.692906
2500 4 -0.25 -7.614005 -6.817242 -5.369132 -2.995211 -0.689922
2500 4.5 -0.25 -8.517534 -7.308242 -5.518233 -3.064574 -0.732102
2500 5 -0.25 -9.393778 -7.735753 -5.65708 -3.023487 -0.604437
2500 5.5 -0.25 -10.162364 -7.838414 -5.498994 -2.775208 -0.276668
2500 6 -0.25 -10.162364 -7.838414 -5.498994 -2.775208 -0.276668
fxt44 commented 2 years ago

what is the read that parses and stores the #Teff line?

rjfarmer commented 2 years ago
character(len=*),dimension(:),pointer,intent(inout) :: col_names
character(len=strlen) :: tmp

read(IO_UBV,fmt=*,iostat=ios) tmp,tmp,tmp,col_names(1:n_colors)

We don't even care what someone calls metallicity we only check for three strings at the start of the line.

fxt44 commented 2 years ago

for debugging maybe try

character(len=strlen) :: tmp,tmp2,tmp3 read(IO_UBV,fmt=*,iostat=ios) tmp,tmp2,tmp3,col_names(1:n_colors)

and a formatted character read

read(IO_UBV,fmt='(a)',iostat=ios) tmp,tmp,tmp,col_names(1:n_colors)

rjfarmer commented 2 years ago

Switching to multiple tmps doesn't help while using a formatted read does at least keeps us reading the same file but we end up reading entire rows into each variable instead of just the column.

Doing something like a fmt='(8(a10))' and configuring the header to have the right string lengths does work (at the expense that now everything must be the same width).

So it seems its unformatted reads of [Fe/H] is causing gfortran to read a different file.

For now we can try can catch the square brackets usage and warn people not to do that.

fxt44 commented 2 years ago

fmt='(a)' should mean read each string, with spaces as the delimiter, however long the string is declared. this gets away from a fixed format read. i'm also going to test this issue offline with ifort and gfortran.

a more formal way is to read the entire string and parse it oneself. slightly more overhead but this approach allows all kinds of errant entries to be caught before processing further (e.g., non printing control characters).

rjfarmer commented 2 years ago

Gfortran bug report https://gcc.gnu.org/bugzilla/show_bug.cgi?id=104939

It seems like the / is breaking the read statement

rjfarmer commented 2 years ago

Not a bug. Turns out / in unformatted reads tells Fortran to stop reading. Because no one would ever want to put a / in a string.

We'll need to turn this into a formatted read and either build the format string or manually split it.

fxt44 commented 2 years ago

i didn't know that! but then i don't do unformatted reads of character strings. i always use fmt='(a)' or format (a).