Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Chris follett #13

Open
wants to merge 2 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 3 additions & 0 deletions DataSet From DARWIN/AddRrsCorrection/README.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@
# Add RrsCorrection.jl will open the DARWIN Dataset from the oceanclustering repository (stored here in this folder ) and apply
# the zenith angle correction to it

10 changes: 10 additions & 0 deletions DataSet From DARWIN/AddRrsCorrection/RrsCorrection.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,10 @@
using Plots, Distributions, NetCDF, NCDatasets, Interpolations

filename="HalfDegree.nc"
varnames=["Rirr412", "Rirr443", "Rirr490", "Rirr510", "Rirr555", "Rirr670"]
for n in 1:length(varnames)
Rirr=ncread(filename, varnames[n])
tmp=Rirr/3
Rrs0=(0.52*tmp)./(1.0 .-1.7*tmp)
ncwrite(Rrs0,"HalfDegree.nc",varnames[n])
end
12 changes: 6 additions & 6 deletions DataSet From DARWIN/Main_MakeGridDarwin2.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,10 +3,10 @@

include("MakeGridDarwin2.jl")

if ~isfile("MasterHalf.nc")
#if ~isfile("MasterHalf.nc")
MakeMasterHalf()
end
realnames=["Latitude","Longitude","Month","THETA", "SALT", "Chl", "Rirr412", "Rirr443", "Rirr490", "Rirr510", "Rirr555","Rirr670","MXLDEPTH","Bottom Depth", "Wind Speed", "TKE","PAR","Euphotic Depth"]
for n in [.5 1 2 4]
MakeDarwinGrid(n,realnames)
end
#end
# realnames=["Latitude","Longitude","Month","THETA", "SALT", "Chl", "Rirr412", "Rirr443", "Rirr490", "Rirr510", "Rirr555","Rirr670","MXLDEPTH","Bottom Depth", "Wind Speed", "TKE","PAR","Euphotic Depth"]
# for n in [.5 1 2 4]
# MakeDarwinGrid(n,realnames)
# end
19 changes: 12 additions & 7 deletions DataSet From DARWIN/MakeGridDarwin2.jl
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,8 @@ function MakeMasterHalf()
#This is the starting folder for the Climatology
#name="C:\\Users\\folle\\Dropbox (MIT)\\201805-CBIOMES-climatology\\interp"
name="C:\\Users\\folle\\Documents\\CBIOMES\\interp"
#name="/nfs/cnhlab001/cnh/cbiomes/201805-CBIOMES-climatology/interp/"

#This is the filename for the .5 degree monthly climatology
filename="HalfDegree.nc"
#These are the variable names which we would like on the grid
Expand Down Expand Up @@ -121,7 +123,10 @@ function MakeMasterHalf()
OCCIbands=["Rirr412", "Rirr443", "Rirr490", "Rirr510", "Rirr555",
"Rirr670"]
#export OCCIbands
if ~isfile("InterpedWavebands.nc")
if isfile("InterpedWavebands3.nc")
rm("InterpedWavebands3.nc")
end
if ~isfile("InterpedWavebands3.nc")
for n in 1:6
nccreate("InterpedWavebands3.nc", OCCIbands[n],"Lon",collect(-179.75:.5:179.75),"Lat",collect(-89.75:.5:89.75),"Month",collect(1:12))
ncwrite(wbcci[:,:,:,n],"InterpedWavebands3.nc", OCCIbands[n])
Expand Down Expand Up @@ -191,7 +196,7 @@ function MakeMasterHalf()
if isfile("EuphoticDepth.nc")
rm("EuphoticDepth.nc")
end
nccreate("EuphoticDepth.nc", "Euphotic Depth","Lon",collect(.5:.5:360),"Lat",collect(-89.75:.5:89.75),"Month",collect(1:12))
nccreate("EuphoticDepth.nc", "Euphotic Depth","Lon",collect(-179.75:.5:179.75),"Lat",collect(-89.75:.5:89.75),"Month",collect(1:12))
ncwrite(Eudepth,"EuphoticDepth.nc", "Euphotic Depth")
#lon lat depth t

Expand All @@ -207,8 +212,8 @@ function MakeMasterHalf()
#We have everything except PAR, Ekman pumping, Kd, and Euphotic Depth.
vars=["THETA", "SALT", "Chl", "Rirr412", "Rirr443", "Rirr490", "Rirr510", "Rirr555",
"Rirr670","MXLDEPTH","Bottom Depth", "EXFwspee", "GGL90TKE","PAR","Euphotic Depth"]
realnames=["THETA", "SALT", "Chl", "Rirr412", "Rirr443", "Rirr490", "Rirr510", "Rirr555",
"Rirr670","MXLDEPTH","Bottom Depth", "Wind Speed", "TKE","PAR","Euphotic Depth"]
realnames=["SST", "SALT", "Chl", "Rrs412", "Rrs443", "Rrs490", "Rrs510", "Rrs555",
"Rrs670","MLD","Bathymetry", "wind", "EKE","PAR","Euphotic Depth"]
#export vars
Tg=Array{Float64,4}(undef,720,360,12,length(vars))
for n in [1 2 3]
Expand Down Expand Up @@ -255,7 +260,7 @@ function MakeDarwinGrid(deg,names)
for n in 1:length(names)
if sum(k.==names[n])>=1
vartmp=ncread("MasterHalf.nc",names[n])
nccreate("DarwinGrid"*string(deg)*".nc", names[n],"Londim",collect(.5:deg:360),"Latdim",collect(-89.75:deg:89.75),"Monthdim",collect(1:12))
nccreate("DarwinGrid"*string(deg)*".nc", names[n],"Londim",collect(-179.75:deg:179.75),"Latdim",collect(-89.75:deg:89.75),"Monthdim",collect(1:12))
#println(size(vartmp[1:Int(deg/.5):end,1:Int(deg/.5):end,:]))
#println(size(vartmp))
tmp=vartmp[1:Int(deg/.5):end,1:Int(deg/.5):end,:];
Expand All @@ -268,15 +273,15 @@ function MakeDarwinGrid(deg,names)
ncwrite(tmp,"DarwinGrid"*string(deg)*".nc", names[n])
end
if sum(k.==names[n])==0
tmp=zeros(length(collect(.5:deg:360)), length(collect(-89.75:deg:89.75)),12)
tmp=zeros(length(collect(-179.75:deg:179.75)), length(collect(-89.75:deg:89.75)),12)
if names[n]=="Latitude"
Lat=collect(-89.75:deg:89.75);
for m in 1:length(Lat)
#fill!(HalfDegree[:,n,:,1],Lat[n])
tmp2=fill(Lat[m],size(tmp[:,m,:]))
tmp[:,m,:]=tmp2
end
nccreate("DarwinGrid"*string(deg)*".nc", names[n],"Londim",collect(.5:deg:360),"Latdim",collect(-89.75:deg:89.75),"Monthdim",collect(1:12))
nccreate("DarwinGrid"*string(deg)*".nc", names[n],"Londim",collect(-179.75:deg:179.75),"Latdim",collect(-89.75:deg:89.75),"Monthdim",collect(1:12))
ncwrite(tmp,"DarwinGrid"*string(deg)*".nc", names[n])
sztmp=size(tmp)
end
Expand Down
22 changes: 11 additions & 11 deletions DataSet From DARWIN/README.md
Original file line number Diff line number Diff line change
Expand Up @@ -14,19 +14,19 @@ This folder contains code written in Julia to generate gridded data for the foll
- "Lat" (degrees)
- "Lon" (degrees)
- "Month" (month)
- "THETA" (degrees C)
- "SST" (degrees C)
- "SALT" (psu)
- "Chl" (mg/l)
- "Rirr412" (irradiance reflectance for waveband at wavelength)
- "Rirr443" (irradiance reflectance for waveband at wavelength)
- "Rirr490" (irradiance reflectance for waveband at wavelength)
- "Rirr510" (irradiance reflectance for waveband at wavelength)
- "Rirr555" (irradiance reflectance for waveband at wavelength)
- "Rirr670" (irradiance reflectance for waveband at wavelength)
- "MXLDEPTH" (depth, meters)
- "Bottom Depth" (depth, meters)
- "Wind Speed" (m/s)
- "TKE" (m^2/s^2) is GGL90 sub-grid turbulent kinetic energy
- "Rrs412" (irradiance reflectance for waveband at wavelength)
- "Rrs443" (irradiance reflectance for waveband at wavelength)
- "Rrs490" (irradiance reflectance for waveband at wavelength)
- "Rrs510" (irradiance reflectance for waveband at wavelength)
- "Rrs555" (irradiance reflectance for waveband at wavelength)
- "Rrs670" (irradiance reflectance for waveband at wavelength)
- "MLD" (depth, meters)
- "Bathymetry" (depth, meters)
- "wind" (m/s)
- "EKE" (m^2/s^2) is GGL90 sub-grid turbulent kinetic energy
- "PAR" ()
- "Euphotic Depth" Depth (m) where light is 1% of the surface values

Expand Down
Binary file added Try Jackson/J17.nc
Binary file not shown.
92 changes: 92 additions & 0 deletions Try Jackson/JacksonModelRemote.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,92 @@
#This script will attempt to apply the Jackson classifier to the DARWIN reflectance data

using Plots, Distributions, NetCDF, NCDatasets, Interpolations


##This is Gael's Code, taken from 05-OceanColourClassifiers.ipynb##

function fcm(M,Sinv,Rrs)
f=Array{Any,1}(undef,length(M))
for ii=1:length(M)
X=vec(Rrs)-M[ii]
Z=transpose(X)*Sinv[ii]*X
f[ii]=ccdf(Chisq(6),Z)
end
f
end

function getJackson(file)
#Jackson et al 2017:
tmpM = ncread("J17.nc", "cluster_means")
tmpSinv = ncread("J17.nc", "inverse_covariance")
wv_cci=[412, 443, 490, 510, 555, 670]
M=Array{Any,1}(undef,14)
Sinv=Array{Any,1}(undef,14)
for ii=1:length(M)
M[ii]=vec(tmpM[ii,:])
Sinv[ii]=tmpSinv[1:6,1:6,ii]
end

J17=Dict("M" => M, "Sinv" => Sinv, "S" => inv.(Sinv))
#plot(wv_cci,M,w=3); xlabel!("nm"); ylabel!("Rrs")


ds = Dataset(file)

Rrs_412=ds["Rrs412"]
Rrs_443=ds["Rrs443"]
Rrs_490=ds["Rrs490"]
Rrs_510=ds["Rrs510"]
Rrs_555=ds["Rrs555"]
Rrs_670=ds["Rrs670"]

Rrs_412=Rrs_412[1:4:720,1:4:360,1:4:12]
Rrs_443=Rrs_443[1:4:720,1:4:360,1:4:12]
Rrs_490=Rrs_490[1:4:720,1:4:360,1:4:12]
Rrs_510=Rrs_510[1:4:720,1:4:360,1:4:12]
Rrs_555=Rrs_555[1:4:720,1:4:360,1:4:12]
Rrs_670=Rrs_670[1:4:720,1:4:360,1:4:12]

tmp=fill(false,size(Rrs_412))
for ii in eachindex(Rrs_412)
!ismissing(Rrs_412[ii]) ? tmp[ii]=!ismissing(Rrs_443[ii].*Rrs_490[ii].*Rrs_510[ii].*Rrs_555[ii].*Rrs_670[ii]) : nothing
end
ii=findall(tmp)
sz=size(Rrs_412)
#export ii, sz, Rrs_412, Rrs_443, Rrs_490, Rrs_510, Rrs_555, Rrs_670
mbrshp=Array{Float64,4}(undef,(sz[1],sz[2],sz[3],14));
for jj=1:length(ii);
kk=ii[jj]
Rrs_tmp=[Rrs_412[kk] Rrs_443[kk] Rrs_490[kk] Rrs_510[kk] Rrs_555[kk] Rrs_670[kk]]
mbrshp[kk[1],kk[2],kk[3],:]=fcm(J17["M"],J17["Sinv"],Rrs_tmp)
println(kk)
end
return mbrshp
end
#path="https://rsg.pml.ac.uk/shared_files/brj/CBIOMES_ecoregions/ver_0_2_5"
path="https://rsg.pml.ac.uk/shared_files/brj/CBIOMES_ecoregions/ver_0_2_5/"
filegeo="gridded_geospatial_montly_clim_360_720_ver_0_2.nc"
filedarwin="gridded_darwin_montly_clim_360_720_ver_0_2_2.nc"
#run(`wget $path/$file`)
if ~isfile(filegeo)
download(string(path, filegeo), filegeo)
end
if ~isfile(filedarwin)
download(string(path, filedarwin), filedarwin)
end

#varnames=["Rirr412", "Rirr443", "Rirr490", "Rirr510", "Rirr555", "Rirr670"]
varnames=["Rrs412", "Rrs443", "Rrs490", "Rrs510", "Rrs555", "Rrs670"]

for n in 1:length(varnames)
Rirr=ncread(filedarwin, varnames[n])
tmp=Rirr/3
Rrs0=(0.52*tmp)./(1.0 .-1.7*tmp)
ncwrite(Rrs0,filedarwin,varnames[n])
#renameVar(filedarwin,varnames[n],newnames[n])
end

DarwinMembers=getJackson(filedarwin)
GeoMembers=getJackson(filegeo)
export DarwinMembers
export GeoMembers
63 changes: 63 additions & 0 deletions Try Jackson/M09.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,63 @@

#Moore et al 2009 mean vectors:
M=Array{Any,1}(undef,8)
M[1]=vec([0.0234 0.0192 0.0129 0.0075 0.0031 0.0002])
M[2]=vec([0.0162 0.0141 0.0112 0.0073 0.0034 0.0002])
M[3]=vec([0.0107 0.0098 0.0092 0.0070 0.0039 0.0003])
M[4]=vec([0.0065 0.0064 0.0070 0.0064 0.0048 0.0006])
M[5]=vec([0.0033 0.0034 0.0042 0.0042 0.0043 0.0009])
M[6]=vec([0.0064 0.0074 0.0105 0.0116 0.0140 0.0041])
M[7]=vec([0.0121 0.0140 0.0192 0.0204 0.0231 0.0084])
M[8]=vec([0.0184 0.0230 0.0333 0.0359 0.0409 0.0137])

#Moore et al 2009 covariance matrices:
S=Array{Any,1}(undef,8)
S[1]=[0.00000959 0.00000556 0.00000138 -0.00000034 -0.00000024 -0.00000003;
0.00000556 0.00000493 0.00000193 0.00000060 0.00000023 0.00000001;
0.00000138 0.00000193 0.00000282 0.00000223 0.00000119 0.00000007;
-0.00000034 0.00000060 0.00000223 0.00000232 0.00000119 0.00000007;
-0.00000024 0.00000023 0.00000119 0.00000119 0.00000071 0.00000005;
-0.00000003 0.00000001 0.00000007 0.00000007 0.00000005 0.00000001]
S[2]=[0.00000346 0.00000186 -0.00000011 -0.00000060 -0.00000062 -0.00000005;
0.00000186 0.00000228 0.00000086 0.00000033 -0.00000007 0.00000003;
-0.00000011 0.00000086 0.00000231 0.00000221 0.00000145 0.00000023;
-0.00000060 0.00000033 0.00000221 0.00000266 0.00000191 0.00000034;
-0.00000062 -0.00000007 0.00000145 0.00000191 0.00000175 0.00000041;
-0.00000005 0.00000003 0.00000023 0.00000034 0.00000041 0.00000021]
S[3]=[0.00000241 0.00000144 0.00000035 -0.00000031 -0.00000063 -0.00000006;
0.00000144 0.00000138 0.00000076 0.00000015 -0.00000021 -0.00000001;
0.00000035 0.00000076 0.00000161 0.00000156 0.00000121 0.00000016;
-0.00000031 0.00000015 0.00000156 0.00000227 0.00000209 0.00000031;
-0.00000063 -0.00000021 0.00000121 0.00000209 0.00000225 0.00000037;
-0.00000006 -0.00000001 0.00000016 0.00000031 0.00000037 0.00000013]
S[4]=[0.00000166 0.00000091 0.00000034 -0.00000009 -0.00000080 -0.00000025;
0.00000091 0.00000097 0.00000071 0.00000025 -0.00000041 -0.00000015;
0.00000034 0.00000071 0.00000118 0.00000103 0.00000072 0.00000003;
-0.00000009 0.00000025 0.00000103 0.00000137 0.00000162 0.00000025;
-0.00000080 -0.00000041 0.00000072 0.00000162 0.00000290 0.00000065;
-0.00000025 -0.00000015 0.00000003 0.00000025 0.00000065 0.00000050]
S[5]=[0.00000178 0.00000132 0.00000104 0.00000081 0.00000018 -0.00000014;
0.00000132 0.00000127 0.00000121 0.00000099 0.00000034 -0.00000010;
0.00000104 0.00000121 0.00000150 0.00000142 0.00000110 0.00000013;
0.00000081 0.00000099 0.00000142 0.00000158 0.00000177 0.00000042;
0.00000018 0.00000034 0.00000110 0.00000177 0.00000351 0.00000131;
-0.00000014 -0.00000010 0.00000013 0.00000042 0.00000131 0.00000081]
S[6]=[0.00000715 0.00000586 0.00000409 0.00000292 0.00000005 -0.00000075;
0.00000586 0.00000589 0.00000520 0.00000398 0.00000027 -0.00000114;
0.00000409 0.00000520 0.00000634 0.00000541 0.00000188 -0.00000097;
0.00000292 0.00000398 0.00000541 0.00000528 0.00000392 0.00000070;
0.00000005 0.00000027 0.00000188 0.00000392 0.00000995 0.00000657;
-0.00000075 -0.00000114 -0.00000097 0.00000070 0.00000657 0.00000819]
S[7]=[0.00002625 0.00001981 0.00001058 0.00000544 -0.00000654 -0.00001122;
0.00001981 0.00001745 0.00001314 0.00000822 -0.00000431 -0.00001228;
0.00001058 0.00001314 0.00001629 0.00001226 0.00000035 -0.00001311;
0.00000544 0.00000822 0.00001226 0.00001170 0.00000742 -0.00000500;
-0.00000654 -0.00000431 0.00000035 0.00000742 0.00002241 0.00001782;
-0.00001122 -0.00001228 -0.00001311 -0.00000500 0.00001782 0.00003987]
S[8]=[0.00001186 0.00001134 0.00001139 0.00000919 0.00000395 -0.00000186;
0.00001134 0.00001484 0.00002034 0.00001907 0.00001531 0.00000087;
0.00001139 0.00002034 0.00003467 0.00003546 0.00003555 0.00000708;
0.00000919 0.00001907 0.00003546 0.00003907 0.00004604 0.00001733;
0.00000395 0.00001531 0.00003555 0.00004604 0.00007306 0.00004953;
-0.00000186 0.00000087 0.00000708 0.00001733 0.00004953 0.00006542];

7 changes: 7 additions & 0 deletions Try Jackson/README.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,7 @@
# TryJackson

This folder contains code that will download the half degree data grids for DARWIN and GEO and attempt to apply the OC-CCI water class algorithm to it.

The important code is JacksonModelRemote.jl. In the current version, the Rirr to Rrs correction is applied to the DARWIN dataset still.


Binary file added Try Jackson/cci_Rrs_490.bin
Binary file not shown.
Binary file added Try Jackson/drwn3_Rirr.bin
Binary file not shown.
Binary file added Try Jackson/lat.bin
Binary file not shown.
Binary file added Try Jackson/lon.bin
Binary file not shown.
85 changes: 85 additions & 0 deletions Try Jackson/tryJackson.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,85 @@
#This script will attempt to apply the Jackson classifier to the DARWIN reflectance data

using Plots, Distributions, NetCDF, NCDatasets, Interpolations


##This is Gael's Code, taken from 05-OceanColourClassifiers.ipynb##

function fcm(M,Sinv,Rrs)
f=Array{Any,1}(undef,length(M))
for ii=1:length(M)
X=vec(Rrs)-M[ii]
Z=transpose(X)*Sinv[ii]*X
f[ii]=ccdf(Chisq(6),Z)
end
f
end

#Jackson et al 2017:
tmpM = ncread("J17.nc", "cluster_means")
tmpSinv = ncread("J17.nc", "inverse_covariance")
wv_cci=[412, 443, 490, 510, 555, 670]
M=Array{Any,1}(undef,14)
Sinv=Array{Any,1}(undef,14)
for ii=1:length(M)
M[ii]=vec(tmpM[ii,:])
Sinv[ii]=tmpSinv[1:6,1:6,ii]
end

J17=Dict("M" => M, "Sinv" => Sinv, "S" => inv.(Sinv))
plot(wv_cci,M,w=3); xlabel!("nm"); ylabel!("Rrs")

#path="https://rsg.pml.ac.uk/shared_files/brj/CBIOMES_ecoregions/ver_0_2_5"
path="https://https://rsg.pml.ac.uk/shared_files/brj/CBIOMES_ecoregions/ver_0_2_5/"
file="gridded_geospatial_montly_clim_360_720_ver_0_2.nc"
file="gridded_darwin_montly_clim_360_720_ver_0_2_2.nc"
#run(`wget $path/$file`)
download(string(path, file), file)
#Using NCDatasets
#Dataset(file)
#fil="C:\Users\folle\Dropbox (MIT)\CBIOMES_TransitionZone\Chris Follett\StatsGroup\ocean_clustering\DataSet From DARWIN\AddRrsCorrection\HalfDegree.nc"
#fil="C:\\Users\\folle\\Dropbox (MIT)\\CBIOMES_TransitionZone\\Chris Follett\\StatsGroup\\ocean_clustering\\DataSet From DARWIN\\AddRrsCorrection\\HalfDegree.nc"
#dir0="../samples/"
#fil="MasterHalf.nc"
##
varnames=["Rirr412", "Rirr443", "Rirr490", "Rirr510", "Rirr555", "Rirr670"]
for n in 1:length(varnames)
Rirr=ncread(filename, varnames[n])
tmp=Rirr/3
Rrs0=(0.52*tmp)./(1.0 .-1.7*tmp)
ncwrite(Rrs0,"HalfDegree.nc",varnames[n])
end
##

ds = Dataset(file)

Rrs_412=ds["Rirr412"]
Rrs_443=ds["Rirr443"]
Rrs_490=ds["Rirr490"]
Rrs_510=ds["Rirr510"]
Rrs_555=ds["Rirr555"]
Rrs_670=ds["Rirr670"]

Rrs_412=Rrs_412[1:4:720,1:4:360,1:4:12]
Rrs_443=Rrs_443[1:4:720,1:4:360,1:4:12]
Rrs_490=Rrs_490[1:4:720,1:4:360,1:4:12]
Rrs_510=Rrs_510[1:4:720,1:4:360,1:4:12]
Rrs_555=Rrs_555[1:4:720,1:4:360,1:4:12]
Rrs_670=Rrs_670[1:4:720,1:4:360,1:4:12]

tmp=fill(false,size(Rrs_412))
for ii in eachindex(Rrs_412)
!ismissing(Rrs_412[ii]) ? tmp[ii]=!ismissing(Rrs_443[ii].*Rrs_490[ii].*Rrs_510[ii].*Rrs_555[ii].*Rrs_670[ii]) : nothing
end
ii=findall(tmp)
sz=size(Rrs_412)
export ii, sz, Rrs_412, Rrs_443, Rrs_490, Rrs_510, Rrs_555, Rrs_670
mbrshp=Array{Float64,4}(undef,(sz[1],sz[2],sz[3],14));
for jj=1:length(ii);
kk=ii[jj]
Rrs_tmp=[Rrs_412[kk] Rrs_443[kk] Rrs_490[kk] Rrs_510[kk] Rrs_555[kk] Rrs_670[kk]]
mbrshp[kk[1],kk[2],kk[3],:]=fcm(J17["M"],J17["Sinv"],Rrs_tmp)
println(kk)
end
export mbrshp
heatmap(mbrshp[:,:,1,10])