geotran.m
by
Webmaster Legos
—
last modified
Jul 05, 2012 12:37 PM
geotran.m
—
Objective-C source code,
4 kB (4,375 bytes)
File contents
% function % KEY: computes geostrophic transport after geovel output % AUTHOR : A.Ganachaud (ganacho@gulf.mit.edu) , May 97 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %ISOPYCNALS ca sigint= [22.00 26.40 26.80 27.10 27.30 27.50 27.70, ... 36.87 36.94 36.98 37.02, ... 45.81 45.85 45.87 45.895 45.91 45.925 48.00]'; sigipref=[ 1 1 1 1 1 1 1, ... 2 2 2 2,... 3 3 3 3 3 3 3]'; pref=[0 2000 4000]; %DB (Meters for the model case) %REFERENCE LEVELS sigid sigid=[12 12 8]; %AS IN ALISON'S %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% p_ref_lcd=0; p_load=1; Meth={'slope0';'slope0';'slope1';'pfit';'polyfit';'hpolyfit'}; Marker={'b-o' ;'r-^' ;'g-s' ;'m-x' ;'c-d' ;'k-h'}; Secid='a48n'; isec=3; Secid='a24n'; isec=1; Secid='flst'; isec=1; Secid='a36n'; isec=2; datadir=['/data1/ganacho/Hdata/' Secid '/Geovelm/']; %datadir=['/data1/ganacho/Hdata/' 'a24n' '/Geovelm/']; %for flst for imethod=1:6 if imethod==1 p_trig=-1; %-1 for zero velocity (to use with slope0 only) else p_trig=1; end if p_load secid=[Secid Meth{imethod}]; getpdat sd=sw_dist(Slat,Slon,'km'); sdd=[0;cumsum(sd)]; sdp=0.5*(sdd(1:Nstat-1)+sdd(2:Nstat)); end % FIND THE REFERENCE LEVEL DEPTH [rlpres]=find_sig_interface(sigint(sigid(isec)),1, ... pref(sigipref(sigid(isec))),ptemp,psali,Pres,Maxdp(:,1),Pbotp); if p_ref_lcd disp('SET REF. TO THE LAST COMMON DEPTH IF IT WAS UNDERNEATH') ishdp = [1:Nstat-1; 2:Nstat];isw = find( diff(Botp) < 0 ); ishdp([1 2],isw) = ishdp([2 1],isw); limitdep=Botp(ishdp(1,:)); gis_at_bot=find(rlpres>limitdep); rlpres(gis_at_bot)=limitdep(gis_at_bot); else disp('REF TO DEEPEST STATION') end %GET REFERENCE LEVEL VELOCITY VELOCITY gvel_rl = getsigprop(Pres(:,1),gvel,Pbotp,rlpres); % GET RELATIVE VELOCITY Ndep=size(gvel,1); gvel_rel=gvel-ones(Ndep,1)*(gvel_rl'); gidry=find(isnan(gvel)); gvel(gidry)=zeros(size(gidry)); gvel_rel(gidry)=zeros(size(gidry)); A0=mk_A0(p_trig,MPres,Pres,Botp,Slat,Slon,ptemp,Maxd,Maxdp,Npair); Ci=ones(size(ptemp)); %ro/1e3; Ci(gidry)=zeros(size(gidry)); %gvelint=(1e3*sd').*trapz(Pres, gvel_rel/100)/1e6; %(in Sv) geotrans=A0.*Ci.*gvel_rel/100/1e6; gvelint=sum(geotrans); %integrate from East end tottran{imethod}=[fliplr(cumsum(fliplr(gvelint)))]; %tottran=cumsum(gvelint); Tnga{imethod}=sum(gvelint) % if imethod==1 %plots bottom topography clf; a1=axes('position',[.13 .7 .775 .25 ],'box','on'); plot(sdp,-rlpres,'-.');hold on ns=Nstat;ds=[0;cumsum(sw_dist(Slat,Slon,'km'))]; dp=0.5*(ds(1:ns-1)+ds(2:ns));d1=min(-100*ceil(Botp/100)); fill([ds(1);ds(ns);ds(ns:-1:1)], ... [d1;d1;-Botp(ns:-1:1)],[.8, .8, .8]); grid on;axis([0 max(ds) d1 0 ]);set(gca,'xticklabel',''); title([Secid(1:4)]) a2=axes('position',[.13 .11 .775 .55 ],'box','on');hold on grid on end pl(imethod)=plot(dp,tottran{imethod},Marker{imethod},'markersize',3); end %imethod set(gca,'xlim',[0 max(ds)]) legend(pl,'null velocity','constant vel','max slope = 1','plane fit',... 'polynomial','horizontal') title('');ylabel('Sv');xlabel('distance (km)') setlargefig if 0 %FOR REFERENCE AT THE BOTTOM %print -deps /data4/ganacho/SW25/FIGURES/fig_thesis_116.eps %a36N% %print -deps /data4/ganacho/SW25/FIGURES/fig_thesis_124.eps %a48N% %print -deps /data4/ganacho/SW25/FIGURES/fig_thesis_125.eps %a24N% %print -deps /data4/ganacho/SW25/FIGURES/fig_thesis_126.eps %flst% %FOR RERERENCE AT LCD %print -deps /data4/ganacho/SW25/FIGURES/fig_thesis_127.eps %a24N% %print -deps /data4/ganacho/SW25/FIGURES/fig_thesis_128.eps %flst% %print -deps /data4/ganacho/SW25/FIGURES/fig_thesis_129.eps %a36N% %print -deps /data4/ganacho/SW25/FIGURES/fig_thesis_130.eps %a48N% %title([secid(1:4) ' relative transport']) %print -depsc /data4/ganacho/SW25/FIGURES/fig_talk_002.epsc surf(Plon,-Pres,gvel_rel) set(gca,'ylim',[-4000 0]) set(gca,'zlim',[-200 200]) xlabel('longitude');ylabel('depth');zlabel('velocity'); title([secid ' geostrophic velocity']) signature set(gcf,'papero','land') set(gcf, 'Paperposition', [1 1 8 6]) end