debug

sed_cores/options.py

 extract_adists_data = True
-adists_bin_data_dir = "./adists_res"
+adists_bin_data_dir = "."
 adists_extracted_data_dir = "./adists_out"
 adists_exe = "adists"
 output_accumulation_to_csv = True

sed_cores/sediment_cores.py

@@ -5,9 +5,9 @@ except:
     print("Error: sys not found")
     exit()
 try:
-    import importlib as imp
+    import imp
 except:
-    print("Error: importlib not found")
+    print("Error: imp not found")
     exit()
 try:
     import getopt
@@ -120,9 +120,12 @@ if not os.path.isfile(options_path):
 else:
     option_file = options_path
 
-imp.machinery.SourceFileLoader("options", option_file)
+print("we are in",os.getcwd())
+print("option file is ",option_file)
+imp.load_source("options", option_file)
 from options import extract_adists_data
-from options import adists_data_dir
+from options import adists_bin_data_dir
+from options import adists_extracted_data_dir
 from options import adists_exe
 from options import output_accumulation_to_csv
 from options import output_accumulation_to_png
@@ -142,6 +145,8 @@ values = dict()
 nb_pk = len(list_pk)
 nb_files = len(list_files)
 nb_sediments = len(list_sediments)
+cwd = os.getcwd()
+print("we are in",os.getcwd())
 
 #===========#
 # read data #
@@ -161,19 +166,20 @@ if extract_adists_data:
 
 if extract_adists_data:
     if (not os.path.isdir(adists_extracted_data_dir)):
-        print("Error: "+adists_extracted_data_dir+" folder not found")
-        exit()
+        os.mkdir(adists_extracted_data_dir)
     os.chdir(adists_extracted_data_dir)
     print("Extracting data from Adis-TS to CSV files...")
+    print("we are in",os.getcwd())
 
     output = ""
     if not verbose: output = " > /dev/null"
     for bed in list_beds:
         for file in list_files:
             for i,pk in enumerate(list_pk):
+                print("Extracting data from "+adists_bin_data_dir+"/"+file+".bin")
                 os.system( adists+" -b "+adists_bin_data_dir+"/"+file+".bin "+bed+"dt "+str(pk)+" "+file+"_"+bed+str(i)+output)
 
-    os.chdir("..")
+    os.chdir(cwd)
 
 print("Reading CSV files...")
 
@@ -212,7 +218,7 @@ for bed in list_beds:
                 for t,row in enumerate(csv_reader):
                     values.get(bed)[j][i][t] = float(row[1])
             f.close()
-os.chdir("..")
+os.chdir(cwd)
 if verbose: print("  Reading done")
 
 #========================#
@@ -247,13 +253,13 @@ for bed in list_beds:
         val_acc[pk][0] = val[0][pk][0]
         val_acc_tot[pk][0] = val[0][pk][0]
         for t in range(nb_timesteps-1):
+            val_acc[pk][t+1] = val[0][pk][t+1] - val[0][pk][t]
             if val[0][pk][t+1]>val[0][pk][t]:
-                val_acc[pk][t+1] = val[0][pk][t+1] - val[0][pk][t]
                 layers[pk] += 1
+                val_acc_tot[pk][t+1] = val_acc_tot[pk][t] + val_acc[pk][t+1]
             else:
-                val_acc[pk][t+1] = 0
-            val_acc_tot[pk][t+1] = val_acc_tot[pk][t] + val_acc[pk][t+1]
-os.chdir("..")
+                val_acc_tot[pk][t+1] = val_acc_tot[pk][t]
+os.chdir(cwd)
 if verbose: print("  Computing accumulation done")
 
 if output_accumulation_to_csv:
@@ -265,9 +271,9 @@ if output_accumulation_to_csv:
     os.chdir("./accumulation")
 
     if def_numpy:
-        row_data=np.zeros((nb_timesteps,3),float)
+        row_data=np.zeros((nb_timesteps,4),float)
     else:
-        row_data = create_float_list([nb_timesteps,3])
+        row_data = create_float_list([nb_timesteps,4])
     for i in range(nb_timesteps):
         row_data[i][0]=list_timesteps[i]
 
@@ -280,16 +286,17 @@ if output_accumulation_to_csv:
             for t in range(nb_timesteps):
                 row_data[t][1]=val_acc[pk][t]
                 row_data[t][2]=val_acc_tot[pk][t]
+                row_data[t][3]=val[0][pk][t]
             if layers[pk]>0:
                 if verbose: print("  Writing values_accumulated_"+list_pk[pk]+"_"+bed+".csv")
                 with open("values_accumulated_"+list_pk[pk]+"_"+bed+".csv", "w") as f:
                     writer = csv.writer(f)
-                    writer.writerow(["timestep","val_acc_instant","val_acc_tot"])
+                    writer.writerow(["timestep","instant_accumulation","positive accumulation","total_accumulation"])
                     for t in range(nb_timesteps):
                         writer.writerow(row_data[t][:])
                 f.close()
     if verbose: print("  Writing accumulation done")
-    os.chdir("..")
+    os.chdir(cwd)
 
 if def_matplotlib and output_accumulation_to_png:
 
@@ -301,6 +308,7 @@ if def_matplotlib and output_accumulation_to_png:
         val_acc = values_accumulated.get(bed)
         val_acc_tot = values_accumulated_tot.get(bed)
         layers = nb_layers.get(bed)
+        val = values.get(bed)
         for pk in range(nb_pk):
             if layers[pk]>0:
                 if verbose: print("  Plotting values_accumulated_"+list_pk[pk]+"_"+bed+".png")
@@ -309,15 +317,16 @@ if def_matplotlib and output_accumulation_to_png:
                 ax.set_xlabel("time")
                 ax.set_ylabel("height")
                 #ax.set_title("Sediment accumulated at pk"+list_pk[pk]+" in bed "+bed)
-                plt.plot(list_timesteps[:],val_acc[pk][:],label="instant accumulated")
-                plt.plot(list_timesteps[:],val_acc_tot[pk][:],label="total accumulated")
+                plt.plot(list_timesteps[:],val_acc[pk][:],label="instant accumulation")
+                plt.plot(list_timesteps[:],val_acc_tot[pk][:],label="positive accumulation")
+                plt.plot(list_timesteps[:],val[0][pk][:],label="total accumulation")
 
                 plt.legend()
                 plt.title(label="Sediment accumulated at pk"+list_pk[pk]+" in bed "+bed+"\n", fontsize = 20)
                 fig.savefig("values_accumulated_"+list_pk[pk]+"_"+bed+".png")
                 plt.close(fig)
     if verbose: print("  Plotting accumulation done")
-    os.chdir("..")
+    os.chdir(cwd)
 
 ####################
 # % of sed classes #
@@ -361,7 +370,7 @@ for bed in list_beds:
                     val_temp1[s] = sum(val[s+j][pk][t] for j in range(1,nb_files,nb_sediments))
                     if val_temp2>0.0 and val_temp1[s]>0.0 :
                         layers[bed][pk][l-1][s+2] = val_temp1[s]/val_temp2
-os.chdir("..")
+os.chdir(cwd)
 if verbose: print("Computing sediment cores done")
 
 if output_cores_to_csv:
@@ -383,7 +392,7 @@ if output_cores_to_csv:
                         writer.writerow(layers[bed][pk][i][:])
                 f.close()
     if verbose: print("  Writing sediment cores done")
-    os.chdir("..")
+    os.chdir(cwd)
 
 if def_matplotlib and output_cores_to_png:
     print("Output sediment cores as png...")
@@ -422,7 +431,7 @@ if def_matplotlib and output_cores_to_png:
                 fig.savefig("sed_core_"+list_pk[pk]+"_"+bed+".png")
                 plt.close(fig)
     if verbose: print("  Plotting sediment cores done")
-    os.chdir("..")
+    os.chdir(cwd)
 
 print("Done !")