function [normpow, pows, powss, coordinates, srfmed] = bright_area_analysis(geomlist, radlist, tolerance)
%bright_area_analysis: Re-analyzes the MARSIS data from Orosei et al. 2018
%to determing the median normalized basal echo power at each different
%frequency within the bright area of the SPLD
%Inputs:
% geomlist:  structure array containing geometry file names, created using
%            dir command.
% radlist:   structure array containing radargram file names, created using
%            dir command. must be same length as geomlist.
% tolerance: radius around central point to include in bright area, in
%            meters

pows                                = [];
powss                               = [];
normpow                             = [];
coordinates                         = [];
mars                                = referenceEllipsoid('mars'); %create reference ellipse object
center                              = [-81.1 193]; %set center lat and lon


for ii = 1:length(geomlist) %loop through data files
    ii                              %for debugging purposes
    geom                            = importdata(geomlist(ii).name); %import geometry file
    rgram                           = importdata(radlist(ii).name); %import radargram file
    radlist(ii).name                %for debugging purposes
    srfmed(ii)                      = 20.*log10( median( max( rgram(1:420,:) ) ) ); %median surface power from entire track used for normalization, following Orosei et al 2018
    lonlat                          = geom(:,5:6); %grab lat and lon
    latlon                          = circshift(lonlat,1,2); %flip columns so lat is first
    [arclen, az]                    = distance(latlon,center,mars); %find distance from center for every radargram trace
    coordinates                     = vertcat(coordinates,latlon(arclen >= tolerance,:)); %grab lat and lon from bright area
    subgram                         = 20.*log10(rgram(:,arclen>=tolerance)); %grab the section of the radargram corresponding to bright area (within tolerance)
    pows                            = horzcat(pows,max(subgram(1:420,:))); %extract surface power THIS IS HARDCODED FOR OROSEI'S MARSIS DATA
    powss                           = horzcat(powss,max(subgram(420:end,:))); %extract subsurface power THIS IS HARDCODED FOR OROSEI'S MARSIS DATA
    normpow                         = horzcat(normpow,max(subgram(420:end,:)) - srfmed(ii));
end

%normpow                             = powss - pows; %normalize subsurface echo to surface echo

end

