Visualisation_of_the_7_ECMWF_models.m

clear%clears workspace
clc %clears command window
fprintf("Do you need color blind support?")
cb_sup=input("\n1-yes\n2-no\ninput:");%user input for colorblind support
%loops while input is invalid
while cb_sup ~= 1 && cb_sup~=2
      cb_sup=input("\ninvalid input\ninput again:");%user input for colorblind support
      if cb_sup == 2
          break
      end
end
cb_input=0;%default value of choice for type of blind-color affection
if cb_sup==1 %checks if supports has been asked for by user
    %prints list of blind-color variants
    fprintf("affection:\n")
    fprintf("1-Protanopia or deuteranopia\n")
    fprintf("2-Tritanopia\n")
    fprintf("0-exit\n")
    cb_input=input("cb input:");%user input for variance or 0-exit in case user did not want colorblind support
    while cb_input ~= 1 && cb_input ~= 2 && cb_input ~= 0%loops while input is invalid
           cb_input=input("cb input again:");%user input for getting out of loop
    end
end 
%% prints the 3 processes which can be started and option for exiting program; The whole program is a single process; however it helps to refer to these options as processes
fprintf("type in:\n")
fprintf(" 1 - visualise one of the seven models at a given hour,\n") 
fprintf(" 2 - visualise the 2 same type figures from 2 climate models at a given hour,\n")
fprintf(" 3 - visualise the 2 same type figures at 2 given hours, \n")
fprintf(" and\n 0 to exit.\n")
ok=input("input:");%user input for choice
while ok<0 || ok>3%loops while user input is invalid
    if ok==0%if ok is 0 then it means user has selected to exit
        fprintf("over\n")
        break
    end
    clc%clears command window
    %prints the 3 processes which can be started and option for exiting program
    fprintf("error!\n")
    fprintf("type in:\n")
    fprintf(" 1 - visualise one of the seven models at a given hour,\n") 
    fprintf(" 2 - visualise the 2 same type figures from 2 climate models at a given hour,\n")
    fprintf(" 3 - visualise the 2 same type figures at 2 given hours, \n")
    fprintf(" and\n 0 to exit.\n")
    ok=input("input:");%user input for choice
end
while ok==1 || ok==2 || ok==3%loops while user still wants to start one of the 3 processes above
    
    NetCDF_list=["Individual_Models\W_fr-meteofrance,MODEL,SILAM+FORECAST+SURFACE+O3+0H24H_C_LFPW_20180701000000.nc",...
            "Individual_Models\W_fr-meteofrance,MODEL,MATCH+FORECAST+SURFACE+O3+0H24H_C_LFPW_20180701000000.nc",...
            "Individual_Models\W_fr-meteofrance,MODEL,ENSEMBLE+FORECAST+SURFACE+O3+0H24H_C_LFPW_20180701000000.nc",...
            "Individual_Models\W_fr-meteofrance,MODEL,CHIMERE+FORECAST+SURFACE+O3+0H24H_C_LFPW_20180701000000.nc",...
            "Individual_Models\W_fr-meteofrance,MODEL,EURAD+FORECAST+SURFACE+O3+0H24H_C_LFPW_20180701000000.nc",...
            "Individual_Models\W_fr-meteofrance,MODEL,MOCAGE+FORECAST+SURFACE+O3+0H24H_C_LFPW_20180701000000.nc",...
            "Individual_Models\W_fr-meteofrance,MODEL,EMEP+FORECAST+SURFACE+O3+0H24H_C_LFPW_20180701000000.nc"];
    meteoModelTitleList=["Meteorologic Model: Silam","Meteorologic Model:Match","Meteorologic Model:Ensemble",...
        "Meteorologic Model:Chimere","Meteorologic Model:Eurad","Meteorologic Model:Mokage","Meteorologic Model:Emep"];
    if ok==0%if ok is 0 then it means user has selected to exit
       fprintf("over\n")
       break
    end
    if ok==1%if ok is 1 then the process corresponding to this number will start
       %models and corresponding number
       fprintf("1=Meteorologic Model Silam\n")
       fprintf("2=Meteorologic Model Match\n")
       fprintf("3=Meteorologic Model Ensemble\n")
       fprintf("4=Meteorologic Model Chimere\n")
       fprintf("5=Meteorologic Model Eurad\n")
       fprintf("6=Meteorologic Model Mocage\n")
       fprintf("7=Meteorologic Model Emep\n")
       m1=input("file:");%user input for model
       while m1<1 || m1>7%loops while usre input is invalid
            fprintf("error \n")
            m1=input("file:");%user input
       end
       %user input for hour at which data was recorded
       h=input("input an integer between 0 and 26 representing 1 of the 25 hours:");
       while h>25 || h<1%loops while user input is invalid
          fprintf("input again\n") 
          %user input for hour at which data was recorded
          h=input("input an integer between 0 and 26 representing 1 of the 25 hours:");
       end 
       for i = 1:length(NetCDF_list)%loops through vector containing names of climate model file
          if i==m1%if input is equal to one of files then
                pnc1 = NetCDF_list(i);%assign name to pnc1
                t1=meteoModelTitleList(i);%assign title of model to t1
          end
       end
          
          figure(1)%creates figure1
          clf%cleares fig
          
          if cb_input == 1%checks for color-blind support
                  colormap Hsv;%colormap colors change 
          else
              if cb_input == 2
                  colormap Jet;%colormap colors change 
              else 
                  colormap; %default colormap
              end
          end
          colorbar;%displays colorbar
          y = double(ncread(pnc1,'longitude'));%creates y variable

          x = double(ncread(pnc1,'latitude'));%creates x variable
          time = ncread(pnc1,'time');%
          ozone = ncread(pnc1,'unknown');%ozone variable containing 4-d array

          z = ozone(:,:,h);%ozone values at a given hour
          title(t1)%titile of climate model
          mesh(x,y,z)%plot mesh
          figure(2)%surm figure
          clf%clears figure
          colorbar;%displays colorbar
          if cb_input == 1%checks for color-blind support
                  colormap Hsv;%colormap colors change 
          else
              if cb_input == 2
                  colormap Jet;%colormap colors change 
              else 
                  colormap ;%default colormap
              end
          end
          title(t1)%displays matching model title
          [x,y] = meshgrid(x, y);%returns 2-D grid coordinates based on the coordinates contained in vectors x and y
          worldmap('Europe'); % set the part of the earth to show
          load coastlines     % loads coaslines
          plotm(coastlat,coastlon)%plots the lat and long of coastlines
          setm(gca,'ffacecolor','#0072BD')%sets ocean color to blue
          %makes land the color green
          land = shaperead('landareas', 'UseGeoCoords', true);
          geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
          %adds wordlakes to the map
          lakes = shaperead('worldlakes', 'UseGeoCoords', true);
          geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
          %adds worldrivers to the map
          rivers = shaperead('worldrivers', 'UseGeoCoords', true);
          geoshow(rivers, 'Color', 'blue')
          %adds worldcities to the map
          cities = shaperead('worldcities', 'UseGeoCoords', true);
          geoshow(cities, 'Marker', 'V', 'Color', 'black')
          %plot the data
          surfm(x,y,z, 'EdgeColor', 'none','FaceAlpha', 0.6)% edge colour outlines the edges, 'FaceAlpha', sets the transparency
           
          figure(3)%contourfm figure
          clf%clear figure

          title(t1)%title of climate model
          colorbar;%displays colorbar
          if cb_input == 1%checks for color blind input
                  colormap Hsv;%colormap colors change 
          else
              if cb_input == 2
                  colormap Jet;%colormap colors change 
              else 
                  colormap ;%default colormap
              end
          end
          worldmap('Europe'); % set the part of the earth to show
          load coastlines     % loads coaslines
          plotm(coastlat,coastlon)%plots the lat and long of coastlines
          setm(gca,'ffacecolor','#0072BD')%sets ocean color to blue
          %makes land the color green
          land = shaperead('landareas', 'UseGeoCoords', true);
          geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
          %adds wordlakes to the map
          lakes = shaperead('worldlakes', 'UseGeoCoords', true);
          geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
          %adds worldrivers to the map
          rivers = shaperead('worldrivers', 'UseGeoCoords', true);
          geoshow(rivers, 'Color', 'blue')
          %adds worldcities to the map
          cities = shaperead('worldcities', 'UseGeoCoords', true);
          geoshow(cities, 'Marker', 'V') 
          NumContours = 5;
          % display the data
          contourfm(x, y, z, NumContours,'--')
          %sets the visibility of the various parts of the plot so the land, cities, rivers, lake etc shows through.
          Plots = findobj(gca,'Type','Axes');%locates graphic object of type axes
          Plots.SortMethod = 'depth';%sorts plots by depth
    end
 
    if ok==2%if ok is 2 then the process corresponding to this number will start
       %user input for type of figure      
       fnum=input("choose kind of fig:1-mesh 2-map 3-countour_map || 0-go back:");
       while fnum<0 || fnum>3 %loops while user input is invalid
           fprintf("input again\n")
           fnum=input("choose kind of fig:1-mesh 2-map 3-countour_map || 0-go back:");%user input
           if fnum==0%if 0 is inputted then exit
             fprintf("exit\n")
             break
           end
       end
       
       while fnum>=1 && fnum<=3%loops while valid input
            %user input for hour at which data was recorded
            h=input("input an integer between 0 and 26 representing 1 of the 25 hours:");
            %loops while user is invalid
            while h>25 || h<1
               fprintf("input again\n") 
               h=input("input an integer between 0 and 26 representing 1 of the 25 hours:");%user input
            end
            %prints list of climate model names
            fprintf("1=Meteorologic Model Silam\n")
            fprintf("2=Meteorologic Model Match\n")
            fprintf("3=Meteorologic Model Ensemble\n")
            fprintf("4=Meteorologic Model Chimere\n")
            fprintf("5=Meteorologic Model Eurad\n")
            fprintf("6=Meteorologic Model Mocage\n")
            fprintf("7=Meteorologic Model Emep\n")
            m1=input("f1:");%user input for file
            
            while m1<1 || m1>7%loops while invalid user input
                fprintf("error \n")
                m1=input("f1:");%user input
            end
            m2=input("f2:");%second user input for 2nd file
            
            while m2<1 || m2>7%loops while invalid user input
                fprintf("error \n")
                m2=input("f2:");%user input
            end
            for i = 1:length(NetCDF_list)%loops through vector containing names of climate model file
                if i==m1%checks if file1 exists in list of climate model names
                    pnc1=NetCDF_list(i);%assigns name to pnc1
                    t1=meteoModelTitleList(i);%assigns title to t1
                end
                if i==m2%checks if file2 exists in list of climate model names
                    pnc2=NetCDF_list(i);%assigns name to pnc1
                    t2=meteoModelTitleList(i);%assigns title to t1
                end
            end    
            %checks if type of figure is a mesh and plots the data
            if fnum==1
               figure(1)%creates figure
               clf
               colorbar;
               %color blind support check
               if cb_input == 1
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               %creates x1,y1 and z1 variables
               y1 = double(ncread(pnc1,'longitude'));
               x1 = double(ncread(pnc1,'latitude'));
               time1 = ncread(pnc1,'time');
               ozone1 = ncread(pnc1,'unknown');

               z1 = ozone1(:,:,h);
               title(t1)
               mesh(x1,y1,z1)%plot data
               

               figure(2)%create 2nd figure
               clf
               colorbar;
               %color blind support check
               if cb_input == 1
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               %creates x2,y2 and z2 variables
               y2 = double(ncread(pnc2,'longitude'));
               zlabel('ozone')
               x2 = double(ncread(pnc2,'latitude'));
               time2 = ncread(pnc2,'time');
               ozone2 = ncread(pnc2,'unknown');
               z2 = ozone2(:,:,h);
               title(t2)
               mesh(x2,y2,z2)%plot data
            end
            %checks if type of figure is surfm
            if fnum==2
               figure(3)%creates figure
               clf
               colorbar;
               %color blind support check
               if cb_input == 1
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               title(t1)
               %creates x1,y1 and z1 variables
               y1 = double(ncread(pnc1,'longitude'));
               zlabel('ozone')
               x1 = double(ncread(pnc1,'latitude'));
               time1 = ncread(pnc1,'time');
               ozone1 = ncread(pnc1,'unknown');

               z1 = ozone1(:,:,h);
               [x1,y1] = meshgrid(x1, y1);%returns 2-D grid coordinates based on the coordinates contained in vectors x and y
               worldmap('Europe'); % set the part of the earth to show
               load coastlines
               plotm(coastlat,coastlon)
               %adds landares,lakes,rivers and cities
               land = shaperead('landareas', 'UseGeoCoords', true);
               geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
               lakes = shaperead('worldlakes', 'UseGeoCoords', true);
               geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
               rivers = shaperead('worldrivers', 'UseGeoCoords', true);
               geoshow(rivers, 'Color', 'blue')
               cities = shaperead('worldcities', 'UseGeoCoords', true);
               geoshow(cities, 'Marker', 'V', 'Color', 'black')
               setm(gca,'ffacecolor','#0072BD')%set color of water bodies to blue
               %plot the data
               surfm(x1,y1,z1, 'EdgeColor', 'none','FaceAlpha', 0.6)% edge colour outlines the edges, 'FaceAlpha', sets the transparency

               figure(4)
               clf
               colorbar; 
               %color blind support check
               if cb_input == 1
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default map
                  end
               end
               title(t2)
               %creates x2,y2 and z2 variables
               y2 = double(ncread(pnc2,'longitude'));
               zlabel('ozone')
               x2 = double(ncread(pnc2,'latitude'));
               time2 = ncread(pnc2,'time');
               ozone2 = ncread(pnc2,'unknown');
               z2 = ozone2(:,:,h);
               [x2,y2] = meshgrid(x2, y2);%returns 2-D grid coordinates based on the coordinates contained in vectors x and y
               
               worldmap('Europe'); % set the part of the earth to show
               load coastlines
               plotm(coastlat,coastlon)
               %adds landares,lakes,rivers and cities
               land = shaperead('landareas', 'UseGeoCoords', true);
               geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
               lakes = shaperead('worldlakes', 'UseGeoCoords', true);
               geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
               rivers = shaperead('worldrivers', 'UseGeoCoords', true);
               geoshow(rivers, 'Color', 'blue')
               cities = shaperead('worldcities', 'UseGeoCoords', true);
               geoshow(cities, 'Marker', 'V', 'Color', 'black')
               setm(gca,'ffacecolor','#0072BD')%set color of water body to blue
               %plot the data
               surfm(x2,y2,z2, 'EdgeColor', 'none','FaceAlpha', 0.6)% edge colour outlines the edges, 'FaceAlpha', sets the transparency
            end
            if fnum == 3%checks if type of figure is contourfm
               %creates x1,y1 and z1 variables
               y1 = double(ncread(pnc1,'longitude'));
               x1 = double(ncread(pnc1,'latitude'));
               time1 = ncread(pnc1,'time');
               ozone1 = ncread(pnc1,'unknown');
               [x1,y1] = meshgrid(x1, y1);%returns 2-D grid coordinates based on the coordinates contained in vectors x and y
               z1 = ozone1(:,:,h);
               figure(5) 
               clf
               colorbar;
               if cb_input == 1%color blind support check
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               title(t1)

               worldmap('Europe'); % set the part of the earth to show
               load coastlines
               plotm(coastlat,coastlon)
               %adds landares,lakes,rivers and cities
               setm(gca,'ffacecolor','#0072BD')
               land = shaperead('landareas', 'UseGeoCoords', true);
               geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
               lakes = shaperead('worldlakes', 'UseGeoCoords', true);
               geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
               rivers = shaperead('worldrivers', 'UseGeoCoords', true);
               geoshow(rivers, 'Color', 'blue')
               cities = shaperead('worldcities', 'UseGeoCoords', true);
               geoshow(cities, 'Marker', 'V', 'Color', 'black')

               NumContours1 = 6;
               contourfm(x1, y1, z1, NumContours1,'--')%sets the visibility of the various parts of the plot so the land, cities, rivers, lake etc shows through.
               Plots = findobj(gca,'Type','Axes');%locates graphic object of type axes
               Plots.SortMethod = 'depth';%sorts plots by depth
               %creates x2,y2 and z2 variables
               y2 = double(ncread(pnc2,'longitude'));
               x2 = double(ncread(pnc2,'latitude'));
               time2 = ncread(pnc2,'time');
               ozone2 = ncread(pnc2,'unknown');
               z2 = ozone2(:,:,h);
               [x2,y2] = meshgrid(x2, y2);%returns 2-D grid coordinates based on the coordinates contained in vectors x and y
                             
               figure(6) 
               clf
               colorbar;
               if cb_input == 1%color blind support check
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               title(t2)

               worldmap('Europe'); % set the part of the earth to show
               load coastlines
               plotm(coastlat,coastlon)
               %adds landares,lakes,rivers and cities
               setm(gca,'ffacecolor','#0072BD')
               land = shaperead('landareas', 'UseGeoCoords', true);
               geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
               lakes = shaperead('worldlakes', 'UseGeoCoords', true);
               geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
               rivers = shaperead('worldrivers', 'UseGeoCoords', true);
               geoshow(rivers, 'Color', 'blue')
               cities = shaperead('worldcities', 'UseGeoCoords', true);
               geoshow(cities, 'Marker', 'V', 'Color', 'black')

               NumContours2 = 5;
               contourfm(x2, y2, z2, NumContours2,'--')%sets the visibility of the various parts of the plot so the land, cities, rivers, lake etc shows through.
               Plots = findobj(gca,'Type','Axes');%locates graphic object of type axes
               Plots.SortMethod = 'depth';%sorts plots by depth
               
            end
            fnum=input("choose kind of fig:1-mesh 2-map 3-countour_map || 0-go back:");%user input for continuing process or go back to previous menu
            while fnum<0 || fnum>3 %loops while user input is invalid
                fprintf("input again\n")
                fnum=input("choose kind of fig:1-mesh 2-map 3-countour_map || 0-go back:");%user input
                if fnum==0%if input is 0 then return to previous menu
                  fprintf("exit\n")
                  break
                end
            end
        end
       
    end
    if ok==3%3rd process that can be started
       
       clc%clears command window
       %user input for type of figure 
       fnum=input("choose kind of fig:1-mesh 2-map 3-countour_map || 0-go back:");
       while fnum<0 || fnum>3 %loops while invalid user input
           fprintf("input again\n")
           fnum=input("choose kind of fig:1-mesh 2-map 3-countour_map || 0-go back:");%user input
           if fnum==0%option for returning to previous menu
             clc
             fprintf("exit\n")
             break
           end
       end
       while fnum>=1 && fnum<=3
           %input for first hour 
           h1=input("input an integer between 0 and 26 representing 1 of the 25 hours:");
           while h1>25 || h1<1%loops while invalid input
             fprintf("input again\n") 
             h1=input("input an integer between 0 and 26 representing 1 of the 25 hours:");%user input
           end
           %input for second hour
           h2=input("input a second integer between 0 and 26 representing 1 of the 25 hours:");
           while h2>25 || h2<1%loops while invalid input
             clc
             fprintf("input again\n") 
             h2=input("input an integer between 0 and 26 representing 1 of the 25 hours:");%user input
           end
           if fnum==1%option for mesh figure
               clc%clears command window
               fprintf("1=Meteorologic Model Silam\n")
               fprintf("2=Meteorologic Model Match\n")
               fprintf("3=Meteorologic Model Ensemble\n")
               fprintf("4=Meteorologic Model Chimere\n")
               fprintf("5=Meteorologic Model Eurad\n")
               fprintf("6=Meteorologic Model Mocage\n")
               fprintf("7=Meteorologic Model Emep\n")
               %input number coresponding to model 
               m1=input("File:");

               while m1<1 || m1>7%loops while invalid user input
                   clc
                   fprintf("error \n")
                   fprintf("1=Meteorologic Model Silam\n")
                   fprintf("2=Meteorologic Model Match\n")
                   fprintf("3=Meteorologic Model Ensemble\n")
                   fprintf("4=Meteorologic Model Chimere\n")
                   fprintf("5=Meteorologic Model Eurad\n")
                   fprintf("6=Meteorologic Model Mocage\n")
                   fprintf("7=Meteorologic Model Emep\n")

                   m1=input("File:");%user input
               end
               for i = 1:length(NetCDF_list)%loops through NetCDF_list to find inputted file
                  if i==m1
                        pnc1=NetCDF_list(i);
                        t1=meteoModelTitleList(i);
                  end
               end
               %creates x,y variables
               y = double(ncread(pnc1,'longitude'));
               zlabel('ozone')
               x = double(ncread(pnc1,'latitude'));
               time = ncread(pnc1,'time');
               ozone = ncread(pnc1,'unknown');%ozone variable
               %creates z1 and z2 variables
               z1 = ozone(:,:,h1);
               z2 = ozone(:,:,h2);
               
               figure(1)
               clf
               colorbar;
               if cb_input == 1%color-blind support check
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%colormap default
                  end
               end
               title(t1)
               mesh(x,y,z1)%plot data
               figure(2)%create figure
               clf
               colorbar;
               if cb_input == 1%color-blind support check
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               title(t1)
               mesh(x,y,z2)%plot data
           end
           if fnum==2%sufm type figures option
               clc
               %prints list of models and corresponding number
               fprintf("1=Meteorologic Model Silam\n")
               fprintf("2=Meteorologic Model Match\n")
               fprintf("3=Meteorologic Model Ensemble\n")
               fprintf("4=Meteorologic Model Chimere\n")
               fprintf("5=Meteorologic Model Eurad\n")
               fprintf("6=Meteorologic Model Mocage\n")
               fprintf("7=Meteorologic Model Emep\n")
               %input file
               m1=input("File:");

               while m1<1 || m1>7%loops while invalid input
                   clc
                   fprintf("error \n")
                   fprintf("1=Meteorologic Model Silam\n")
                   fprintf("2=Meteorologic Model Match\n")
                   fprintf("3=Meteorologic Model Ensemble\n")
                   fprintf("4=Meteorologic Model Chimere\n")
                   fprintf("5=Meteorologic Model Eurad\n")
                   fprintf("6=Meteorologic Model Mocage\n")
                   fprintf("7=Meteorologic Model Emep\n")

                   m1=input("File:");%user input
               end
               for i = 1:length(NetCDF_list)%loops through NetCDF_list to find correct net
                  if i==m1%check if file exists
                        pnc1=NetCDF_list(i);%assign name to pnc1
                        t1=meteoModelTitleList(i);%assign title
                  end
               end
               %create x and y values
               y = double(ncread(pnc1,'longitude'));
               zlabel('ozone')
               x = double(ncread(pnc1,'latitude'));
               time = ncread(pnc1,'time');
               ozone = ncread(pnc1,'unknown');
               %creates z1,z2 values
               z1 = ozone(:,:,h1);                     
               z2 = ozone(:,:,h2);
               [x,y] = meshgrid(x, y);%returns 2-D grid coordinates based on the coordinates contained in vectors x and y
               figure(3)
               
               clf 
               if cb_input == 1%color-blind support check
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               title(t1)
               colorbar;
               worldmap('Europe'); % set the part of the earth to show
               load coastlines
               plotm(coastlat,coastlon)
               %adds landares,lakes,rivers and cities
               land = shaperead('landareas', 'UseGeoCoords', true);
               geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
               lakes = shaperead('worldlakes', 'UseGeoCoords', true);
               geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
               rivers = shaperead('worldrivers', 'UseGeoCoords', true);
               geoshow(rivers, 'Color', 'blue')
               cities = shaperead('worldcities', 'UseGeoCoords', true);
               geoshow(cities, 'Marker', 'V', 'Color', 'black')
               setm(gca,'ffacecolor','#0072BD')
               %Plot the data
               surfm(x,y,z1, 'EdgeColor', 'none','FaceAlpha', 0.6)
               % edge colour outlines the edges, 'FaceAlpha', sets the transparency

               figure(4)
               clf
               if cb_input == 1%color-blind support check
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               title(t1)
               colorbar;
              
               worldmap('Europe'); % set the part of the earth to show
               load coastlines
               plotm(coastlat,coastlon)
               %adds landares,lakes,rivers and cities
               land = shaperead('landareas', 'UseGeoCoords', true);
               geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
               lakes = shaperead('worldlakes', 'UseGeoCoords', true);
               geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
               rivers = shaperead('worldrivers', 'UseGeoCoords', true);
               geoshow(rivers, 'Color', 'blue')
               cities = shaperead('worldcities', 'UseGeoCoords', true);
               geoshow(cities, 'Marker', 'V', 'Color', 'black')
               setm(gca,'ffacecolor','#0072BD')
               %Plot the data
               surfm(x,y,z2, 'EdgeColor', 'none','FaceAlpha', 0.6) 
               % edge colour outlines the edges, 'FaceAlpha', sets the transparency
           end
           if fnum==3%option for choosing contourfm type figure
               clc
               %prints list of models and corresponding number
               fprintf("1=Meteorologic Model Silam\n")
               fprintf("2=Meteorologic Model Match\n")
               fprintf("3=Meteorologic Model Ensemble\n")
               fprintf("4=Meteorologic Model Chimere\n")
               fprintf("5=Meteorologic Model Eurad\n")
               fprintf("6=Meteorologic Model Mocage\n")
               fprintf("7=Meteorologic Model Emep\n")
               
               m1=input("File:");%user input for file

               while m1<1 || m1>7%loops while invalid user input
                   clc
                   %prints list of models and corresponding number
                   fprintf("error \n")
                   fprintf("1=Meteorologic Model Silam\n")
                   fprintf("2=Meteorologic Model Match\n")
                   fprintf("3=Meteorologic Model Ensemble\n")
                   fprintf("4=Meteorologic Model Chimere\n")
                   fprintf("5=Meteorologic Model Eurad\n")
                   fprintf("6=Meteorologic Model Mocage\n")
                   fprintf("7=Meteorologic Model Emep\n")

                   m1=input("File:");%user input
               end
               

               for i = 1:length(NetCDF_list)%loops through vector containing names of climate model file
                  if i==m1%checks if num inputted corresponds to a model
                        pnc1=NetCDF_list(i);%assigns name to pnc1
                        t1=meteoModelTitleList(i);%assigns corresponding title
                  end
               end
               %creates x and y values
               y = double(ncread(pnc1,'longitude'));
               zlabel('ozone')
               x = double(ncread(pnc1,'latitude'));
               time = ncread(pnc1,'time');
               ozone = ncread(pnc1,'unknown');
               %creates z1 and z2 values
               z1 = ozone(:,:,h1);
               z2 = ozone(:,:,h2);

               [x,y] = meshgrid(x, y);%returns 2-D grid coordinates based on the coordinates contained in vectors x and y
               
               figure(5) 
               clf
               colorbar;
               title(t1)
               if cb_input == 1%color-blind support check
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end
               worldmap('Europe'); % set the part of the earth to show
               load coastlines
               plotm(coastlat,coastlon)
               %adds landares,lakes,rivers and cities
               setm(gca,'ffacecolor','#0072BD')
               land = shaperead('landareas', 'UseGeoCoords', true);
               geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
               lakes = shaperead('worldlakes', 'UseGeoCoords', true);
               geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
               rivers = shaperead('worldrivers', 'UseGeoCoords', true);
               geoshow(rivers, 'Color', 'blue')
               cities = shaperead('worldcities', 'UseGeoCoords', true);
               geoshow(cities, 'Marker', 'V', 'Color', 'black')
               NumContours1=5;
               contourfm(x, y, z1, NumContours1,'--')%sets the visibility of the various parts of the plot so the land, cities, rivers, lake etc shows through.
               Plots = findobj(gca,'Type','Axes');%locates graphic object of type axes
               Plots.SortMethod = 'depth';%sorts plots by depth
               
               figure(6) 
               clf
               colorbar;
               title(t1)
               if cb_input == 1%color-blind support check
                  colormap Hsv;%colormap colors change
               else
                  if cb_input == 2
                      colormap Jet;%colormap colors change
                  else 
                      colormap ;%default colormap
                  end
               end   
               worldmap('Europe'); % set the part of the earth to show
               load coastlines
               plotm(coastlat,coastlon)
               setm(gca,'ffacecolor','#0072BD')
               land = shaperead('landareas', 'UseGeoCoords', true);
               geoshow(gca, land, 'FaceColor', [0.5 0.7 0.5])
               lakes = shaperead('worldlakes', 'UseGeoCoords', true);
               geoshow(lakes, 'FaceColor', 'blue', 'EdgeColor', 'cyan')
               rivers = shaperead('worldrivers', 'UseGeoCoords', true);
               geoshow(rivers, 'Color', 'blue')
               cities = shaperead('worldcities', 'UseGeoCoords', true);
               geoshow(cities, 'Marker', 'V', 'Color', 'black')
               NumContours2 = 5;
               contourfm(x, y, z2, NumContours2,'LineStyle',':')%sets the visibility of the various parts of the plot so the land, cities, rivers, lake etc shows through.
               Plots = findobj(gca,'Type','Axes');%locates graphic object of type axes
               Plots.SortMethod = 'depth';%sorts plots by depth
           end
           fnum=input("choose kind of fig:1-mesh 2-map 3-countour_map || 0-go back:");%nput for repeating process
           while fnum<0 || fnum>3 %loops while invalid input
               fprintf("input again\n")
               fnum=input("choose kind of fig:1-mesh 2-map 3-countour_map || 0-go back:");%type of figure input
               if fnum==0%option for returning to previous menu
                 fprintf("exit\n")
                 break
               end
           end
       end
       
    end
    clc%clears command window
    %prints list of options
    fprintf("type in:\n")
    fprintf(" 1 - visualise one of the seven models at a given hour,\n") 
    fprintf(" 2 - visualise the 2 same type figures from 2 climate models at a given hour,\n")
    fprintf(" 3 - visualise the 2 same type figures at 2 given hours, \n")
    fprintf(" and\n 0 to exit.\n")
    ok=input("input:");%user input for option
    while ok<0 || ok>3%loops while invalid user input
        if ok==0%breaks if option inputted is 0
            fprintf("over\n")
            break
        end
        clc
        %prints list of options
        fprintf("type in:\n")
        fprintf(" 1 - visualise one of the seven models at a given hour,\n") 
        fprintf(" 2 - visualise the 2 same type figures from 2 climate models at a given hour,\n")
        fprintf(" 3 - visualise the 2 same type figures at 2 given hours, \n")
        fprintf(" and\n 0 to exit.\n")
        ok=input("input:");%user input
    end
end
clear
close all
%end of program