tadd new figure with climate record - cosmo - front and backend for Markov-Chain Monte Carlo inversion of cosmogenic nuclide concentrations

git clone git://src.adamsgaard.dk/cosmo

Log

Files

Refs

README

LICENSE

---

commit e25b3d107641ec4057698d8e444d95e46a6c59dd

parent a0b34d9dea120b697868fbd2aef1004cb4b6497e

Author: Anders Damsgaard 

Date:   Fri, 23 Oct 2015 13:42:13 +0200

add new figure with climate record

Diffstat:
  M index.php                           |      36 ++++++++++++++++++++++++++++++-
  M matlab/file_scanner_mcmc_starter.m  |       3 ++-
  M matlab/generate_plots.m             |     107 ++++++++++++++++++++++++++-----

3 files changed, 127 insertions(+), 19 deletions(-)
---

diff --git a/index.php b/index.php

t@@ -76,7 +76,7 @@ if (isset($_GET['wait_id']) && !empty($_GET['wait_id'])) {
                   

                   
-                    Fourth figure
+                    Model parameter values
                   
                   

                       
The histograms show the distribution of (a)
t@@ -105,6 +105,40 @@ if (isset($_GET['wait_id']) && !empty($_GET['wait_id'])) {
               
             
 
+            

+              

+                

+                  

+                  
+                    Erosion history
+                  
+                  

+                      

+                      (a) Two-stage glacial-interglacial model based on the
+                      selected climate record.
+                      
+                      (b) Timing and duration of glacial and interglacial
+                      periods are defined by a threshold value (d18Othreshold)
+                      that is applied to the climate record. 
+                      
+                      
+                  
+                  

+                  Link to PNG
+                  Link to PDF
+                  Link to FIG
+                  
+                
+              
+            
+
 

diff --git a/matlab/file_scanner_mcmc_starter.m b/matlab/file_scanner_mcmc_starter.m

t@@ -106,7 +106,8 @@ while 1
         
         % generate plots
         %CompareWalks2(Ss, save_file)
-        generate_plots(Ss, save_file, graphics_formats, show_figures);
+        generate_plots(Ss, record, matlab_scripts_folder, ...
+            save_file, graphics_formats, show_figures);
         
         % close all figures
         close all;

diff --git a/matlab/generate_plots.m b/matlab/generate_plots.m

t@@ -1,4 +1,5 @@
-function generate_plots(Ss, save_file, formats, show_figures)
+function generate_plots(Ss, record, matlab_scripts_folder,...
+    save_file, formats, show_figures)
 
 %% Copied from m_pakke2014maj11/CompareWalks2.m
 % Generates and saves all relevant figures
t@@ -161,7 +162,7 @@ for i1=1:M
   end
 end
 
-% Histogram plots for all four parameters, one subplot per walker, fh(4)
+%% Histogram plots for all four parameters, one subplot per walker, fh(4)
 fh = [fh;figure('visible', show_figures)];
 for i1 = 1:M
   for iwalk=1:min(4,Nwalkers)
t@@ -180,8 +181,24 @@ for i1 = 1:M
   end
 end
 
+%% Titles for Mads' plots
+if titles
+    %Putting in titles over figure 2:4
+    figure(fh(2)); set(fh(2), 'Visible', show_figures)
+    subplot(5,4,1)
+    title(['Density cross-plots A. Result file = ',save_file])
+
+    figure(fh(3)); set(fh(3), 'Visible', show_figures)
+    subplot(5,4,1)
+    title(['Density cross-plots B. Result file = ',save_file])
+
+    figure(fh(4)); set(fh(4), 'Visible', show_figures)
+    subplot(M,Nwalkers,1)
+    %title(['Histograms. Result file =',save_file],'interp','none')
+    title('Distribution of model parameter values')
+end
 
-% Plot one histogram per model parameter per walker
+%% Plot one histogram per model parameter per walker
 fh = [fh;figure('visible', show_figures)];
 for i1 = 1:M % for each model parameter
   for iwalk=1:min(4,Nwalkers) % for each walker
t@@ -224,23 +241,79 @@ for i1 = 1:M % for each model parameter
   end
 end
 
-if titles
-    %Putting in titles over figure 2:4
-    figure(fh(2)); set(fh(2), 'Visible', show_figures)
-    subplot(5,4,1)
-    title(['Density cross-plots A. Result file = ',save_file])
 
-    figure(fh(3)); set(fh(3), 'Visible', show_figures)
-    subplot(5,4,1)
-    title(['Density cross-plots B. Result file = ',save_file])
 
-    figure(fh(4)); set(fh(4), 'Visible', show_figures)
-    subplot(M,Nwalkers,1)
-    %title(['Histograms. Result file =',save_file],'interp','none')
-    title('Distribution of model parameter values')
+
+%% Plot d18O curve, from PlotSmoothZachos.m
+fh = [fh;figure('visible', show_figures)];
+
+if strcmp(record, 'rec_5kyr')
+    d18O_filename = 'lisiecki_triinterp_2p6Ma_5ky.mat'; %  zachos_triinterp_2p6Ma
+elseif strcmp(record, 'rec_20kyr')
+    d18O_filename = 'lisiecki_triinterp_2p6Ma_20ky.mat'; %  zachos_triinterp_2p6Ma
+elseif strcmp(record, 'rec_30kyr')
+    d18O_filename = 'lisiecki_triinterp_2p6Ma_30ky.mat'; %  zachos_triinterp_2p6Ma
+else
+    error(['record ' record ' not understood']);
 end
 
-% position figure windows at certain coordinates on the screen
+load([matlab_scripts_folder, d18O_filename])
+xs = ti;
+dt = diff(ti(1:2));
+ys = d18O_triang;
+%colpos = [0.5,0.5,1];
+colpos = [0,0,1];
+%colneg = [0.5,1,0.5];
+colneg = [1,0,0];
+midvalue = 3.75;
+axh(1)=subplot(3,1,1);
+[~,~,~,i_cross]=fill_red_blue(xs,ys,colpos,colneg,midvalue);
+ylabel('\delta18O')
+axis tight
+axis([-0.1,2.7,3.0,5.2])
+%axis([-0.05,0.3,3.0,5.2])
+xlim([-0.1,2.7])
+%xlim([-0.05,0.3])
+hold on
+%plot(A(:,1)/1000,A(:,2),'.-m')
+axis ij
+xs_cross = (i_cross-1)*dt;
+xs_cross = [0;xs_cross];
+%title([fn,'.  minrad = ',num2str(minrad),' Ma'],'interp','none')
+
+xs_cross(2) = 0.011;
+
+axh(2)=subplot(3,1,2);
+stairs(xs_cross,(1+-1*(-1).^(1:length(xs_cross)))/2,'b','linewidth',1.5);
+hold on
+start1 = [xs_cross(end),2.7];
+start2 = [1,1];
+plot(start1,start2,'b','linewidth',1.5);
+%title('Exposure. 0 = glaciated, 1 = not glaciated')
+axis([-0.1,2.7,-0.5,1.5])
+%axis([-0.05,0.3,-0.5,1.5])
+xlabel('Age BP [Ma]')
+
+axh(2)=subplot(3,1,3);
+stairs(xs_cross,(1+-1*(-1).^(1:length(xs_cross)))/2,'r','linewidth',1.5);
+hold on
+plot(start1,start2,'r','linewidth',1.5);
+%title('Exposure. 0 = glaciated, 1 = not glaciated')
+axis([-0.1,2.7,-1,2])
+%axis([-0.05,0.3,-1,2])
+xlabel('Age BP [Ma]')
+
+
+hold on
+d18Oth = midvalue;
+ErateInt=1e-6; ErateGla=1e-7;
+[~,~] = ExtractHistory2(ti,d18O_triang,d18Oth,ErateInt,ErateGla);
+
+linkaxes(axh,'x')
+% stairs(-tStarts,relExpos+2,'r')
+
+
+%% position figure windows at certain coordinates on the screen
 if strcmp(show_figures, 'on')
     figpos1 = [6         474        1910         504];
     figpos2 =[    12 94 645 887];
t@@ -253,7 +326,7 @@ if strcmp(show_figures, 'on')
     figure(fh(1))
 end
 
-% save all figures
+%% save all figures
 for i=1:length(fh)
     figure_save_multiformat(fh(i), ...
         strcat(save_file, '-', num2str(i)), ...