diff --git a/GIV_GUI_initialize.m b/GIV_GUI_initialize.m index f9392ef..4063c1c 100644 --- a/GIV_GUI_initialize.m +++ b/GIV_GUI_initialize.m @@ -1,3 +1,4 @@ +<<<<<<< HEAD function GIV_GUI_initialize() % % @@ -1254,3 +1255,1260 @@ function GIV_GUI_initialize() end +======= +function GIV_GUI_initialize() +% +% +%This function starts off the main image velocity calculation code. It will +%open a graphical user interface where the different input values can be +%entered, and then run the functions to calcuate velocities for each image +%pair inputted. +% +% +%Note this can take a long time (hours) for very large datasets, it is +%sometimes worth testing it with smaller datasets to begin with (and/or +%using low temporal oversampling values). +% +%Most of this very long file is simply setting the location of the boxes, +%etc. Portions were generated automatically using MATLAB's app builder. You +%can tweak it to improve layout if you like, but I would not recommend +%changing it very much. + + + + + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %% GLACIER IMAGE VELOCIMETRY (GIV) %% +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +%Code written by Max Van Wyk de Vries @ University of Minnesota +%Credit to Ben Popken and Andrew Wickert for portions of the toolbox. +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +%Portions of this toolbox are based on a number of codes written by +%previous authors, including matPIV, IMGRAFT, PIVLAB, M_Map and more. +%Credit and thanks are due to the authors of these toolboxes, and for +%sharing their codes online. See the user manual for a full list of third +%party codes used here. Accordingly, you are free to share, edit and +%add to this GIV code. Please give us credit if you do, and share your code +%with the same conditions as this. + +% Read the associated paper here: +% https://doi.org/10.5194/tc-2020-204 +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %Version 0.7, Autumn 2020% +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %Feel free to contact me at vanwy048@umn.edu% + + + +%% Make sure to add functions to matlab path so that they are accessible +if ~isdeployed +addpath(genpath(pwd)); +end +%Create array to hold all the components in + +GUIinputs = []; + + +%% This portion of the code generated the GUI and input boxes. It is not very interesting. + +% Create UIFigure and hide until all components are created + GUIFigure = uifigure('Visible', 'off'); + GUIFigure.Position = [100 100 656 736]; + GUIFigure.Name = 'GLACIER IMAGE VELOCIMETRY'; + + % Create TabGroup + TabGroup = uitabgroup(GUIFigure); + TabGroup.Position = [1 194 660 543]; + + % Create REQUIREDINPUTSTab + REQUIREDINPUTSTab = uitab(TabGroup); + REQUIREDINPUTSTab.Title = ' REQUIRED INPUTS '; + + % Create PathtoimagesfolderEditFieldLabel + PathtoimagesfolderEditFieldLabel = uilabel(REQUIREDINPUTSTab); + PathtoimagesfolderEditFieldLabel.HorizontalAlignment = 'right'; + PathtoimagesfolderEditFieldLabel.FontName = 'Cambria'; + PathtoimagesfolderEditFieldLabel.Position = [209 305 113 22]; + PathtoimagesfolderEditFieldLabel.Text = 'Path to images folder'; + + % Create PathtoimagesfolderEditField + PathtoimagesfolderEditField = uieditfield(REQUIREDINPUTSTab, 'text'); + PathtoimagesfolderEditField.Position = [337 305 100 22]; + PathtoimagesfolderEditField.Value = 'C:/EXAMPLE/PATH/TO/IMAGES'; + + % Create Image + Image = uiimage(REQUIREDINPUTSTab); + Image.Position = [146 355 347 154]; + Image.ImageSource = 'logo.png'; + + % Create BASICINPUTSREQUIREDLabel + BASICINPUTSREQUIREDLabel = uilabel(REQUIREDINPUTSTab); + BASICINPUTSREQUIREDLabel.FontName = 'Cambria'; + BASICINPUTSREQUIREDLabel.FontSize = 15; + BASICINPUTSREQUIREDLabel.FontWeight = 'bold'; + BASICINPUTSREQUIREDLabel.Position = [236 344 191 22]; + BASICINPUTSREQUIREDLabel.Text = 'BASIC INPUTS (REQUIRED)'; + + + % Create MinimumLatitudeLabel + MinimumLatitudeLabel = uilabel(REQUIREDINPUTSTab); + MinimumLatitudeLabel.HorizontalAlignment = 'right'; + MinimumLatitudeLabel.FontName = 'Cambria'; + MinimumLatitudeLabel.Position = [213 271 109 22]; + MinimumLatitudeLabel.Text = 'Minimum Latitude '; + + % Create MinimumLatitudeEditField + MinimumLatitudeEditField = uieditfield(REQUIREDINPUTSTab, 'numeric'); + MinimumLatitudeEditField.Position = [337 271 100 22]; + MinimumLatitudeEditField.Value = -50.9640; + + % Create MaximumLatitudeEditFieldLabel + MaximumLatitudeEditFieldLabel = uilabel(REQUIREDINPUTSTab); + MaximumLatitudeEditFieldLabel.HorizontalAlignment = 'right'; + MaximumLatitudeEditFieldLabel.FontName = 'Cambria'; + MaximumLatitudeEditFieldLabel.Position = [212 237 110 22]; + MaximumLatitudeEditFieldLabel.Text = 'Maximum Latitude '; + + % Create MaximumLatitudeEditField + MaximumLatitudeEditField = uieditfield(REQUIREDINPUTSTab, 'numeric'); + MaximumLatitudeEditField.Position = [337 237 100 22]; + MaximumLatitudeEditField.Value = -50.8920; + + % Create MinimumLongitudeEditFieldLabel + MinimumLongitudeEditFieldLabel = uilabel(REQUIREDINPUTSTab); + MinimumLongitudeEditFieldLabel.HorizontalAlignment = 'right'; + MinimumLongitudeEditFieldLabel.FontName = 'Cambria'; + MinimumLongitudeEditFieldLabel.Position = [207 205 115 22]; + MinimumLongitudeEditFieldLabel.Text = 'Minimum Longitude '; + + % Create MinimumLongitudeEditField + MinimumLongitudeEditField = uieditfield(REQUIREDINPUTSTab, 'numeric'); + MinimumLongitudeEditField.Position = [337 205 100 22]; + MinimumLongitudeEditField.Value = -73.7460; + + % Create MaximumLongitudeEditFieldLabel + MaximumLongitudeEditFieldLabel = uilabel(REQUIREDINPUTSTab); + MaximumLongitudeEditFieldLabel.HorizontalAlignment = 'right'; + MaximumLongitudeEditFieldLabel.FontName = 'Cambria'; + MaximumLongitudeEditFieldLabel.Position = [208 170 114 22]; + MaximumLongitudeEditFieldLabel.Text = 'Maximum Longitude '; + + % Create MaximumLongitudeEditField + MaximumLongitudeEditField = uieditfield(REQUIREDINPUTSTab, 'numeric'); + MaximumLongitudeEditField.Position = [337 170 100 22]; + MaximumLongitudeEditField.Value = -73.4730; + + + % Create ParralelizeCodeLabel + geotiffinputLabel = uilabel(REQUIREDINPUTSTab); + geotiffinputLabel.FontName = 'Cambria'; + geotiffinputLabel.Position = [515 250 110 22]; + geotiffinputLabel.Text = 'geoTiff input?'; + + % Create geotiffinput + geotiffinputSwitch = uiswitch(REQUIREDINPUTSTab, 'slider'); + geotiffinputSwitch.Items = {'No', 'Yes'}; + geotiffinputSwitch.FontName = 'Cambria'; + geotiffinputSwitch.Position = [525 220 110 22]; + geotiffinputSwitch.Value = 'Yes'; + + + % Create TimeoversamplingfactorEditFieldLabel + TimeoversamplingfactorEditFieldLabel = uilabel(REQUIREDINPUTSTab); + TimeoversamplingfactorEditFieldLabel.HorizontalAlignment = 'right'; + TimeoversamplingfactorEditFieldLabel.FontName = 'Cambria'; + TimeoversamplingfactorEditFieldLabel.Position = [185 134 137 22]; + TimeoversamplingfactorEditFieldLabel.Text = 'Time oversampling factor'; + + % Create TimeoversamplingfactorEditField + TimeoversamplingfactorEditField = uieditfield(REQUIREDINPUTSTab, 'numeric'); + TimeoversamplingfactorEditField.Position = [337 134 100 22]; + TimeoversamplingfactorEditField.Value = 1; + + % Create ParralelizecodeSwitch + ParralelizecodeSwitch = uiswitch(REQUIREDINPUTSTab, 'slider'); + ParralelizecodeSwitch.Items = {'No', 'Yes'}; + ParralelizecodeSwitch.FontName = 'Cambria'; + ParralelizecodeSwitch.Position = [309 15 45 20]; + ParralelizecodeSwitch.Value = 'Yes'; + + + % Create FilenametosaveasEditFieldLabel + FilenametosaveasEditFieldLabel = uilabel(REQUIREDINPUTSTab); + FilenametosaveasEditFieldLabel.HorizontalAlignment = 'right'; + FilenametosaveasEditFieldLabel.FontName = 'Cambria'; + FilenametosaveasEditFieldLabel.Position = [215 99 107 22]; + FilenametosaveasEditFieldLabel.Text = 'File name to save as'; + + % Create FilenametosaveasEditField + FilenametosaveasEditField = uieditfield(REQUIREDINPUTSTab, 'text'); + FilenametosaveasEditField.Position = [337 99 100 22]; + FilenametosaveasEditField.Value = 'Glacier Velocities 1'; + + + + % Create ParralelizeCodeLabel + ParralelizeCodeLabel = uilabel(REQUIREDINPUTSTab); + ParralelizeCodeLabel.Position = [282 48 101 22]; + ParralelizeCodeLabel.FontName = 'Cambria'; + ParralelizeCodeLabel.Text = 'Parralelize Code?'; + + % Create SelectButton + SelectButton = uibutton(REQUIREDINPUTSTab, 'push'); + SelectButton.BackgroundColor = [0.8 0.8 0.8]; + SelectButton.FontName = 'Cambria'; + SelectButton.FontWeight = 'bold'; + SelectButton.Position = [517 305 100 23]; + SelectButton.Text = 'Select'; + SelectButton.ButtonPushedFcn = @(btn,event) plotButtonPushed1(btn,GUIinputs); + + + % Create orLabel + orLabel = uilabel(REQUIREDINPUTSTab); + orLabel.FontName = 'Cambria'; + orLabel.Position = [468 305 25 22]; + orLabel.Text = 'or'; + + % Create AdvancedInputsTab + AdvancedInputsTab = uitab(TabGroup); + AdvancedInputsTab.Title = ' Advanced Inputs '; + + % Create Image_2 + Image_2 = uiimage(AdvancedInputsTab); + Image_2.Position = [146 355 347 154]; + Image_2.ImageSource = 'logo.png'; + + % Create TabGroup2 + TabGroup2 = uitabgroup(AdvancedInputsTab); + TabGroup2.Position = [1 1 659 368]; + + % Create TemplatematchingTab + TemplatematchingTab = uitab(TabGroup2); + TemplatematchingTab.Title = ' Template matching '; + + % Create Switch + Switch = uiswitch(TemplatematchingTab, 'slider'); + Switch.Items = {'Single Pass', 'Multipass'}; + Switch.FontName = 'Cambria'; + Switch.Position = [308 311 45 20]; + Switch.Value = 'Multipass'; + + + + % Create SignaltonoiseratioEditFieldLabel + SignaltonoiseratioEditFieldLabel = uilabel(TemplatematchingTab); + SignaltonoiseratioEditFieldLabel.HorizontalAlignment = 'right'; + SignaltonoiseratioEditFieldLabel.FontName = 'Cambria'; + SignaltonoiseratioEditFieldLabel.Position = [208 215 106 22]; + SignaltonoiseratioEditFieldLabel.Text = 'Signal to noise ratio'; + + % Create SignaltonoiseratioEditField + SignaltonoiseratioEditField = uieditfield(TemplatematchingTab, 'numeric'); + SignaltonoiseratioEditField.Position = [329 215 100 22]; + SignaltonoiseratioEditField.Value = 5; + + + % Create PeakratioEditFieldLabel + PeakratioEditFieldLabel = uilabel(TemplatematchingTab); + PeakratioEditFieldLabel.HorizontalAlignment = 'right'; + PeakratioEditFieldLabel.FontName = 'Cambria'; + PeakratioEditFieldLabel.Position = [208 170 106 22]; + PeakratioEditFieldLabel.Text = 'Peak ratio'; + + % Create PeakratioEditFieldField + PeakratioEditFieldField = uieditfield(TemplatematchingTab, 'numeric'); + PeakratioEditFieldField.Position = [329 170 100 22]; + PeakratioEditFieldField.Value = 1.3; + + + + + +% % Create MultipassoptionsLabel +% MultipassoptionsLabel = uilabel(TemplatematchingTab); +% MultipassoptionsLabel.FontName = 'Cambria'; +% MultipassoptionsLabel.Position = [447 259 97 22]; +% MultipassoptionsLabel.Text = 'Multipass options'; +% +% % Create WindowmatchoverlapEditFieldLabel +% WindowmatchoverlapEditFieldLabel = uilabel(TemplatematchingTab); +% WindowmatchoverlapEditFieldLabel.HorizontalAlignment = 'right'; +% WindowmatchoverlapEditFieldLabel.FontName = 'Cambria'; +% WindowmatchoverlapEditFieldLabel.Position = [387 210 124 22]; +% WindowmatchoverlapEditFieldLabel.Text = 'Window match overlap'; +% +% % Create WindowmatchoverlapEditField +% WindowmatchoverlapEditField = uieditfield(TemplatematchingTab, 'numeric'); +% WindowmatchoverlapEditField.Position = [526 210 100 22]; +% WindowmatchoverlapEditField.Value = 0.5; +% +% % Create MatchwindowsizeEditFieldLabel +% MatchwindowsizeEditFieldLabel = uilabel(TemplatematchingTab); +% MatchwindowsizeEditFieldLabel.HorizontalAlignment = 'right'; +% MatchwindowsizeEditFieldLabel.FontName = 'Cambria'; +% MatchwindowsizeEditFieldLabel.Position = [408 170 103 22]; +% MatchwindowsizeEditFieldLabel.Text = 'Match window size'; +% +% % Create MatchwindowsizeEditField +% MatchwindowsizeEditField = uieditfield(TemplatematchingTab, 'numeric'); +% MatchwindowsizeEditField.FontName = 'Cambria'; +% MatchwindowsizeEditField.Position = [526 170 100 22]; +% MatchwindowsizeEditField.Value = 64; +% +% +% +% +% % Create SinglepassoptionsLabel +% SinglepassoptionsLabel = uilabel(TemplatematchingTab); +% SinglepassoptionsLabel.FontName = 'Cambria'; +% SinglepassoptionsLabel.Position = [89 259 103 22]; +% SinglepassoptionsLabel.Text = 'Single pass options'; + + % Create IdealresolutionofoutputdataEditFieldLabel + IdealresolutionofoutputdataEditFieldLabel = uilabel(TemplatematchingTab); + IdealresolutionofoutputdataEditFieldLabel.HorizontalAlignment = 'right'; + IdealresolutionofoutputdataEditFieldLabel.FontName = 'Cambria'; + IdealresolutionofoutputdataEditFieldLabel.Position = [152 259 161 22]; + IdealresolutionofoutputdataEditFieldLabel.Text = 'Ideal resolution of output data'; + + % Create IdealresolutionofoutputdataEditField + IdealresolutionofoutputdataEditField = uieditfield(TemplatematchingTab, 'numeric'); + IdealresolutionofoutputdataEditField.Position = [328 259 100 22]; + IdealresolutionofoutputdataEditField.Value = 50; + + +% % Create SearchwindowSizeEditFieldLabel +% SearchwindowSizeEditFieldLabel = uilabel(TemplatematchingTab); +% SearchwindowSizeEditFieldLabel.HorizontalAlignment = 'right'; +% SearchwindowSizeEditFieldLabel.FontName = 'Cambria'; +% SearchwindowSizeEditFieldLabel.Position = [77 170 107 22]; +% SearchwindowSizeEditFieldLabel.Text = 'Search window Size'; +% +% % Create SearchwindowSizeEditField +% SearchwindowSizeEditField = uieditfield(TemplatematchingTab, 'numeric'); +% SearchwindowSizeEditField.Position = [199 170 100 22]; +% SearchwindowSizeEditField.Value = 30; +% +% % Create MinimumsearchareaEditFieldLabel +% MinimumsearchareaEditFieldLabel = uilabel(TemplatematchingTab); +% MinimumsearchareaEditFieldLabel.HorizontalAlignment = 'right'; +% MinimumsearchareaEditFieldLabel.FontName = 'Cambria'; +% MinimumsearchareaEditFieldLabel.Position = [68 132 116 22]; +% MinimumsearchareaEditFieldLabel.Text = 'Minimum search area'; +% +% % Create MinimumsearchareaEditField +% MinimumsearchareaEditField = uieditfield(TemplatematchingTab, 'numeric'); +% MinimumsearchareaEditField.Position = [199 132 100 22]; +% MinimumsearchareaEditField.Value = 50; + + + % Create DateoptionsTab + DateoptionsTab = uitab(TabGroup2); + DateoptionsTab.Title = ' Date options '; + + % Create MinimumYearEditFieldLabel + MinimumYearEditFieldLabel = uilabel(DateoptionsTab); + MinimumYearEditFieldLabel.HorizontalAlignment = 'right'; + MinimumYearEditFieldLabel.FontName = 'Cambria'; + MinimumYearEditFieldLabel.Position = [245 288 85 22]; + MinimumYearEditFieldLabel.Text = 'Minimum Year '; + + % Create MinimumYearEditField + MinimumYearEditField = uieditfield(DateoptionsTab, 'numeric'); + MinimumYearEditField.Position = [345 288 100 22]; + MinimumYearEditField.Value = 1900; + + % Create MaximumYearEditFieldLabel + MaximumYearEditFieldLabel = uilabel(DateoptionsTab); + MaximumYearEditFieldLabel.HorizontalAlignment = 'right'; + MaximumYearEditFieldLabel.FontName = 'Cambria'; + MaximumYearEditFieldLabel.Position = [243 254 87 22]; + MaximumYearEditFieldLabel.Text = 'Maximum Year '; + + % Create MaximumYearEditField + MaximumYearEditField = uieditfield(DateoptionsTab, 'numeric'); + MaximumYearEditField.Position = [345 254 100 22]; + MaximumYearEditField.Value = 2050; + + + % Create MinimumMonthEditFieldLabel + MinimumMonthEditFieldLabel = uilabel(DateoptionsTab); + MinimumMonthEditFieldLabel.HorizontalAlignment = 'right'; + MinimumMonthEditFieldLabel.FontName = 'Cambria'; + MinimumMonthEditFieldLabel.Position = [239 222 91 22]; + MinimumMonthEditFieldLabel.Text = 'Minimum Month'; + + % Create MinimumMonthEditField + MinimumMonthEditField = uieditfield(DateoptionsTab, 'numeric'); + MinimumMonthEditField.Position = [345 222 100 22]; + MinimumMonthEditField.Value = 1; + + % Create MaximumMonthEditFieldLabel + MaximumMonthEditFieldLabel = uilabel(DateoptionsTab); + MaximumMonthEditFieldLabel.HorizontalAlignment = 'right'; + MaximumMonthEditFieldLabel.FontName = 'Cambria'; + MaximumMonthEditFieldLabel.Position = [237 187 93 22]; + MaximumMonthEditFieldLabel.Text = 'Maximum Month'; + + % Create MaximumMonthEditField + MaximumMonthEditField = uieditfield(DateoptionsTab, 'numeric'); + MaximumMonthEditField.Position = [345 187 100 22]; + MaximumMonthEditField.Value = 12; + + + + % Create MinimumDayEditFieldLabel + MinimumDayEditFieldLabel = uilabel(DateoptionsTab); + MinimumDayEditFieldLabel.HorizontalAlignment = 'right'; + MinimumDayEditFieldLabel.FontName = 'Cambria'; + MinimumDayEditFieldLabel.Position = [253 152 77 22]; + MinimumDayEditFieldLabel.Text = 'Minimum Day'; + + % Create MinimumDayEditField + MinimumDayEditField = uieditfield(DateoptionsTab, 'numeric'); + MinimumDayEditField.Position = [345 152 100 22]; + MinimumDayEditField.Value = 1; + + % Create MaximumDayEditFieldLabel + MaximumDayEditFieldLabel = uilabel(DateoptionsTab); + MaximumDayEditFieldLabel.HorizontalAlignment = 'right'; + MaximumDayEditFieldLabel.FontName = 'Cambria'; + MaximumDayEditFieldLabel.Position = [251 117 79 22]; + MaximumDayEditFieldLabel.Text = 'Maximum Day'; + + % Create MaximumDayEditField + MaximumDayEditField = uieditfield(DateoptionsTab, 'numeric'); + MaximumDayEditField.Position = [345 117 100 22]; + MaximumDayEditField.Value = 31; + + + % Create MinimumintervalforimagepairsEditFieldLabel + MinimumintervalforimagepairsEditFieldLabel = uilabel(DateoptionsTab); + MinimumintervalforimagepairsEditFieldLabel.HorizontalAlignment = 'right'; + MinimumintervalforimagepairsEditFieldLabel.FontName = 'Cambria'; + MinimumintervalforimagepairsEditFieldLabel.Position = [153 70 177 22]; + MinimumintervalforimagepairsEditFieldLabel.Text = 'Minimum interval for image pairs'; + + % Create MinimumintervalforimagepairsEditField + MinimumintervalforimagepairsEditField = uieditfield(DateoptionsTab, 'numeric'); + MinimumintervalforimagepairsEditField.Position = [345 70 100 22]; + MinimumintervalforimagepairsEditField.Value = 0.019; + + + % Create MaximumintervalforimagepairsEditFieldLabel + MaximumintervalforimagepairsEditFieldLabel = uilabel(DateoptionsTab); + MaximumintervalforimagepairsEditFieldLabel.HorizontalAlignment = 'right'; + MaximumintervalforimagepairsEditFieldLabel.FontName = 'Cambria'; + MaximumintervalforimagepairsEditFieldLabel.Position = [151 35 179 22]; + MaximumintervalforimagepairsEditFieldLabel.Text = 'Maximum interval for image pairs'; + + % Create MaximumintervalforimagepairsEditField + MaximumintervalforimagepairsEditField = uieditfield(DateoptionsTab, 'numeric'); + MaximumintervalforimagepairsEditField.FontName = 'Cambria'; + MaximumintervalforimagepairsEditField.Position = [345 35 100 22]; + MaximumintervalforimagepairsEditField.Value = 0.75; + + + % Create ImagefilteringTab + ImagefilteringTab = uitab(TabGroup2); + ImagefilteringTab.Title = ' Image filtering '; + + % Create Switch_2 + Switch_2 = uiswitch(ImagefilteringTab, 'slider'); + Switch_2.Items = {'OFF', 'ON'}; + Switch_2.FontName = 'Cambria'; + Switch_2.Position = [200 293 45 20]; + Switch_2.Value = 'OFF'; + + % Create ContrastLimitedHistogramEqualisationLabel + ContrastLimitedHistogramEqualisationLabel = uilabel(ImagefilteringTab); + ContrastLimitedHistogramEqualisationLabel.Position = [123 319 220 22]; + ContrastLimitedHistogramEqualisationLabel.Text = 'Contrast Limited Histogram Equalisation'; + + % Create Switch_3 + Switch_3 = uiswitch(ImagefilteringTab, 'slider'); + Switch_3.Items = {'OFF', 'ON'}; + Switch_3.FontName = 'Cambria'; + Switch_3.Position = [201 240 45 20]; + Switch_3.Value = 'OFF'; + + % Create HighpassFilterLabel + HighpassFilterLabel = uilabel(ImagefilteringTab); + HighpassFilterLabel.Position = [180 266 85 22]; + HighpassFilterLabel.Text = 'Highpass Filter'; + + + % Create Switch_4 + Switch_4 = uiswitch(ImagefilteringTab, 'slider'); + Switch_4.Items = {'OFF', 'ON'}; + Switch_4.FontName = 'Cambria'; + Switch_4.Position = [202 191 45 20]; + Switch_4.Value = 'OFF'; + + % Create IntensityCapLabel + IntensityCapLabel = uilabel(ImagefilteringTab); + IntensityCapLabel.Position = [187 215 75 22]; + IntensityCapLabel.Text = 'Intensity Cap'; + + + % Create Switch_5 + Switch_5 = uiswitch(ImagefilteringTab, 'slider'); + Switch_5.Items = {'OFF', 'ON'}; + Switch_5.FontName = 'Cambria'; + Switch_5.Position = [201 130 45 20]; + Switch_5.Value = 'ON'; + + % Create OrientationFilterLabel + OrientationFilterLabel = uilabel(ImagefilteringTab); + OrientationFilterLabel.Position = [178 154 94 22]; + OrientationFilterLabel.Text = 'Orientation Filter'; + + + + % Create CLAHEsizeEditFieldLabel + CLAHEsizeEditFieldLabel = uilabel(ImagefilteringTab); + CLAHEsizeEditFieldLabel.HorizontalAlignment = 'right'; + CLAHEsizeEditFieldLabel.Position = [384 298 70 22]; + CLAHEsizeEditFieldLabel.Text = 'CLAHE size'; + + % Create CLAHEsizeEditField + CLAHEsizeEditField = uieditfield(ImagefilteringTab, 'numeric'); + CLAHEsizeEditField.Position = [469 298 100 22]; + CLAHEsizeEditField.Value = 10; + + + + + % Create HighpasssizeEditFieldLabel + HighpasssizeEditFieldLabel = uilabel(ImagefilteringTab); + HighpasssizeEditFieldLabel.HorizontalAlignment = 'right'; + HighpasssizeEditFieldLabel.Position = [375 245 80 22]; + HighpasssizeEditFieldLabel.Text = 'Highpass size'; + + % Create HighpasssizeEditField + HighpasssizeEditField = uieditfield(ImagefilteringTab, 'numeric'); + HighpasssizeEditField.Position = [470 245 100 22]; + HighpasssizeEditField.Value = 10; + + + + % Create Switch_10 + Switch_10 = uiswitch(ImagefilteringTab, 'slider'); + Switch_10.Items = {'OFF', 'ON'}; + Switch_10.FontName = 'Cambria'; + Switch_10.Position = [202 70 45 20]; + Switch_10.Value = 'OFF'; + + % Create Sobel Filter + Sobel_Filter = uilabel(ImagefilteringTab); + Sobel_Filter.Position = [187 94 75 22]; + Sobel_Filter.Text = 'Sobel Filter'; + + + % Create Switch_11 + Switch_11 = uiswitch(ImagefilteringTab, 'slider'); + Switch_11.Items = {'OFF', 'ON'}; + Switch_11.FontName = 'Cambria'; + Switch_11.Position = [201 9 45 20]; + Switch_11.Value = 'OFF'; + + % Create OrientationFilterLabel_2 + Laplacian_filter = uilabel(ImagefilteringTab); + Laplacian_filter.Position = [178 33 94 22]; + Laplacian_filter.Text = 'Laplacian filter'; + + + % Create SavingTab + SavingTab = uitab(TabGroup2); + SavingTab.Title = ' Saving '; + + % Create Switch_6 + Switch_6 = uiswitch(SavingTab, 'slider'); + Switch_6.Items = {'OFF', 'ON'}; + Switch_6.FontName = 'Cambria'; + Switch_6.Position = [298 275 45 20]; + Switch_6.Value = 'ON'; + + % Create SaveMATLABdataarraysLabel + SaveMATLABdataarraysLabel = uilabel(SavingTab); + SaveMATLABdataarraysLabel.Position = [254 302 153 22]; + SaveMATLABdataarraysLabel.Text = 'Save MATLAB data arrays?'; + + + % Create Switch_7 + Switch_7 = uiswitch(SavingTab, 'slider'); + Switch_7.Items = {'OFF', 'ON'}; + Switch_7.FontName = 'Cambria'; + Switch_7.Position = [298 198 45 20]; + Switch_7.Value = 'ON'; + + % Create Switch_8 + Switch_8 = uiswitch(SavingTab, 'slider'); + Switch_8.Items = {'OFF', 'ON'}; + Switch_8.FontName = 'Cambria'; + Switch_8.Position = [299 19 45 20]; + Switch_8.Value = 'ON'; + + % Create SavegeoreferencedtifimagesLabel + SavegeoreferencedtifimagesLabel = uilabel(SavingTab); + SavegeoreferencedtifimagesLabel.Position = [233 46 178 22]; + SavegeoreferencedtifimagesLabel.Text = 'Save georeferenced .tif images?'; + + % Create SaveimagesofkeyvelocitiesLabel + SaveimagesofkeyvelocitiesLabel = uilabel(SavingTab); + SaveimagesofkeyvelocitiesLabel.Position = [237 229 169 22]; + SaveimagesofkeyvelocitiesLabel.Text = 'Save images of key velocities?'; + + + % Create FormatofimagestosaveEditFieldLabel + FormatofimagestosaveEditFieldLabel = uilabel(SavingTab); + FormatofimagestosaveEditFieldLabel.HorizontalAlignment = 'right'; + FormatofimagestosaveEditFieldLabel.Position = [147 119 141 22]; + FormatofimagestosaveEditFieldLabel.Text = 'Format of images to save'; + + % Create FormatofimagestosaveEditField + FormatofimagestosaveEditField = uieditfield(SavingTab, 'text'); + FormatofimagestosaveEditField.Position = [303 119 100 22]; + FormatofimagestosaveEditField.Value = 'png'; + + + % Create OtherTab + OtherTab = uitab(TabGroup2); + OtherTab.Title = ' Other '; + + % Create MaximumVelocityEditFieldLabel + MaximumVelocityEditFieldLabel = uilabel(OtherTab); + MaximumVelocityEditFieldLabel.HorizontalAlignment = 'right'; + MaximumVelocityEditFieldLabel.FontName = 'Cambria'; + MaximumVelocityEditFieldLabel.Position = [222 300 100 22]; + MaximumVelocityEditFieldLabel.Text = 'Maximum Velocity'; + + % Create MaximumVelocityEditField + MaximumVelocityEditField = uieditfield(OtherTab, 'numeric'); + MaximumVelocityEditField.Position = [337 300 100 22]; + MaximumVelocityEditField.Value = 2500; + + % Create Excludedangle1minimumEditFieldLabel + Excludedangle1minimumEditFieldLabel = uilabel(OtherTab); + Excludedangle1minimumEditFieldLabel.HorizontalAlignment = 'right'; + Excludedangle1minimumEditFieldLabel.FontName = 'Cambria'; + Excludedangle1minimumEditFieldLabel.Position = [58 210 144 22]; + Excludedangle1minimumEditFieldLabel.Text = 'Excluded angle 1 minimum'; + + % Create Excludedangle1minimumEditField + Excludedangle1minimumEditField = uieditfield(OtherTab, 'numeric'); + Excludedangle1minimumEditField.Position = [217 210 100 22]; + + % Create Excludedangle2minimumEditFieldLabel + Excludedangle2minimumEditFieldLabel = uilabel(OtherTab); + Excludedangle2minimumEditFieldLabel.HorizontalAlignment = 'right'; + Excludedangle2minimumEditFieldLabel.FontName = 'Cambria'; + Excludedangle2minimumEditFieldLabel.Position = [58 174 144 22]; + Excludedangle2minimumEditFieldLabel.Text = 'Excluded angle 2 minimum'; + + % Create Excludedangle2minimumEditField + Excludedangle2minimumEditField = uieditfield(OtherTab, 'numeric'); + Excludedangle2minimumEditField.Position = [217 174 100 22]; + Excludedangle2minimumEditField.Value = 360; + + % Create Excludedangle2maximumEditFieldLabel + Excludedangle2maximumEditFieldLabel = uilabel(OtherTab); + Excludedangle2maximumEditFieldLabel.HorizontalAlignment = 'right'; + Excludedangle2maximumEditFieldLabel.FontName = 'Cambria'; + Excludedangle2maximumEditFieldLabel.Position = [353 174 146 22]; + Excludedangle2maximumEditFieldLabel.Text = 'Excluded angle 2 maximum'; + + % Create Excludedangle2maximumEditField + Excludedangle2maximumEditField = uieditfield(OtherTab, 'numeric'); + Excludedangle2maximumEditField.Position = [514 174 100 22]; + Excludedangle2maximumEditField.Value = 360; + + % Create Excludedangle1maximumEditFieldLabel + Excludedangle1maximumEditFieldLabel = uilabel(OtherTab); + Excludedangle1maximumEditFieldLabel.HorizontalAlignment = 'right'; + Excludedangle1maximumEditFieldLabel.FontName = 'Cambria'; + Excludedangle1maximumEditFieldLabel.Position = [353 210 146 22]; + Excludedangle1maximumEditFieldLabel.Text = 'Excluded angle 1 maximum'; + + % Create Excludedangle1maximumEditField + Excludedangle1maximumEditField = uieditfield(OtherTab, 'numeric'); + Excludedangle1maximumEditField.Position = [514 210 100 22]; + + +% % % % Create MinimumangleEditFieldLabel +% % % MinimumangleEditFieldLabel = uilabel(OtherTab); +% % % MinimumangleEditFieldLabel.HorizontalAlignment = 'right'; +% % % MinimumangleEditFieldLabel.FontName = 'Cambria'; +% % % MinimumangleEditFieldLabel.Position = [230 187 85 22]; +% % % MinimumangleEditFieldLabel.Text = 'Minimum angle'; +% % % +% % % % Create MinimumangleEditField +% % % MinimumangleEditField = uieditfield(OtherTab, 'numeric'); +% % % MinimumangleEditField.Position = [330 187 100 22]; +% % % MinimumangleEditField.Value = 0; +% % % +% % % +% % % % Create MaximumAngleEditFieldLabel +% % % MaximumAngleEditFieldLabel = uilabel(OtherTab); +% % % MaximumAngleEditFieldLabel.HorizontalAlignment = 'right'; +% % % MaximumAngleEditFieldLabel.FontName = 'Cambria'; +% % % MaximumAngleEditFieldLabel.Position = [227 151 88 22]; +% % % MaximumAngleEditFieldLabel.Text = 'Maximum Angle'; +% % % +% % % % Create MaximumAngleEditField +% % % MaximumAngleEditField = uieditfield(OtherTab, 'numeric'); +% % % MaximumAngleEditField.Position = [330 151 100 22]; +% % % MaximumAngleEditField.Value = 360; + + % Create Switch_9 + Switch_9 = uiswitch(OtherTab, 'slider'); + Switch_9.Items = {'NO', 'YES'}; + Switch_9.FontName = 'Cambria'; + Switch_9.Position = [316 243 45 20]; + Switch_9.Value = 'NO'; + + + % Create FilterbasedonflowdirectionLabel + FilterbasedonflowdirectionLabel = uilabel(OtherTab); + FilterbasedonflowdirectionLabel.Position = [272 270 158 22]; + FilterbasedonflowdirectionLabel.Text = 'Filter based on flow direction'; + + + + % Create SmoothingofcompositearrayofallvelocitiesDropDownLabel + SmoothingofcompositearrayofallvelocitiesDropDownLabel = uilabel(OtherTab); + SmoothingofcompositearrayofallvelocitiesDropDownLabel.HorizontalAlignment = 'right'; + SmoothingofcompositearrayofallvelocitiesDropDownLabel.Position = [103 18 248 22]; + SmoothingofcompositearrayofallvelocitiesDropDownLabel.Text = 'Smoothing of composite array of all velocities'; + + % Create SmoothingofcompositearrayofallvelocitiesDropDown + SmoothingofcompositearrayofallvelocitiesDropDown = uidropdown(OtherTab); + SmoothingofcompositearrayofallvelocitiesDropDown.Items = {'No Smoothing', 'Smoothing in time', 'Smoothing in time and space'}; + SmoothingofcompositearrayofallvelocitiesDropDown.Position = [365 18 192 22]; + SmoothingofcompositearrayofallvelocitiesDropDown.Value = 'Smoothing in time and space'; + + % Create NumberofiterationsformonthlyvelocitiesEditFieldLabel + NumberofiterationsformonthlyvelocitiesEditFieldLabel = uilabel(OtherTab); + NumberofiterationsformonthlyvelocitiesEditFieldLabel.HorizontalAlignment = 'right'; + NumberofiterationsformonthlyvelocitiesEditFieldLabel.Position = [134 59 231 22]; + NumberofiterationsformonthlyvelocitiesEditFieldLabel.Text = 'Number of iterations for monthly velocities'; + + % Create NumberofiterationsformonthlyvelocitiesEditField + NumberofiterationsformonthlyvelocitiesEditField = uieditfield(OtherTab, 'numeric'); + NumberofiterationsformonthlyvelocitiesEditField.Position = [380 59 100 22]; + NumberofiterationsformonthlyvelocitiesEditField.Value = 0; + + % Create Switch_12 + Switch_12 = uiswitch(OtherTab, 'slider'); + Switch_12.Items = {'NO', 'YES'}; + Switch_12.FontName = 'Cambria'; + Switch_12.Position = [321 102 45 20]; + Switch_12.Value = 'NO'; + + % Create NormalizetovelocityofastableregionLabel + NormalizetovelocityofastableregionLabel = uilabel(OtherTab); + NormalizetovelocityofastableregionLabel.Position = [245 130 212 22]; + NormalizetovelocityofastableregionLabel.Text = 'Normalize to velocity of a stable region'; + + + % Create GIVAGLACIERVELOCITYCALCULATIONTOOLBOXBYMAXVANWYKDEVRIESLabel + GIVAGLACIERVELOCITYCALCULATIONTOOLBOXBYMAXVANWYKDEVRIESLabel = uilabel(GUIFigure); + GIVAGLACIERVELOCITYCALCULATIONTOOLBOXBYMAXVANWYKDEVRIESLabel.FontName = 'Cambria'; + GIVAGLACIERVELOCITYCALCULATIONTOOLBOXBYMAXVANWYKDEVRIESLabel.FontAngle = 'italic'; + GIVAGLACIERVELOCITYCALCULATIONTOOLBOXBYMAXVANWYKDEVRIESLabel.Position = [116 26 447 22]; + GIVAGLACIERVELOCITYCALCULATIONTOOLBOXBYMAXVANWYKDEVRIESLabel.Text = 'GIV: A GLACIER VELOCITY TOOLBOX BY VAN WYK DE VRIES AND WICKERT'; + + % Create WWWGIVGLACIERCOMCONTACTMEATVANWY048UMNEDULabel + WWWGIVGLACIERCOMCONTACTMEATVANWY048UMNEDULabel = uilabel(GUIFigure); + WWWGIVGLACIERCOMCONTACTMEATVANWY048UMNEDULabel.FontName = 'Cambria'; + WWWGIVGLACIERCOMCONTACTMEATVANWY048UMNEDULabel.FontSize = 10; + WWWGIVGLACIERCOMCONTACTMEATVANWY048UMNEDULabel.FontAngle = 'italic'; + WWWGIVGLACIERCOMCONTACTMEATVANWY048UMNEDULabel.Position = [184 6 295 22]; + WWWGIVGLACIERCOMCONTACTMEATVANWY048UMNEDULabel.Text = ' CONTACT ME AT VANWY048@UMN.EDU '; + + % Create CALCULATEVELOCITIESButton + CALCULATEVELOCITIESButton = uibutton(GUIFigure, 'push'); + CALCULATEVELOCITIESButton.BackgroundColor = [0.302 0.7451 0.9333]; + CALCULATEVELOCITIESButton.FontName = 'Cambria'; + CALCULATEVELOCITIESButton.FontSize = 25; + CALCULATEVELOCITIESButton.FontWeight = 'bold'; + CALCULATEVELOCITIESButton.FontAngle = 'italic'; + CALCULATEVELOCITIESButton.FontColor = [0 0.149 1]; + CALCULATEVELOCITIESButton.Position = [186 59 310 71]; + CALCULATEVELOCITIESButton.Text = 'CALCULATE VELOCITIES'; + CALCULATEVELOCITIESButton.ButtonPushedFcn = @(btn,event) plotButtonPushed2 (btn); + + % Create LoadsetupButton + LoadsetupButton = uibutton(GUIFigure, 'push'); + LoadsetupButton.FontName = 'Cambria'; + LoadsetupButton.FontWeight = 'bold'; + LoadsetupButton.Position = [26 77 100 34]; + LoadsetupButton.Text = 'Load setup'; + LoadsetupButton.ButtonPushedFcn = @(btn,event) plotButtonPushed3(btn); + + % Create SavesetupButton + SavesetupButton = uibutton(GUIFigure, 'push'); + SavesetupButton.BackgroundColor = [0.9412 0.9412 0.9412]; + SavesetupButton.FontName = 'Cambria'; + SavesetupButton.FontWeight = 'bold'; + SavesetupButton.Position = [536 77 100 34]; + SavesetupButton.Text = 'Save setup'; + SavesetupButton.ButtonPushedFcn = @(btn,event) plotButtonPushed4 (btn); + + % Create AnalyseImagePairsButton + AnalyseImagePairsButton = uibutton(GUIFigure, 'push'); + AnalyseImagePairsButton.BackgroundColor = [0.302 0.7451 0.9333]; + AnalyseImagePairsButton.FontName = 'Cambria'; + AnalyseImagePairsButton.FontSize = 18; + AnalyseImagePairsButton.FontWeight = 'bold'; + AnalyseImagePairsButton.FontAngle = 'italic'; + AnalyseImagePairsButton.Position = [248 142 185 40]; + AnalyseImagePairsButton.Text = 'Analyse Image Pairs'; + AnalyseImagePairsButton.ButtonPushedFcn = @(btn,event) plotButtonPushed5 (btn); + + + + % Show the figure after all components are created + GUIFigure.Visible = 'on'; + + + %% This portion makes the buttons call the relevant functions. + + %For ease, functions are all wrapped in a single GIV_GUI_main, + %which is called here. + + + % This button press allows you to select the input folder. + function [btn,event,GUIinputs]= plotButtonPushed1(btn,GUIinputs) + a = uigetdir('.'); + PathtoimagesfolderEditField.Value = a; + GUIinputs.folder = a; + figure(GUIFigure) + end + + %This button runs GIV + function [btn,event]= plotButtonPushed2(btn,GUIinputs) + + %Display message to user + logo = imread('GIV_LOGO_SMALL.png'); + msgbox({'RUNNING...YOU MAY CLOSE THE INPUTS BOX. IT MAY TAKE A FEW SECONDS TO CLOSE'},... + 'GIV is running','custom',logo); + + %Combine input into a MATLAB struct array. These are more human + %readable than a simple array, and easier to modify in the + %future. + GUIinputs.folder = PathtoimagesfolderEditField.Value; + GUIinputs.isgeotiff = geotiffinputSwitch.Value; + GUIinputs.minlat = MinimumLatitudeEditField.Value; + GUIinputs.maxlat = MaximumLatitudeEditField.Value; + GUIinputs.minlon = MinimumLongitudeEditField.Value; + GUIinputs.maxlon = MaximumLongitudeEditField.Value; + GUIinputs.temporaloversampling = TimeoversamplingfactorEditField.Value; + GUIinputs.parralelize = ParralelizecodeSwitch.Value; + GUIinputs.name = FilenametosaveasEditField.Value; + + if strcmpi(Switch.Value, 'Multipass') + GUIinputs.numpass = 'Multi'; + else + GUIinputs.numpass = 'Single'; + end + + GUIinputs.snr = SignaltonoiseratioEditField.Value; + GUIinputs.pkr = PeakratioEditFieldField.Value; + GUIinputs.windowoverlap = 0.5; + GUIinputs.idealresolution = IdealresolutionofoutputdataEditField.Value; + GUIinputs.searchwindowsize = 30; + GUIinputs.minsearcharea = 50; + GUIinputs.minyear = MinimumYearEditField.Value; + GUIinputs.maxyear = MaximumYearEditField.Value; + GUIinputs.minmonth = MinimumMonthEditField.Value; + GUIinputs.maxmonth = MaximumMonthEditField.Value; + GUIinputs.minday = MinimumDayEditField.Value; + GUIinputs.maxday = MaximumDayEditField.Value; + GUIinputs.mininterval = MinimumintervalforimagepairsEditField.Value; + GUIinputs.maxinterval = MaximumintervalforimagepairsEditField.Value; + + if strcmpi(Switch_2.Value, 'ON') + GUIinputs.CLAHE = 1; + else + GUIinputs.CLAHE = 0; + end + + if strcmpi(Switch_3.Value, 'ON') + GUIinputs.hipass = 1; + else + GUIinputs.hipass = 0; + end + + if strcmpi(Switch_4.Value, 'ON') + GUIinputs.intenscap = 1; + else + GUIinputs.intenscap = 0; + end + + if strcmpi(Switch_5.Value, 'ON') + GUIinputs.NAOF = 1; + else + GUIinputs.NAOF = 0; + end + + GUIinputs.CLAHEsize = CLAHEsizeEditField.Value; + GUIinputs.hipasssize = HighpasssizeEditField.Value; + + if strcmpi(Switch_10.Value, 'ON') + GUIinputs.sobel = 1; + else + GUIinputs.sobel = 0; + end + + if strcmpi(Switch_11.Value, 'ON') + GUIinputs.laplacian = 1; + else + GUIinputs.laplacian = 0; + end + + if strcmpi(Switch_6.Value, 'ON') + GUIinputs.savearrays = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savearrays = 'No'; + end + + if strcmpi(Switch_7.Value, 'ON') + GUIinputs.savekeyvel = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savekeyvel = 'No'; + end + + if strcmpi(Switch_8.Value, 'ON') + GUIinputs.savegeotiff = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savegeotiff = 'No'; + end + + GUIinputs.imageformat = FormatofimagestosaveEditField.Value; + GUIinputs.maxvel = MaximumVelocityEditField.Value; + GUIinputs.excudedangle1.min = Excludedangle1minimumEditField.Value; + GUIinputs.excudedangle1.max = Excludedangle1maximumEditField.Value; + GUIinputs.excudedangle2.min = Excludedangle2minimumEditField.Value; + GUIinputs.excudedangle2.max = Excludedangle2maximumEditField.Value; + + if strcmpi(Switch_12.Value, 'YES') + GUIinputs.stable = 'Yes'; %'Yes' or 'No' + else + GUIinputs.stable = 'No'; + end + + if strcmpi(Switch_9.Value, 'YES') + GUIinputs.excludeangle = 'Yes'; %'Yes' or 'No' + else + GUIinputs.excludeangle = 'No'; + end + + if strcmpi(SmoothingofcompositearrayofallvelocitiesDropDown.Value, 'Smoothing in time and space') + GUIinputs.finalsmooth = 'Time and Space'; + elseif strcmpi(SmoothingofcompositearrayofallvelocitiesDropDown.Value, 'Smoothing in time') + GUIinputs.finalsmooth = 'Time'; + else + GUIinputs.finalsmooth = 'None'; + end + + GUIinputs.nummonthiter = NumberofiterationsformonthlyvelocitiesEditField.Value; + + %Rename this array + inputs = GUIinputs; + + %Call the main wrapper function + GIV_GUI_main(inputs) + end + + %Load a previously saved inputs file. + function [btn,event]= plotButtonPushed3(btn) + + %Prompt user to select file + [fname,load_file] = uigetfile('input files', 'Select an input file you previously saved.'); + load(strcat(load_file,'/',fname),'inputs'); + + %Display message + logo = imread('GIV_LOGO_SMALL.png'); + msgbox({'RUNNING...YOU MAY CLOSE THE INPUTS BOX. IT MAY TAKE A FEW SECONDS TO CLOSE'},... + 'GIV is running','custom',logo); + + %Call main wrapper function + GIV_GUI_main(inputs) + end + + %Save an input setup + function [btn,event]= plotButtonPushed4(btn,GUIinputs) + + %Combine input into a MATLAB struct array. These are more human + %readable than a simple array, and easier to modify in the + %future. + GUIinputs.folder = PathtoimagesfolderEditField.Value; + GUIinputs.isgeotiff = geotiffinputSwitch.Value; + GUIinputs.minlat = MinimumLatitudeEditField.Value; + GUIinputs.maxlat = MaximumLatitudeEditField.Value; + GUIinputs.minlon = MinimumLongitudeEditField.Value; + GUIinputs.maxlon = MaximumLongitudeEditField.Value; + GUIinputs.temporaloversampling = TimeoversamplingfactorEditField.Value; + GUIinputs.parralelize = ParralelizecodeSwitch.Value; + GUIinputs.name = FilenametosaveasEditField.Value; + + if strcmpi(Switch.Value, 'Multipass') + GUIinputs.numpass = 'Multi'; + else + GUIinputs.numpass = 'Single'; + end + + GUIinputs.snr = SignaltonoiseratioEditField.Value; + GUIinputs.pkr = PeakratioEditFieldField.Value; + GUIinputs.windowoverlap = 0.5; + GUIinputs.idealresolution = IdealresolutionofoutputdataEditField.Value; + GUIinputs.searchwindowsize = 30; + GUIinputs.minsearcharea = 50; + GUIinputs.minyear = MinimumYearEditField.Value; + GUIinputs.maxyear = MaximumYearEditField.Value; + GUIinputs.minmonth = MinimumMonthEditField.Value; + GUIinputs.maxmonth = MaximumMonthEditField.Value; + GUIinputs.minday = MinimumDayEditField.Value; + GUIinputs.maxday = MaximumDayEditField.Value; + GUIinputs.mininterval = MinimumintervalforimagepairsEditField.Value; + GUIinputs.maxinterval = MaximumintervalforimagepairsEditField.Value; + + if strcmpi(Switch_2.Value, 'ON') + GUIinputs.CLAHE = 1; + else + GUIinputs.CLAHE = 0; + end + + if strcmpi(Switch_3.Value, 'ON') + GUIinputs.hipass = 1; + else + GUIinputs.hipass = 0; + end + + if strcmpi(Switch_4.Value, 'ON') + GUIinputs.intenscap = 1; + else + GUIinputs.intenscap = 0; + end + + if strcmpi(Switch_5.Value, 'ON') + GUIinputs.NAOF = 1; + else + GUIinputs.NAOF = 0; + end + + GUIinputs.CLAHEsize = CLAHEsizeEditField.Value; + GUIinputs.hipasssize = HighpasssizeEditField.Value; + + if strcmpi(Switch_10.Value, 'ON') + GUIinputs.sobel = 1; + else + GUIinputs.sobel = 0; + end + + if strcmpi(Switch_11.Value, 'ON') + GUIinputs.laplacian = 1; + else + GUIinputs.laplacian = 0; + end + + if strcmpi(Switch_6.Value, 'ON') + GUIinputs.savearrays = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savearrays = 'No'; + end + + if strcmpi(Switch_7.Value, 'ON') + GUIinputs.savekeyvel = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savekeyvel = 'No'; + end + + if strcmpi(Switch_8.Value, 'ON') + GUIinputs.savegeotiff = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savegeotiff = 'No'; + end + + GUIinputs.imageformat = FormatofimagestosaveEditField.Value; + GUIinputs.maxvel = MaximumVelocityEditField.Value; + GUIinputs.excudedangle1.min = Excludedangle1minimumEditField.Value; + GUIinputs.excudedangle1.max = Excludedangle1maximumEditField.Value; + GUIinputs.excudedangle2.min = Excludedangle2minimumEditField.Value; + GUIinputs.excudedangle2.max = Excludedangle2maximumEditField.Value; + + if strcmpi(Switch_12.Value, 'YES') + GUIinputs.stable = 'Yes'; %'Yes' or 'No' + else + GUIinputs.stable = 'No'; + end + + if strcmpi(Switch_9.Value, 'YES') + GUIinputs.excludeangle = 'Yes'; %'Yes' or 'No' + else + GUIinputs.excludeangle = 'No'; + end + + if strcmpi(SmoothingofcompositearrayofallvelocitiesDropDown.Value, 'Smoothing in time and space') + GUIinputs.finalsmooth = 'Time and Space'; + elseif strcmpi(SmoothingofcompositearrayofallvelocitiesDropDown.Value, 'Smoothing in time') + GUIinputs.finalsmooth = 'Time'; + else + GUIinputs.finalsmooth = 'None'; + end + + GUIinputs.nummonthiter = NumberofiterationsformonthlyvelocitiesEditField.Value; + + %Rename this array + inputs = GUIinputs; + + %Save the file + filename2 = strcat(inputs.folder,'/input files'); + if ~exist(filename2) + mkdir(filename2) + end + uisave('inputs',strcat(filename2,'/','savefile')); + + end + + %Calculate number of viable images for a given setup. + function [btn,event]= plotButtonPushed5(btn,GUIinputs) + + %Display message to user + logo = imread('GIV_LOGO_SMALL.png'); + msgbox({'CALCULATING NUMBER OF IMAGE PAIRS...YOU MAY CLOSE THE INPUTS BOX. IT MAY TAKE A FEW SECONDS TO CLOSE'},... + 'GIV is running','custom',logo); + + %Combine input into a MATLAB struct array. These are more human + %readable than a simple array, and easier to modify in the + %future. + GUIinputs.folder = PathtoimagesfolderEditField.Value; + GUIinputs.isgeotiff = geotiffinputSwitch.Value; + GUIinputs.minlat = MinimumLatitudeEditField.Value; + GUIinputs.maxlat = MaximumLatitudeEditField.Value; + GUIinputs.minlon = MinimumLongitudeEditField.Value; + GUIinputs.maxlon = MaximumLongitudeEditField.Value; + GUIinputs.temporaloversampling = TimeoversamplingfactorEditField.Value; + GUIinputs.parralelize = ParralelizecodeSwitch.Value; + GUIinputs.name = FilenametosaveasEditField.Value; + + if strcmpi(Switch.Value, 'Multipass') + GUIinputs.numpass = 'Multi'; + else + GUIinputs.numpass = 'Single'; + end + + GUIinputs.snr = SignaltonoiseratioEditField.Value; + GUIinputs.pkr = PeakratioEditFieldField.Value; + GUIinputs.windowoverlap = 0.5; + GUIinputs.idealresolution = IdealresolutionofoutputdataEditField.Value; + GUIinputs.searchwindowsize = 30; + GUIinputs.minsearcharea = 50; + GUIinputs.minyear = MinimumYearEditField.Value; + GUIinputs.maxyear = MaximumYearEditField.Value; + GUIinputs.minmonth = MinimumMonthEditField.Value; + GUIinputs.maxmonth = MaximumMonthEditField.Value; + GUIinputs.minday = MinimumDayEditField.Value; + GUIinputs.maxday = MaximumDayEditField.Value; + GUIinputs.mininterval = MinimumintervalforimagepairsEditField.Value; + GUIinputs.maxinterval = MaximumintervalforimagepairsEditField.Value; + + if strcmpi(Switch_2.Value, 'ON') + GUIinputs.CLAHE = 1; + else + GUIinputs.CLAHE = 0; + end + + if strcmpi(Switch_3.Value, 'ON') + GUIinputs.hipass = 1; + else + GUIinputs.hipass = 0; + end + + if strcmpi(Switch_4.Value, 'ON') + GUIinputs.intenscap = 1; + else + GUIinputs.intenscap = 0; + end + + if strcmpi(Switch_5.Value, 'ON') + GUIinputs.NAOF = 1; + else + GUIinputs.NAOF = 0; + end + + GUIinputs.CLAHEsize = CLAHEsizeEditField.Value; + GUIinputs.hipasssize = HighpasssizeEditField.Value; + + if strcmpi(Switch_10.Value, 'ON') + GUIinputs.sobel = 1; + else + GUIinputs.sobel = 0; + end + + if strcmpi(Switch_11.Value, 'ON') + GUIinputs.laplacian = 1; + else + GUIinputs.laplacian = 0; + end + + if strcmpi(Switch_6.Value, 'ON') + GUIinputs.savearrays = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savearrays = 'No'; + end + + if strcmpi(Switch_7.Value, 'ON') + GUIinputs.savekeyvel = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savekeyvel = 'No'; + end + + if strcmpi(Switch_8.Value, 'ON') + GUIinputs.savegeotiff = 'Yes'; %'Yes' or 'No' + else + GUIinputs.savegeotiff = 'No'; + end + + GUIinputs.imageformat = FormatofimagestosaveEditField.Value; + GUIinputs.maxvel = MaximumVelocityEditField.Value; + GUIinputs.excudedangle1.min = Excludedangle1minimumEditField.Value; + GUIinputs.excudedangle1.max = Excludedangle1maximumEditField.Value; + GUIinputs.excudedangle2.min = Excludedangle2minimumEditField.Value; + GUIinputs.excudedangle2.max = Excludedangle2maximumEditField.Value; + + if strcmpi(Switch_12.Value, 'YES') + GUIinputs.stable = 'Yes'; %'Yes' or 'No' + else + GUIinputs.stable = 'No'; + end + + if strcmpi(Switch_9.Value, 'YES') + GUIinputs.excludeangle = 'Yes'; %'Yes' or 'No' + else + GUIinputs.excludeangle = 'No'; + end + + if strcmpi(SmoothingofcompositearrayofallvelocitiesDropDown.Value, 'Smoothing in time and space') + GUIinputs.finalsmooth = 'Time and Space'; + elseif strcmpi(SmoothingofcompositearrayofallvelocitiesDropDown.Value, 'Smoothing in time') + GUIinputs.finalsmooth = 'Time'; + else + GUIinputs.finalsmooth = 'None'; + end + + GUIinputs.nummonthiter = NumberofiterationsformonthlyvelocitiesEditField.Value; + + %Rename this array + inputs = GUIinputs; + + %Call the viable image calculation script. + GIVruntime(inputs); + end + + +end + +>>>>>>> f7a161863e66aebc8cf45066f0da532c2fb7904e diff --git a/functions/image cross correlation/GIVcore.m b/functions/image cross correlation/GIVcore.m index b20fb9e..fcb64f0 100644 --- a/functions/image cross correlation/GIVcore.m +++ b/functions/image cross correlation/GIVcore.m @@ -95,6 +95,11 @@ inputs.realresolution = ceil(inputs.idealresolution/mean_resolution)*8; % +%Check if cell size is at least 32 (needs to be at least 8*4) +if inputs.realresolution < 32 + inputs.realresolution = 32; +end + %Initialize some parameters newcol1 = {}; %import to external (new) array in order to be parralelized newcol2 = {}; diff --git a/functions/image cross correlation/GIVtrackmultifinal.asv b/functions/image cross correlation/GIVtrackmultifinal.asv deleted file mode 100644 index 5b75198..0000000 --- a/functions/image cross correlation/GIVtrackmultifinal.asv +++ /dev/null @@ -1,231 +0,0 @@ -function [up,vp,SnR2,SnR]=GIVtrackmultifinal(A,B,winsize,overlap,initialdx,initialdy) -% function [x,y,u,v,SnR,PeakHeight]=finalpass(A,B,winsize,overlap,initialdx,initialdy) -% -% Provides the final pass to get the displacements with -% subpixel resolution. -% -% - - -%This function is based upon a multipass solver written by Kristian Sveen -%as part of the matPIV toolbox. It has been adapted for use as part of GIV. -%It is distributed under the terms of the Gnu General Public License -%manager. - - -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - %% GLACIER IMAGE VELOCIMETRY (GIV) %% -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%Code written by Max Van Wyk de Vries @ University of Minnesota -%Credit to Ben Popken and Andrew Wickert for portions of the toolbox. -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%Portions of this toolbox are based on a number of codes written by -%previous authors, including matPIV, IMGRAFT, PIVLAB, M_Map and more. -%Credit and thanks are due to the authors of these toolboxes, and for -%sharing their codes online. See the user manual for a full list of third -%party codes used here. Accordingly, you are free to share, edit and -%add to this GIV code. Please give us credit if you do, and share your code -%with the same conditions as this. - -% Read the associated paper here: -% doi.org/10.5194/tc-15-2115-2021 -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - %Version 1.0, Spring-Summer 2021% -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - %Feel free to contact me at vanwy048@umn.edu% - -if length(winsize)==1 - M=winsize; -else - M=winsize(1); winsize=winsize(2); -end -cj=1; -[sy,sx]=size(A); - -% Allocate space for matrixes -xp=zeros(ceil((size(A,1)-winsize)/((1-overlap)*winsize))+1, ... - ceil((size(A,2)-M)/((1-overlap)*M))+1); -yp=xp; up=xp; vp=xp; SnR2=xp; SnR=xp; Pkh=xp; - -IN=zeros(size(A)); - -%%%%%%%%%%%%%%% MAIN LOOP %%%%%%%%%%%%%%%%%%%%%%%%% -tic -for jj=1:((1-overlap)*winsize):sy-winsize+1 - ci=1; - for ii=1:((1-overlap)*M):sx-M+1 - if IN(jj+winsize/2,ii+M/2)~=1 - if isnan(initialdx(cj,ci)) - initialdx(cj,ci)=0; - end - if isnan(initialdy(cj,ci)) - initialdy(cj,ci)=0; - end - if jj+initialdy(cj,ci)<1 - initialdy(cj,ci)=1-jj; - elseif jj+initialdy(cj,ci)>sy-winsize+1 - initialdy(cj,ci)=sy-winsize+1-jj; - end - if ii+initialdx(cj,ci)<1 - initialdx(cj,ci)=1-ii; - elseif ii+initialdx(cj,ci)>sx-M+1 - initialdx(cj,ci)=sx-M+1-ii; - end - D2=B(jj+initialdy(cj,ci):jj+winsize-1+initialdy(cj,ci),ii+initialdx(cj,ci):ii+M-1+initialdx(cj,ci)); - E=A(jj:jj+winsize-1,ii:ii+M-1); - stad1=std(E(:)); - stad2=std(D2(:)); - if stad1==0, stad1=1; end - if stad2==0, stad2=1; end - E=E-mean(E(:)); - F=D2-mean(D2(:)); - %E(E<0)=0; F(F<0)=0; - - %%%%%%%%%%%%%%%%%%%%%% Calculate the normalized correlation: - R=xcorrelate(E,F)./(winsize*M*stad1*stad2); - %%%%%%%%%%%%%%%%%%%%%% Find the position of the maximal value of R - %%%%%%%%%%%%%%%%%%%%%% _IF_ the standard deviation is NOT NaN. - if all(~isnan(R(:))) && ~all(R(:)==0) %~isnan(stad1) & ~isnan(stad2) - if size(R,1)==(winsize-1) - [max_y1,max_x1]=find(R==max(R(:))); - - else - [max_y1,max_x1]=find(R==max(max(R(0.5*winsize+2:1.5*winsize-3,... - 0.5*M+2:1.5*M-3)))); - end - if length(max_x1)>1 - max_x1=round(sum(max_x1.^2)./sum(max_x1)); - max_y1=round(sum(max_y1.^2)./sum(max_y1)); - end - if max_x1==1, max_x1=2; end - if max_y1==1, max_y1=2; end - - - %Sub-pixel estimator: - % 3-point peak fit using centroid, gaussian (default) - % or parabolic fit -% [x0, y0]=GIVtrackmultipeak(max_x1,max_y1,R(max_y1,max_x1),... -% R(max_y1,max_x1-1),R(max_y1,max_x1+1),... -% R(max_y1-1,max_x1),R(max_y1+1,max_x1),3,[M,winsize]); - [x0, y0]=subpeak(max_x1,max_y1,R); - - x0 = x0-max_x1; y0 = y0-max_y1; - - % Find the signal to Noise ratio - R2=R; - try - %consider changing this from try-catch to a simpler - %distance check. The key here is the distance tot he - %image edge. When peak is close to edge, this NaN - %allocation may fail. - R2(max_y1-3:max_y1+3,max_x1-3:max_x1+3)=NaN; - catch - R2(max_y1-1:max_y1+1,max_x1-1:max_x1+1)=NaN; - end - - if size(R,1)==(winsize-1) - [p2_y2,p2_x2]=find(R2==max(R2(:))); - else - [p2_y2,p2_x2]=find(R2==max(max(R2(0.5*winsize:1.5*winsize-1,0.5*M:1.5*M-1)))); - end - - if length(p2_x2)>=1 - p2_x2=p2_x2(round(length(p2_x2)/2)); - p2_y2=p2_y2(round(length(p2_y2)/2)); - - % signal to noise: - pkr=R(max_y1,max_x1)/R2(p2_y2,p2_x2); - snr=R(max_y1,max_x1)/nanmean(abs(R),'all'); - - else - - % signal to noise: - pkr=NaN; - snr=NaN; - - end - - - %%%%%%%%%%%%%%%%%%%%%% Store the displacements, SnR and Peak Height. - up(cj,ci)=(-x0+initialdx(cj,ci)); - vp(cj,ci)=(-y0+initialdy(cj,ci)); - xp(cj,ci)=(ii+(M/2)-1); - yp(cj,ci)=(jj+(winsize/2)-1); - SnR(cj,ci)=pkr; - SnR2(cj,ci)=snr; - Pkh(cj,ci)=R(max_y1,max_x1); - else - up(cj,ci)=NaN; vp(cj,ci)=NaN; SnR(cj,ci)=NaN; Pkh(cj,ci)=0; - xp(cj,ci)=(ii+(M/2)-1); - yp(cj,ci)=(jj+(winsize/2)-1); - end - ci=ci+1; - else - xp(cj,ci)=(M/2)+ii-1; - yp(cj,ci)=(winsize/2)+jj-1; - up(cj,ci)=NaN; vp(cj,ci)=NaN; - SnR(cj,ci)=NaN; Pkh(cj,ci)=NaN;ci=ci+1; - end - end - cj=cj+1; -end - - -% now we inline the function xcorrelate to shave off some time. - function c = xcorrelate(a,b) - % c = xcorrelate(a,b) - % - % - % Two-dimensional cross-correlation using Fourier transforms. - - %This function is based upon an adaptation of the xcorrf tool written by - %R. Johnson. It has been adapted for use as part of GIV. - - - if nargin<3 - pad='yes'; - end - - - [ma,na] = size(a); - if nargin == 1 - b = a; - end - [mb,nb] = size(b); - % make reverse conjugate of one array - b = conj(b(mb:-1:1,nb:-1:1)); - - if strcmp(pad,'yes') - % use power of 2 transform lengths - mf = 2^nextpow2(ma+mb); - nf = 2^nextpow2(na+nb); - at = fft2(b,mf,nf); - bt = fft2(a,mf,nf); - elseif strcmp(pad,'no') - at = fft2(b); - bt = fft2(a); - else - disp('Wrong input to xcorrelate'); return - end - - % multiply transforms then inverse transform - c = ifft2(at.*bt); - % make real output for real input - if ~any(any(imag(a))) && ~any(any(imag(b))) - c = real(c); - end - % trim to standard size - - if strcmp(pad,'yes') - c(ma+mb:mf,:) = []; - c(:,na+nb:nf) = []; - elseif strcmp(pad,'no') - c=fftshift(c(1:end-1,1:end-1)); - - c(ma+mb:mf,:) = []; - c(:,na+nb:nf) = []; - end - - - end -end \ No newline at end of file diff --git a/functions/image cross correlation/GIVtrackmultifinal.m b/functions/image cross correlation/GIVtrackmultifinal.m index 59e96b5..02e9fa8 100644 --- a/functions/image cross correlation/GIVtrackmultifinal.m +++ b/functions/image cross correlation/GIVtrackmultifinal.m @@ -106,7 +106,7 @@ % or parabolic fit [x0, y0]=GIVtrackmultipeak(max_x1,max_y1,R(max_y1,max_x1),... R(max_y1,max_x1-1),R(max_y1,max_x1+1),... - R(max_y1-1,max_x1),R(max_y1+1,max_x1),3,[M,winsize]); + R(max_y1-1,max_x1),R(max_y1+1,max_x1),2,[M,winsize]); % Find the signal to Noise ratio diff --git a/functions/image cross correlation/GIVtrackmultifirst.m b/functions/image cross correlation/GIVtrackmultifirst.m index 08beab4..19cf94d 100644 --- a/functions/image cross correlation/GIVtrackmultifirst.m +++ b/functions/image cross correlation/GIVtrackmultifirst.m @@ -87,6 +87,7 @@ else [max_y1,max_x1]=find(R==max(max(R(0.5*winsize+2:1.5*winsize-3,0.5*M+2:1.5*M-3)))); end + if length(max_x1)>1 max_x1=round(sum(max_x1.*(1:length(max_x1))')./sum(max_x1)); diff --git a/functions/image processing and filters/cropmask.m b/functions/image processing and filters/cropmask.m index deb1f7f..14b61d6 100644 --- a/functions/image processing and filters/cropmask.m +++ b/functions/image processing and filters/cropmask.m @@ -93,12 +93,12 @@ rightsizex = min(size_matrix(:,2)); if size(mask_0_1,1) ~= rightsizey || size(mask_0_1,2) ~=rightsizex - mask_0_1 = interp2(mask_0_1, linspace(1, size(mask_0_1,2), rightsizex).', linspace(1, size(mask_0_1,1), rightsizey)); + mask_0_1 = interp2(double(mask_0_1), linspace(1, size(mask_0_1,2), rightsizex).', linspace(1, size(mask_0_1,1), rightsizey)); end for i = 2:size(images,1) if size(images{i,3},1) ~= rightsizey || size(images{i,3},2) ~=rightsizex - images{i,3} = interp2(images{i,3}, (linspace(1, size(images{i,3},2), rightsizex).'), linspace(1, size(images{i,3},1), rightsizey)); + images{i,3} = interp2(double(images{i,3}), (linspace(1, size(images{i,3},2), rightsizex).'), linspace(1, size(images{i,3},1), rightsizey)); end end end diff --git a/functions/post-processing/fd_outlier_filter.m b/functions/post-processing/fd_outlier_filter.m new file mode 100644 index 0000000..3aced1c --- /dev/null +++ b/functions/post-processing/fd_outlier_filter.m @@ -0,0 +1,48 @@ +function [out_limits_fd] = fd_outlier_filter(full_fd,mean_fd,std_fd,num_stand_dev) +%Input = multi-array of flow directions +% mean flow direction +% stv of flow direction +%Output = binary array with location of outliers +% +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %% GLACIER IMAGE VELOCIMETRY (GIV) %% +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %Code written by Max Van Wyk de Vries @ University of Minnesota +% %Credit to Ben Popken and Andrew Wickert for portions of the toolbox. +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %Portions of this toolbox are based on a number of codes written by +% %previous authors, including matPIV, IMGRAFT, PIVLAB, M_Map and more. +% %Credit and thanks are due to the authors of these toolboxes, and for +% %sharing their codes online. See the user manual for a full list of third +% %party codes used here. Accordingly, you are free to share, edit and +% %add to this GIV code. Please give us credit if you do, and share your code +% %with the same conditions as this. +% +% % Read the associated paper here: +% % doi.org/10.5194/tc-15-2115-2021 +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %Version 1.0, Spring-Summer 2021% +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %Feel free to contact me at vanwy048@umn.edu% + + + +%Calculate the radial distance between each measurement and the mean flow +%direction + +fd_dist = min(abs(full_fd-mean_fd),360-abs(full_fd-mean_fd)); + +% Create stv outlier detection matrix + +fd_detection = (repmat(std_fd,size(full_fd,1),1)*num_stand_dev); + + +%Identify regions where the radial distance is larger than XX times the stv + +out_limits_fd = double(fd_dist>fd_detection); + +%Remove any pixels with a stv*num_stand_dev larger than 180 (in this case, +%there are no outliers + +out_limits_fd = out_limits_fd-double(fd_detection>=180); + diff --git a/functions/post-processing/filtall.asv b/functions/post-processing/filtall.asv new file mode 100644 index 0000000..f2fccdc --- /dev/null +++ b/functions/post-processing/filtall.asv @@ -0,0 +1,416 @@ +function [images, images_stack]=filtall(images,inputs) +% %% Take newly created velocity maps and filter based upon entire time series +% % +% %Aim here is to use the information from the entire dataset to filter +% %individual velocity maps- including the knowledge that glacier velocity +% %will be smoothly varying through time (no jumps) and will only vary at a +% %relatively slow rate in general. +% +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %% GLACIER IMAGE VELOCIMETRY (GIV) %% +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %Code written by Max Van Wyk de Vries @ University of Minnesota +% %Credit to Ben Popken and Andrew Wickert for portions of the toolbox. +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %Portions of this toolbox are based on a number of codes written by +% %previous authors, including matPIV, IMGRAFT, PIVLAB, M_Map and more. +% %Credit and thanks are due to the authors of these toolboxes, and for +% %sharing their codes online. See the user manual for a full list of third +% %party codes used here. Accordingly, you are free to share, edit and +% %add to this GIV code. Please give us credit if you do, and share your code +% %with the same conditions as this. +% +% % Read the associated paper here: +% % doi.org/10.5194/tc-15-2115-2021 +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %Version 1.0, Spring-Summer 2021% +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% %Feel free to contact me at vanwy048@umn.edu% + +%% calculate velocity statistics + +%Initialize +meta_dum = 0; +array_pos = 0; +emptycount_inner = 0; +emptycount_outer = 0; +number_in_range = 0; +number_with_values = (size(images,2)-6)/2; + +%loop through all images +for time_loop = 1:number_with_values + position1 = 6+time_loop*2; + for inner_loop = 2:inputs.numimages-time_loop + if ~isempty(images{inner_loop,position1}) + number_in_range = number_in_range + 1; + end + end +end + +%Create matrix to store all x and y velocity components +full_u = NaN(number_in_range,inputs.sizevel(1)*inputs.sizevel(2)+2); +full_v = NaN(number_in_range,inputs.sizevel(1)*inputs.sizevel(2)+2); + +%loop through all +for time_loop = 1:number_with_values + position1 = time_loop + 6 + array_pos; + position2 = time_loop + 7 + array_pos; + for inner_loop = 2:inputs.numimages-time_loop + if ~isempty(images{inner_loop,position1}) + %calculate time interval + time_gap = (images{inner_loop+time_loop,5}-images{inner_loop,5}); + % Here calculates a median date for the interval that we will use + date_current = round(images{inner_loop+time_loop,4}-(time_gap/2)); + full_u(inner_loop+meta_dum-1-emptycount_inner-emptycount_outer,1)= date_current; + full_v(inner_loop+meta_dum-1-emptycount_inner-emptycount_outer,1)= date_current; + + %Recreate x-y velocity from V and fd + [u,v] = Vtoxy(images{inner_loop,position1},images{inner_loop,position2}); + + %Linearize each velocity map + lin_u = reshape(u,[1 inputs.sizevel(1)*inputs.sizevel(2)]); + full_u(inner_loop+meta_dum-1-emptycount_inner-emptycount_outer,3:inputs.sizevel(1)*inputs.sizevel(2)+2)= lin_u; + + %Linearize each flow direction map + lin_v = reshape(v,[1 inputs.sizevel(1)*inputs.sizevel(2)]); + full_v(inner_loop+meta_dum-1-emptycount_inner-emptycount_outer,3:inputs.sizevel(1)*inputs.sizevel(2)+2)= lin_v; + else + emptycount_inner = emptycount_inner + 1; + end + end + + emptycount_outer = emptycount_outer + emptycount_inner; + emptycount_inner = 0; + meta_dum = meta_dum + inputs.numimages-time_loop-1; + array_pos = array_pos + 1; +end + +%Create a column with sequential numbers in order for the initial order to +%be recovered (we are going to sort by date to filter). Add to second +%column so that sortrows function can sort by date (looks at first column). +numbers_order=(1:size(full_u,1))'; +full_u(:,2) = numbers_order; +full_v(:,2) = numbers_order; + +% Sort based on date. +full_u = sortrows(full_u, 1); +full_v = sortrows(full_v, 1); + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +%strip the numbering and date list off the main array and store separately. +only_numlist = full_u(:,2); %list of numbers for re-ordering + +%delete first two columns to leave just data itself +full_u(:,1:2)=[]; +full_v(:,1:2)=[]; + +%% remove outliers in stack +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% If only a single pair, no averaging is possible. Skip this step and +% create 'average' maps anyway + +if size(full_u,1) == 1 + + %Mean + [mean_vel,mean_fd] = xytoV_basic(full_u,full_v); + + ordered_vel = mean_vel; + ordered_fd = mean_fd; + + mean_vel=reshape(mean_vel,inputs.sizevel); + mean_fd=reshape(mean_fd,inputs.sizevel); + + %Median + median_vel=mean_vel; + median_fd=mean_fd; + + %Min + min_vel=mean_vel; + min_fd=mean_fd; + + %Max + max_vel=mean_vel; + max_fd=mean_fd; + + %std + std_vel = zeros(size(mean_fd,1),size(mean_fd,2)); + std_fd = std_vel; + + + perc_error_vel = std_vel; + + mean_u = full_u; + + mean_v = full_v; + + median_u = full_u; + + median_v = full_v; + + std_u = std_vel; + + std_v = std_vel; + + +else % More than one image pair, do averaging + + + num_stand_dev = 10; %Number of standard deviations to consider outliers. Set to 2 by default + + + %%%%%%%%%%%%%%%% First, x and y components + %u (x velocity component). + mean_u = nanmean(full_u); + std_u = nanstd(full_u); + + % find points in the stack that are more than xx standard deviations outside + % of the mean. + low_u = mean_u - num_stand_dev .* std_u; + high_u = mean_u + num_stand_dev .* std_u; + full_low_u = []; + full_high_u = []; + + %faster than repmat + for index = 1:size(full_u,1) %make full arrays of limits + full_low_u(index,:) = low_u; + full_high_u(index,:) = high_u; + end + + %create array of where it is out of limits + out_limits_low_u = double(full_u < full_low_u); + out_limits_high_u = double(full_u > full_high_u); + out_limits_u = double((out_limits_low_u + out_limits_high_u)>=1); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %v (y velocity component). + mean_v = nanmean(full_v); + std_v = nanstd(full_v); + + % find points in the stack that are more than xx standard deviations outside + % of the mean. + low_v = mean_v - num_stand_dev * std_v; + high_v = mean_v + num_stand_dev * std_v; + full_low_v = []; + full_high_v = []; + + %faster than repmat + for index = 1:size(full_v,1) %make full arrays of limits + full_low_v(index,:) = low_v; + full_high_v(index,:) = high_v; + end + + %create array of where it is out of limits + out_limits_low_v = double(full_v < full_low_v); + out_limits_high_v = double(full_v > full_high_v); + out_limits_v = double((out_limits_low_v + out_limits_high_v)>=1); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %remove values more than 2 std away from mean in either u or v + full_u(out_limits_u==1)=NaN; + full_v(out_limits_u==1)=NaN; + full_u(out_limits_v==1)=NaN; + full_v(out_limits_v==1)=NaN; + + + + %%%%%%%%%%%%%%%% Next Velocity and Standard deviation + + %recalculate Vel and flowdir + [full_vel,full_fd]... + =xytoV_basic(full_u,full_v); + + + %Calculate mean and standard deviation using circular statistics + [mean_fd,std_fd] = GIV_circstats(full_fd); + + % find points in the stack that are more than num_stand_dev standard deviations outside + % of the mean. + low_fd = mean_fd - num_stand_dev * std_fd; + templow_fd = mean_fd - num_stand_dev * std_fd; + high_fd = mean_fd + num_stand_dev * std_fd; + temphigh_fd = mean_fd + num_stand_dev * std_fd; + %Find if any are out of 0-360 circle + poslow = low_fd<0; + poshigh = temphigh_fd>360; + %flip where needed + high_fd(poslow)=360+templow_fd(poslow); + low_fd(poslow)=temphigh_fd(poslow); + low_fd(poshigh)=temphigh_fd(poshigh)-360; + high_fd(poshigh)=templow_fd(poshigh); + + %make full arrays of limits + full_low_fd = []; + full_high_fd = []; + for index = 1:size(full_fd,1) + full_low_fd(index,:) = low_fd; + full_high_fd(index,:) = high_fd; + end + + % where is it out of limit in BOTH (i.e. not caused by the boundary) + out_limits_fd = double((double(full_fdfull_high_fd))>=1); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %Secondly, lets do this for the velocities. + mean_vel = nanmean(full_vel); + std_vel = nanstd(full_vel); + + % find points in the stack that are more than 2 standard deviations outside + % of the mean. + low_vel = mean_vel - num_stand_dev * std_vel; + high_vel = mean_vel + num_stand_dev * std_vel; + full_low_vel = []; + full_high_vel = []; + + for index = 1:size(full_vel,1) %make full arrays of limits + full_low_vel(index,:) = low_vel; + full_high_vel(index,:) = high_vel; + end + + out_limits_low_vel = double(full_vel < full_low_vel); + out_limits_high_vel = double(full_vel > full_high_vel); + out_limits_vel = double((out_limits_low_vel + out_limits_high_vel)>=1); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %remove values more than 2 std away from mean in data + full_v(out_limits_fd==1)=NaN; + full_u(out_limits_fd==1)=NaN; + full_v(out_limits_vel==1)=NaN; + full_u(out_limits_vel==1)=NaN; + + + %% Filter data in time and space if requested + + + %TIME FILTER + if strcmpi(inputs.finalsmooth, 'Time') + full_u = nanfill_time(full_u, inputs, 2, 3); + full_v = nanfill_time(full_v, inputs, 2, 3); + %TIME AND SPACE FILTER + elseif strcmpi(inputs.finalsmooth, 'Time and Space') + full_u = nanfill_timeandspace(full_u, inputs, 4, 3); + full_v = nanfill_timeandspace(full_v, inputs, 4, 3); + end + + + %% Calculate statistics in x-y space + + %Mean + mean_u = nanmean(full_u); + mean_v = nanmean(full_v); + + %std + std_u = nanstd(full_u); + std_v = nanstd(full_v); + + %Median + median_u = nanmedian(full_u); + median_v = nanmedian(full_v); + + %local smoothing and NaN filling + mean_u = nanfill(mean_u, 4,3); + mean_v = nanfill(mean_v, 4,3); + std_u = nanfill(std_u, 4,3); + std_v = nanfill(std_v, 4,3); + median_u = nanfill(median_u, 4,3); + median_v = nanfill(median_v, 4,3); + + + + %% Convert to V - fd space + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %Calculate some additional statistics: Minimum, Maximum, Median, 'Percent + %error' + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %Mean + [mean_vel,mean_fd] = xytoV_basic(mean_u,mean_v); + mean_vel=reshape(mean_vel,inputs.sizevel); + mean_fd=reshape(mean_fd,inputs.sizevel); + + %Median + [median_vel,median_fd] = xytoV_basic(median_u,median_v); + median_vel=reshape(median_vel,inputs.sizevel); + median_fd=reshape(median_fd,inputs.sizevel); + + %Min + [min_vel,min_fd] = xytoV_basic(min(mean_u-std_u,mean_u+std_u),min(mean_v-std_v,mean_v+std_v)); + min_vel=reshape(min_vel,inputs.sizevel); + min_fd=reshape(min_fd,inputs.sizevel); + + %Max + [max_vel,max_fd] = xytoV_basic(max(mean_u-std_u,mean_u+std_u),max(mean_v-std_v,mean_v+std_v)); + max_vel=reshape(max_vel,inputs.sizevel); + max_fd=reshape(max_fd,inputs.sizevel); + + %std + std_vel = abs(max_vel-min_vel); + fdtemp(:,:,1) = max_fd-min_fd; fdtemp(:,:,2) = min_fd-max_fd; fdtemp(:,:,3) = max_fd+360-min_fd; fdtemp(:,:,4) = min_fd+360-max_fd; + std_fd = min(fdtemp,[],3); %Smallest angle + std_vel=reshape(std_vel,inputs.sizevel); + std_fd=reshape(std_fd,inputs.sizevel); + + perc_error_vel = 100*std_vel./mean_vel; + + %% Covert back to V, fd stack and order + + + %Add back in numbering scheme, sort rows back into order then remove it + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + [ordered_vel(:,2:inputs.sizevel(1)*inputs.sizevel(2)+1), ... + ordered_fd(:,2:inputs.sizevel(1)*inputs.sizevel(2)+1)]... + =xytoV_basic(full_u,full_v); + ordered_vel(:,1)=only_numlist; + ordered_fd(:,1)=only_numlist; + ordered_vel = sortrows(ordered_vel, 1); + ordered_fd = sortrows(ordered_fd, 1); + ordered_vel(:,1)=[]; + ordered_fd(:,1)=[]; + +end + + +%% Create output array with useful data in it + +%initialize stack +images_stack = {}; + +%Load data into it +%All velocities (for plotting timeseries and querying points) +images_stack{1,1} = 'Full Velocity Array'; images_stack{1,2} = ordered_vel; + +%All flow directions (for plotting timeseries and querying points) +images_stack{2,1} = 'Full Flow Direction Array'; images_stack{2,2} = ordered_fd; + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +%Load additional statistics into output array +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +images_stack{3,1} = 'Mean Velocity'; images_stack{3,2} = mean_vel; +images_stack{5,1} = 'Velocity Standard Deviation'; images_stack{5,2} = std_vel; + +images_stack{4,1} = 'Mean Flow Direction'; images_stack{4,2} = mean_fd; +images_stack{6,1} = 'Flow Direction Standard Deviation'; images_stack{6,2} = std_fd; + +images_stack{7,1} = 'Median Velocity'; images_stack{7,2} = median_vel; +images_stack{8,1} = 'Median Flow Direction'; images_stack{8,2} = median_fd; + +images_stack{9,1} = 'Maximum Velocity'; images_stack{9,2} = max_vel; +images_stack{10,1} = 'Maximum Flow Direction'; images_stack{10,2} = max_fd; + +images_stack{11,1} = 'Minimum Velocity'; images_stack{11,2} = min_vel; +images_stack{12,1} = 'Minimum flow direction'; images_stack{12,2} = min_fd; + +images_stack{13,1} = 'Percentage Error-Variability in Velocity'; images_stack{13,2} = perc_error_vel; + +%Mean and standard deviation of velocity arrays +images_stack{14,1} = 'Mean Velocity x component'; images_stack{14,2} = reshape(mean_u,inputs.sizevel); +images_stack{15,1} = 'Velocity x component Standard Deviation'; images_stack{15,2} = reshape(std_u,inputs.sizevel); + +%Mean and standard deviation of velocity arrays +images_stack{16,1} = 'Mean Velocity y component'; images_stack{16,2} = reshape(mean_v,inputs.sizevel); +images_stack{17,1} = 'Velocity y component Standard Deviation'; images_stack{17,2} = reshape(std_v,inputs.sizevel); + +images_stack{18,1} = 'Median Velocity x component'; images_stack{18,2} = reshape(median_u,inputs.sizevel); +images_stack{19,1} = 'Median Velocity y component'; images_stack{19,2} = reshape(median_v,inputs.sizevel); + + diff --git a/functions/post-processing/filtall.m b/functions/post-processing/filtall.m index b72f758..cfed76e 100644 --- a/functions/post-processing/filtall.m +++ b/functions/post-processing/filtall.m @@ -15,12 +15,12 @@ % %Portions of this toolbox are based on a number of codes written by % %previous authors, including matPIV, IMGRAFT, PIVLAB, M_Map and more. % %Credit and thanks are due to the authors of these toolboxes, and for -% %sharing their codes online. See the user manual for a full list of third +% %sharing their codes online. See the user manual for a full list of third % %party codes used here. Accordingly, you are free to share, edit and -% %add to this GIV code. Please give us credit if you do, and share your code +% %add to this GIV code. Please give us credit if you do, and share your code % %with the same conditions as this. -% -% % Read the associated paper here: +% +% % Read the associated paper here: % % doi.org/10.5194/tc-15-2115-2021 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % %Version 1.0, Spring-Summer 2021% @@ -107,217 +107,245 @@ %% remove outliers in stack %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -num_stand_dev = 2; %Number of standard deviations to consider outliers. Set to 2 by default - - -%%%%%%%%%%%%%%%% First, x and y components -%u (x velocity component). -mean_u = nanmean(full_u); -std_u = nanstd(full_u); - -% find points in the stack that are more than xx standard deviations outside -% of the mean. -low_u = mean_u - num_stand_dev * std_u; -high_u = mean_u + num_stand_dev * std_u; -full_low_u = []; -full_high_u = []; - -%faster than repmat -for index = 1:size(full_u,1) %make full arrays of limits - full_low_u(index,:) = low_u; - full_high_u(index,:) = high_u; -end - -%create array of where it is out of limits -out_limits_low_u = double(full_u < full_low_u); -out_limits_high_u = double(full_u > full_high_u); -out_limits_u = double((out_limits_low_u + out_limits_high_u)>=1); -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%v (y velocity component). -mean_v = nanmean(full_v); -std_v = nanstd(full_v); - -% find points in the stack that are more than xx standard deviations outside -% of the mean. -low_v = mean_v - num_stand_dev * std_v; -high_v = mean_v + num_stand_dev * std_v; -full_low_v = []; -full_high_v = []; - -%faster than repmat -for index = 1:size(full_v,1) %make full arrays of limits - full_low_v(index,:) = low_v; - full_high_v(index,:) = high_v; -end - -%create array of where it is out of limits -out_limits_low_v = double(full_v < full_low_v); -out_limits_high_v = double(full_v > full_high_v); -out_limits_v = double((out_limits_low_v + out_limits_high_v)>=1); -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%remove values more than 2 std away from mean in either u or v -full_u(out_limits_u==1)=NaN; -full_v(out_limits_u==1)=NaN; -full_u(out_limits_v==1)=NaN; -full_v(out_limits_v==1)=NaN; - - - -%%%%%%%%%%%%%%%% Next Velocity and Standard deviation - -%recalculate Vel and flowdir -[full_vel,full_fd]... - =xytoV_basic(full_u,full_v); - - -%Calculate mean and standard deviation using circular statistics -[mean_fd,std_fd] = GIV_circstats(full_fd); - -% find points in the stack that are more than 2 standard deviations outside -% of the mean. -low_fd = mean_fd - 2 * std_fd; -templow_fd = mean_fd - 2 * std_fd; -high_fd = mean_fd + 2 * std_fd; -temphigh_fd = mean_fd + 2 * std_fd; -%Find if any are out of 0-360 circle -poslow = low_fd<0; -poshigh = temphigh_fd>360; -%flip where needed -high_fd(poslow)=360+templow_fd(poslow); -low_fd(poslow)=temphigh_fd(poslow); -low_fd(poshigh)=temphigh_fd(poshigh)-360; -high_fd(poshigh)=templow_fd(poshigh); - -%make full arrays of limits -full_low_fd = []; -full_high_fd = []; -for index = 1:size(full_fd,1) - full_low_fd(index,:) = low_fd; - full_high_fd(index,:) = high_fd; -end - -% where is it out of limit in BOTH (i.e. not caused by the boundary) -out_limits_fd = double((double(full_fdfull_high_fd))>=1); - -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%Secondly, lets do this for the velocities. -mean_vel = nanmean(full_vel); -std_vel = nanstd(full_vel); - -% find points in the stack that are more than 2 standard deviations outside -% of the mean. -low_vel = mean_vel - 2 * std_vel; -high_vel = mean_vel + 2 * std_vel; -full_low_vel = []; -full_high_vel = []; - -for index = 1:size(full_vel,1) %make full arrays of limits - full_low_vel(index,:) = low_vel; - full_high_vel(index,:) = high_vel; -end - -out_limits_low_vel = double(full_vel < full_low_vel); -out_limits_high_vel = double(full_vel > full_high_vel); -out_limits_vel = double((out_limits_low_vel + out_limits_high_vel)>=1); - -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%remove values more than 2 std away from mean in data -full_v(out_limits_fd==1)=NaN; -full_u(out_limits_fd==1)=NaN; -full_v(out_limits_vel==1)=NaN; -full_u(out_limits_vel==1)=NaN; +% If only a single pair, no averaging is possible. Skip this step and +% create 'average' maps anyway - -%% Filter data in time and space if requested - - -%TIME FILTER -if strcmpi(inputs.finalsmooth, 'Time') - full_u = nanfill_time(full_u, inputs, 2, 3); - full_v = nanfill_time(full_v, inputs, 2, 3); - %TIME AND SPACE FILTER -elseif strcmpi(inputs.finalsmooth, 'Time and Space') - full_u = nanfill_timeandspace(full_u, inputs, 4, 3); - full_v = nanfill_timeandspace(full_v, inputs, 4, 3); +if size(full_u,1) == 1 + + %Mean + [mean_vel,mean_fd] = xytoV_basic(full_u,full_v); + + ordered_vel = mean_vel; + ordered_fd = mean_fd; + + mean_vel=reshape(mean_vel,inputs.sizevel); + mean_fd=reshape(mean_fd,inputs.sizevel); + + %Median + median_vel=mean_vel; + median_fd=mean_fd; + + %Min + min_vel=mean_vel; + min_fd=mean_fd; + + %Max + max_vel=mean_vel; + max_fd=mean_fd; + + %std + std_vel = zeros(size(mean_fd,1),size(mean_fd,2)); + std_fd = std_vel; + + + perc_error_vel = std_vel; + + mean_u = full_u; + + mean_v = full_v; + + median_u = full_u; + + median_v = full_v; + + std_u = std_vel; + + std_v = std_vel; + + +else % More than one image pair, do averaging + + + num_stand_dev = 2; %Number of standard deviations to consider outliers. Set to 2 by default + + + %%%%%%%%%%%%%%%% First, x and y components + %u (x velocity component). + mean_u = nanmean(full_u); + std_u = nanstd(full_u); + + % find points in the stack that are more than xx standard deviations outside + % of the mean. + low_u = mean_u - num_stand_dev .* std_u; + high_u = mean_u + num_stand_dev .* std_u; + full_low_u = []; + full_high_u = []; + + %faster than repmat + for index = 1:size(full_u,1) %make full arrays of limits + full_low_u(index,:) = low_u; + full_high_u(index,:) = high_u; + end + + %create array of where it is out of limits + out_limits_low_u = double(full_u < full_low_u); + out_limits_high_u = double(full_u > full_high_u); + out_limits_u = double((out_limits_low_u + out_limits_high_u)>=1); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %v (y velocity component). + mean_v = nanmean(full_v); + std_v = nanstd(full_v); + + % find points in the stack that are more than xx standard deviations outside + % of the mean. + low_v = mean_v - num_stand_dev * std_v; + high_v = mean_v + num_stand_dev * std_v; + full_low_v = []; + full_high_v = []; + + %faster than repmat + for index = 1:size(full_v,1) %make full arrays of limits + full_low_v(index,:) = low_v; + full_high_v(index,:) = high_v; + end + + %create array of where it is out of limits + out_limits_low_v = double(full_v < full_low_v); + out_limits_high_v = double(full_v > full_high_v); + out_limits_v = double((out_limits_low_v + out_limits_high_v)>=1); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %remove values more than 2 std away from mean in either u or v + full_u(out_limits_u==1)=NaN; + full_v(out_limits_u==1)=NaN; + full_u(out_limits_v==1)=NaN; + full_v(out_limits_v==1)=NaN; + + + + %%%%%%%%%%%%%%%% Next Velocity and Standard deviation + + %recalculate Vel and flowdir + [full_vel,full_fd]... + =xytoV_basic(full_u,full_v); + + + %Calculate mean and standard deviation using circular statistics + [mean_fd,std_fd] = GIV_circstats(full_fd); + + %Call seperate function to filter flow directions + [out_limits_fd] = fd_outlier_filter(full_fd,mean_fd,std_fd,num_stand_dev); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %Secondly, lets do this for the velocities. + mean_vel = nanmean(full_vel); + std_vel = nanstd(full_vel); + + % find points in the stack that are more than 2 standard deviations outside + % of the mean. + low_vel = mean_vel - num_stand_dev * std_vel; + high_vel = mean_vel + num_stand_dev * std_vel; + full_low_vel = []; + full_high_vel = []; + + for index = 1:size(full_vel,1) %make full arrays of limits + full_low_vel(index,:) = low_vel; + full_high_vel(index,:) = high_vel; + end + + out_limits_low_vel = double(full_vel < full_low_vel); + out_limits_high_vel = double(full_vel > full_high_vel); + out_limits_vel = double((out_limits_low_vel + out_limits_high_vel)>=1); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %remove values more than 2 std away from mean in data + full_v(out_limits_fd==1)=NaN; + full_u(out_limits_fd==1)=NaN; + full_v(out_limits_vel==1)=NaN; + full_u(out_limits_vel==1)=NaN; + + + %% Filter data in time and space if requested + + + %TIME FILTER + if strcmpi(inputs.finalsmooth, 'Time') + full_u = nanfill_time(full_u, inputs, 2, 3); + full_v = nanfill_time(full_v, inputs, 2, 3); + %TIME AND SPACE FILTER + elseif strcmpi(inputs.finalsmooth, 'Time and Space') + full_u = nanfill_timeandspace(full_u, inputs, 4, 3); + full_v = nanfill_timeandspace(full_v, inputs, 4, 3); + end + + + %% Calculate statistics in x-y space + + %Mean + mean_u = nanmean(full_u); + mean_v = nanmean(full_v); + + %std + std_u = nanstd(full_u); + std_v = nanstd(full_v); + + %Median + median_u = nanmedian(full_u); + median_v = nanmedian(full_v); + + %local smoothing and NaN filling + mean_u = nanfill(mean_u, 4,3); + mean_v = nanfill(mean_v, 4,3); + std_u = nanfill(std_u, 4,3); + std_v = nanfill(std_v, 4,3); + median_u = nanfill(median_u, 4,3); + median_v = nanfill(median_v, 4,3); + + + + %% Convert to V - fd space + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + %Calculate some additional statistics: Minimum, Maximum, Median, 'Percent + %error' + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %Mean + [mean_vel,mean_fd] = xytoV_basic(mean_u,mean_v); + mean_vel=reshape(mean_vel,inputs.sizevel); + mean_fd=reshape(mean_fd,inputs.sizevel); + + %Median + [median_vel,median_fd] = xytoV_basic(median_u,median_v); + median_vel=reshape(median_vel,inputs.sizevel); + median_fd=reshape(median_fd,inputs.sizevel); + + %Min + [min_vel,min_fd] = xytoV_basic(min(mean_u-std_u,mean_u+std_u),min(mean_v-std_v,mean_v+std_v)); + min_vel=reshape(min_vel,inputs.sizevel); + min_fd=reshape(min_fd,inputs.sizevel); + + %Max + [max_vel,max_fd] = xytoV_basic(max(mean_u-std_u,mean_u+std_u),max(mean_v-std_v,mean_v+std_v)); + max_vel=reshape(max_vel,inputs.sizevel); + max_fd=reshape(max_fd,inputs.sizevel); + + %std + std_vel = abs(max_vel-min_vel); + fdtemp(:,:,1) = max_fd-min_fd; fdtemp(:,:,2) = min_fd-max_fd; fdtemp(:,:,3) = max_fd+360-min_fd; fdtemp(:,:,4) = min_fd+360-max_fd; + std_fd = min(fdtemp,[],3); %Smallest angle + std_vel=reshape(std_vel,inputs.sizevel); + std_fd=reshape(std_fd,inputs.sizevel); + + perc_error_vel = 100*std_vel./mean_vel; + + %% Covert back to V, fd stack and order + + + %Add back in numbering scheme, sort rows back into order then remove it + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + [ordered_vel(:,2:inputs.sizevel(1)*inputs.sizevel(2)+1), ... + ordered_fd(:,2:inputs.sizevel(1)*inputs.sizevel(2)+1)]... + =xytoV_basic(full_u,full_v); + ordered_vel(:,1)=only_numlist; + ordered_fd(:,1)=only_numlist; + ordered_vel = sortrows(ordered_vel, 1); + ordered_fd = sortrows(ordered_fd, 1); + ordered_vel(:,1)=[]; + ordered_fd(:,1)=[]; + end -%% Calculate statistics in x-y space - -%Mean -mean_u = nanmean(full_u); -mean_v = nanmean(full_v); - -%std -std_u = nanstd(full_u); -std_v = nanstd(full_v); - -%Median -median_u = nanmedian(full_u); -median_v = nanmedian(full_v); - -%local smoothing and NaN filling -mean_u = nanfill(mean_u, 4,3); -mean_v = nanfill(mean_v, 4,3); -std_u = nanfill(std_u, 4,3); -std_v = nanfill(std_v, 4,3); -median_u = nanfill(median_u, 4,3); -median_v = nanfill(median_v, 4,3); - - - -%% Convert to V - fd space - -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%Calculate some additional statistics: Minimum, Maximum, Median, 'Percent -%error' -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -%Mean -[mean_vel,mean_fd] = xytoV_basic(mean_u,mean_v); -mean_vel=reshape(mean_vel,inputs.sizevel); -mean_fd=reshape(mean_fd,inputs.sizevel); - -%Median -[median_vel,median_fd] = xytoV_basic(median_u,median_v); -median_vel=reshape(median_vel,inputs.sizevel); -median_fd=reshape(median_fd,inputs.sizevel); - -%Min -[min_vel,min_fd] = xytoV_basic(min(mean_u-std_u,mean_u+std_u),min(mean_v-std_v,mean_v+std_v)); -min_vel=reshape(min_vel,inputs.sizevel); -min_fd=reshape(min_fd,inputs.sizevel); - -%Max -[max_vel,max_fd] = xytoV_basic(max(mean_u-std_u,mean_u+std_u),max(mean_v-std_v,mean_v+std_v)); -max_vel=reshape(max_vel,inputs.sizevel); -max_fd=reshape(max_fd,inputs.sizevel); - -%std -std_vel = abs(max_vel-min_vel); -fdtemp(:,:,1) = max_fd-min_fd; fdtemp(:,:,2) = min_fd-max_fd; fdtemp(:,:,3) = max_fd+360-min_fd; fdtemp(:,:,4) = min_fd+360-max_fd; -std_fd = min(fdtemp,[],3); %Smallest angle -std_vel=reshape(std_vel,inputs.sizevel); -std_fd=reshape(std_fd,inputs.sizevel); - -perc_error_vel = 100*std_vel./mean_vel; - -%% Covert back to V, fd stack and order - - -%Add back in numbering scheme, sort rows back into order then remove it -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -[ordered_vel(:,2:inputs.sizevel(1)*inputs.sizevel(2)+1), ... - ordered_fd(:,2:inputs.sizevel(1)*inputs.sizevel(2)+1)]... - =xytoV_basic(full_u,full_v); -ordered_vel(:,1)=only_numlist; -ordered_fd(:,1)=only_numlist; -ordered_vel = sortrows(ordered_vel, 1); -ordered_fd = sortrows(ordered_fd, 1); -ordered_vel(:,1)=[]; -ordered_fd(:,1)=[]; - - %% Create output array with useful data in it %initialize stack