BIBLIOGRAPHIE
Centre Canadien de Télédétection. (1998).
Notions fondamentales de télédétection.
CHAPONNIERE A. et al. 2005. A combined high and low spatial
resolution approach for mapping snow covered area in the Atlas Mountains,
International Journal of remote sensing.
CHAPONNIERE A. (2005). Fonctionnement hydrologique d'un
bassin versant montagneux semi- aride, Cas du bassin versant du Rehraya (Haut
Atlas marocain).Thèse de doctorat, Institut National Agronomique
Paris-Grignon.
DE SOLAN. B. et al (2002). Cartographie de l'enneigement par
télédétection à partir d'images SPOT VEGETATION et
Landsat TM: application a l'Atlas marocain. Séminaire international
`Hydrologie nivale en méditerranée'. Beyrouth, Lebon. p
2.
Duchemin. B. et al. (2002). Normalisation of directional
effects in 10-day global syntheses derived from VEGETATION/SPOT. II Validation
of an operational method on actual data sets. Remote Sensing of Environment
81: 101-113.
JUILLERAT M., 2004, Etude hydrologique comparative de 5 bassins
versants du Haut Atlas marocain, mémoire de maîtrise,
université de franche-comté.
HANICH, L et al. (2003). Snow covers mapping using SPOT
VEGETATION with high resolution data: application in the Moroccan Atlas
mountains. Proceedings of IEEE International Geoscience and Remote Sensing
Symposium (IGARSS). Toulouse, France.
IRD - Valérie Rotival, (2003).Communiqué de presse
sur la conférence Internationale "Hydrologie des régions
méditerranéennes et semi-arides".
Klaus Seidel, Jaroslav Martinec, 2004, «Remote sensing
in snow hydrology: Runoff Modelling, effect of climate change»,
Springer-Paraxis, Chichester, UK, 150 p.
LACOMBE J. (2005/2006), TD d'Initialisation au Logiciel ENVI
4.1, ENSAT, INP (Toulouse).
LISSENS G. et al. (2000), Development of a cloud, snow and cloud
shadow mask for VEGETATION imagery. Proceedings of VEGETATION 2000
Symposium, 3-6 April 2000, Lake Maggiore, Italy.
Musy A. (2004). "Cours d' Hydrologie générale
en ligne." EPFL Lausanne.
ORSTOM. (1976), Hydrologie du bassin du Tensift,
Ministère des travaux publics et communications, Direction de
l'Hydraulique, Division des Ressources en Eau.
SAIDI M. et al. (2006). Effet de la morphologie et de
l'exposition sur les ressources en eau superficielle de part et d'autre du Haut
Atlas (Maroc) ; exemple des bassins versants de l'Ourika et du Marghène.
Bulletin de l'Institut Scientifique, Rabat, section Sciences de la
Terre, n°28.
W. GARETH REES. (2006), Remote sensing of snow and ice,
CRC Press Book, London, 285p.
ANNEXES
Programme d'analyse statistique des images
VGT
%Analyse statistiques :
%MOY,MIN,MAX et NOMBRE D'IMAGES VIDES
%ANGLE DE VISEE « VZA » AU DESSUS DE MARRAKECH %Ecris
par Abdelghani BOUDHAR ___Mai 2006
close all, clear all,clc,fclose all,
chemin1='C: \Documents and
Settings\Irimed\Bureau\images_vgt_extraites\2005- 01-11_2005-06-30\Tensift';
FID1=fopen([chemin1, '\liste.txt'], 'r');
FID2=fopen('statistique_2005.txt', 'w'); fprintf(FID2,'%s \n',' REFLECTANCE_B0
PIEMONT
ATLASIQUE VZA à MARRAKECH');
fprintf(FID2,'%s \n',' date min max moy
std Pixels_vides');
fprintf(FID2, '%s \n','
=' );
for i=1: 360%'eof'
tline1 = fgetl(FID1);
if size(tline1,2)>1
date=tline1 (1:8)
name=tline1 (1:12);
name_ref=strcat(chemin1, '\',name, '_ref');
%ouverture des images VGT avec les 4 bandes
VGT = multibandread(name _ref, [191,280,4], 'int16=>double',0,
'bsq', 'ieeele', {'Band', 'Direct', [1 2 3 4] });
B0=VGT(:, :,1) ;B0=min(B0,2000) ';% elimination des dc
>2000
B2=VGT(:, :,2) ;B2=min(B2,2000) '; B3=VGT(:, :,3)
;B3=min(B3,2000) '; MIR=VGT(:,:,4);MIR=min(MIR,2000)';
%calcule du moyen min, max et l'écart type sur la zone
masquée
;Reflectance en % = CN*100/2000
moyenne _b0=mean2 (B0 (64:101,41:176)) /20; minimum _b0=min (min
(B0 (64:101,41:176))) /20; maximum _b0=max(max(B0 (64:101,41:176))) /20; ecart
_type _b0= std2(B0(64:101,41:176))/20; %nbre de pixel contenant 0 (pas de
données) NB_pixel_vide=sum(sum(B0==0));
%calcul de l'angle VZA à Marrakech en degré=CN*0.5
name_ang=strcat (chemin1, '\',name, '_ang');
VZA= multibandread (name _ang, [24,35,4], 'uint8=>uint8' ,0,
'bsq', 'ieee-
le', { 'Band','Direct', [1] });
ANGlE=VZA(12,21)/2 ;
tline1 = fgetl(FID1);tline1 = fgetl(FID1);%sauter les deux noms:
sm et
ang dans la liste
fprintf(FID2,'%s \n',' ');
fprintf(FID2, '%s \t', '',date);
fprintf(FID2, '%6.2f \t', minimum _b0);fprintf (FID2, '%6.2f
\t',
maximum _b0);fprintf (FID2, '%6.2f \t', moyenne _b0);fprintf
(FID2, '%6.2f \t', ecart_type_b0);
fprintf(FID2, '%6.2f \t', NB_pixel_vide);fprintf(FID2, '%6.2f
\t', ANGlE); end
end
fclose (FID1);
fclose (FID2);
%FIN
Programme de calcul de la surface neigeuse sur le Haut
Atlas et par tranche d'altitude
%Calcul de la couverture neigeuse de tout le Haut Atlas de
Marrakech et par %tranche d'altitude
%A.BOUDHAR ____le 13 juin 2006
close all, clear all,clc,fclose all,
SI100=1000;%dc=1000 pour SI100;
%ouverture image sol nu
SOL_NU_S I0=multibandread ( 'sol
nu _SI0', [191,280,1], 'int16=>double',0, 'bsq', 'ieeele', {
'Band','Direct', [1] });
%ouverture d'MNT
MNT =
multibandread('MNT _atlas _3.214s', [1910,2800,1],
'int16=>double' ,0, 'bsq', 'ieee -le', { 'Band','Direct', [1] });
%MNT réechantillonné à 1/112 degrée
et de même taille que les images vgt
for i=1:191
for j=1:280
MNT(i,j)=mean2(MNT((i*10)-9:i*10,(j*10)-9:j*10));
end
end
MNT _vgt=MNT (1: 191, 1: 280);
clear MNT;
%ouverture de carte des pentes
Pente =
multibandread('slope_3.214s', [1910,2800,1],
'float=>double',0, 'bsq', 'ieee-
le', { 'Band','Direct', [1] });
%pente réechantillonnée à 1/112
degrée et de même taille que les images vgt for i=1:191
for j=1:280
Pente(i,j)=mean2(Pente((i*10)-9:i*10, (j*10)-9:j*10));
end
end
pente _vgt=Pente (1: 191, 1: 280);
clear Pente;
%lecture et puis ouverture de la liste des images vgt une par
une
chemin1='C: \Documents and Settings\Irimed\Bureau\images
_selectionnées _2002- 2005\2002-2003';
%chemin1= [pwd, '\2003'];
FID1=fopen([chemin1, '\liste.txt'], 'r');
%FID2=fopen('Surface _neige _altitude _2004-2005.txt', 'w');
chemin2= [pwd, '\Surface_neigeuse\2002-2003'];
FID2=fopen([chemin2, '\surface _par _altitudes.txt'], 'w');
for i=1:'eof'
tline1 = fgetl(FID1);
if size(tline1,2)>1
date=tline1 (1:8)
name=tline1 (1:12);
name_ref=strcat(chemin1, '\',name, '_ref');
%ouverture des images VGT avec les 4 bandes
VGT = multibandread(name _ref, [191,280,4], 'int16=>double',0,
'bsq', 'ieeele', {'Band', 'Direct', [1 2 3 4] });
B0=VGT (:, : , 1); B2=VGT (:, : , 2); B3=VGT (:, : , 3);
MIR=VGT (:, : , 4);
%COUVERTURE NEIGEUSE en pourcent sur tout le Haut Atlas de
Marrakech surf=strcat ( 'surf_' ,date);
SI=( (B0+B2) /2) -MIR;
MSI=SI-((SOL_NU_SI0).*((SI100-SI) ./(SI100-SOL_NU_SI0)));
pourcent=100*(1-0.492132*exp(-0.6901528*((MSI/100)+1))) .^20; pourcent1=max
(0,pourcent);
pourcent2=min (100,pourcent);
Surface=(pourcent2.*0.84)/100; %1pix===>1/112degré=
0. 99*0 .85=0.84km2
%altitude<1000m ====>pixels sans neige
[L,C]=find(MNT_vgt<1000); for i=1:size(L)
Surface(L(i) ,C(i))=0;
end
% % %correction bruits
% % Application des seuils radiométrique pour q'un pixel
est déclaré neigeuse
for i=1:191
for j=1:280
%
if MIR(i,j)<481 & SI(i,j)>77 %&
((B0(i,j)-MIR(i,j))./(B0(i,j)- MIR(i,j)))*1000>=87&B2(i,j)>=400%&
((B0(i,j)-
B3 (i, j) ) . /B0 (i,j) +B3 (i, j) ) *1000>-773 ;
else end
Surface (i, j ) =Surface (i, j); Surface (i, j ) =0;
end end
%correction de l'effet de pente
surface _corr=(Surface)
./cos((((atan((pente_vgt)/100))*180)/pi)*0.0175);%en
Km2,pi/180=0.0175==>convertion du radian au degrée
% figure;imagesc(pourcent2)
%Surface par tranche d'altitude
Surface _totale _atlas=sum (sum (surface_corr)) fprintf(FID2,'%s
\n',' ');
fprintf(FID2, '%s \n',date, 'Surface_totale en Km2=
');
fprintf(FID2, '%6.2f \n' ,Surface _totale _atlas);
fprintf(FID2, '%s \t','ALTITUDE', 'Surf _Moy',
'Surf_max','Surf_totale');
for i=1:8
switch i
case 1, Alt1=1000;Alt2=1400; case 2, Alt1=1400;Alt2=1800; case
3, Alt1=1800;Alt2=2200; case 4, Alt1=2200;Alt2=2600; case 5,
Alt1=2600;Alt2=3000; case 6, Alt1=3000;Alt2=3400; case 7, Alt1=3400;Alt2=3800;
case 8, Alt1=3800;Alt2=4200;
end
%Sur _totale _tranches=(size (f ind(MNT _vgt>Alt1&MNT
_vgt<Alt2) ) /2) *0.84;
Surface _neige _moy=mean (surface_corr (f ind (MNT _vgt>Alt1
&MNT _vgt<Alt2))); Surface _neige _max=max (surface_corr (f ind (MNT
_vgt>Alt1 &MNT _vgt<Alt2)));
Surface_neige_totale=sum(sum(surface_corr(find(MNT_vgt>Alt1&MNT_vgt<Alt2))))
;
ALTITUDE= (Alt1+Alt2) /2;
% exportation des resultats dans un fichier txt
fprintf(FID2,'%s \n',' ');
fprintf(FID2, '%d \t',ALTITUDE);
fprintf(FID2, '%6.2f \t' ,Surface_neige_moy); fprintf(FID2,
'%6.2f \t' ,Surface_neige_max); fprintf(FID2, '%6.2f \t' ,Surface _neige
_totale);
end
fprintf(FID2,'%s \n',' ');
%ecriture des images de surface dans un fichier bsq lu par ENVI
%multibandwrite (Surface, surf, 'bsq');
multibandwrite (surface_corr, [chemin2, '\',surf], 'bsq');
%multibandwrite (surface_corr, surf, 'bsq');
Hdr1 = fopen(strcat([chemin2,'\',surf,'.hdr']),'w');
fprintf(Hdr1,'%s \n','ENVI');
fprintf(Hdr1, '%s \n', 'description = {'); fprintf(Hdr1, '%s \n',
'Exported from MATLAB}');
fprintf(Hdr1, '%s \n', 'samples = 280'); fprintf(Hdr1, '%s \n',
'lines = 191'); fprintf(Hdr1, '%s \n', 'bands = 1'); fprintf(Hdr1,'%s
\n','header offset = 0'); fprintf(Hdr1,'%s \n','data type = 5') fprintf(Hdr1,
'%s \n', 'interleave = bsq'); fprintf(Hdr1,'%s \n','map info = {');
fprintf(Hdr1, '%s \n', 'Geographic Lat/Lon, 1.0000, 1.0000,
-9.50446429,
32.21446429, 8.9285714000e-003, 8.9285714000e-003, WGS 1984,
units=Degrees }');
fprintf(Hdr1, '%s \n', 'band names = {Surface}');
fclose (Hdr1);
end
end
fclose (FID1);
fclose (FID2);
%FIN
Programme de calcul de la surface neigeuse par
exposition
%Calcul de la couverture neigeuse de tout le Haut Atlas de
Marrakech par %Exposition Nord et Sud
%A.BOUDHAR ____le 13 juin 2006
close all, clear all,clc,fclose all,
SI100=1000;%dc=1000 pour SI100;
%ouverture image sol nu
SOL_NU_S I0=multibandread ( 'sol
nu _SI0', [191,280,1], 'int16=>double',0, 'bsq', 'ieeele', {
'Band', 'Direct', [1] });
%OUVERTURE D'MNT
MNT =
multibandread('MNT _atlas _3.214s', [1910,2800,1],
'int16=>double' ,0, 'bsq', 'ieee -le', { 'Band','Direct', [1] });
%MNT réechantillonné à 1/112 degrée
et de même taille que les images vgt
for i=1:191
for j=1:280 MNT(i,j)=mean2(MNT((i*10)-9:i*10,(j*10)-9:j*10));
end
end
MNT _vgt=MNT (1: 191, 1: 280);
clear MNT;
%OUVERTURE DE CARTE DES PENTES
Pente =
multibandread('slope_3.214s', [1910,2800,1],
'float=>double',0, 'bsq', 'ieee-
le', { 'Band','Direct', [1] });
%pente réechantillonnée à 1/112
degrée et de même taille que les images vgt for i=1:191
for j=1:280
Pente(i,j)=mean2(Pente((i*10)-9:i*10, (j*10)-9:j*10));
end
end
pente _vgt=Pente (1: 191, 1: 280);
clear Pente;
%OUVERTURE DE CARTE DES ORIENTATIONS
Aspect =
multibandread('Aspect _3.2s _vgt', [1910,2800,1],
'float=>double',0, 'bsq', 'ieeele', { 'Band','Direct', [1] });
%carte orientations réechantillonnée à 1/112
degrée et de même taille que
les images vgt for i=1:191
for j=1:280
Aspect _max(i,j)=max(max(Aspect((i*10)-9:i*10,
(j*10)-9:j*10)));
Aspect _min(i,j)=min(min(Aspect((i*10)-9:i*10,
(j*10)-9:j*10)));
end
end
Aspect _max _vgt=Aspect _max(1 :191,1:280);
Aspect _min _vgt=Aspect _min(1 :191,1:280);
clear Aspect; for i=1:191
for j=1:280
if Aspect _max _vgt(i,j)>315 &
Aspect_min_vgt(i,j)<45
Aspect _vgt(i,j)=1;% Versant Nord
else if Aspect_max_vgt(i,j)<225 & Aspect _min
_vgt(i,j)>135 Aspect_vgt(i,j)=2;% Versant Sud
else
Aspect_vgt (i, j) =0;
end
end
end
end
Aspect_vgt;
%lecture et puis ouverture de la liste des images vgt une par une
chemin1='C: \Documents and Settings\Irimed\Bureau\images _selectionnées
_2002- 2005\2004-2005';
%chemin1= [pwd, '\2003'];
FID1=fopen([chemin1, '\liste.txt'], 'r');
chemin2= [pwd, '\Surface_neigeuse\2004-2005'];
FID2=fopen( [chemin2, '\surface_par_oriatation.txt'], 'w');
for i=1:'eof'
tline1 = fgetl(FID1);
if size(tline1,2)>1 date=tline1 (1:8) name=tline1 (1:12);
name_ref=strcat(chemin1, '\',name, '_ref');
%ouverture des images VGT avec les 4 bandes
VGT = multibandread(name _ref, [191,280,4], 'int16=>double',0,
'bsq', 'ieeele', {'Band', 'Direct', [1 2 3 4] });
B0=VGT (:, : , 1); B2=VGT (:, : , 2); B3=VGT (:, : , 3);
MIR=VGT (:, : , 4);
%COUVERTURE NEIGEUSE EN POURCENT SUR TOUT LE HAUT ATLAS DE
MARRAKECH
surf=strcat ( 'surf_' ,date);
SI=( (B0+B2) /2) -MIR;
MSI=SI-((SOL _NU _SI0).*((SI100-SI) ./(SI100-SOL _NU _SI0)));
pourcent=100*(1-0.492132*exp(-0.6901528*((MSI/100)+1))) .^20;
pourcent1=max (0,pourcent);
pourcent2=min (100,pourcent);
Surface=(pourcent2.*0.84)/100; %1pix===>1/112degré=
0.99*0 .85=0.84km2
%altitude<1000m ====>pixels sans neige
[L,C]=find(MNT_vgt<1000);
for i=1:size(L)
Surface(L(i) ,C(i))=0;
end
% % %CORRECTION BRUITS
% % Application des 3 seuils radiométrique pour q'un pixel
est déclaré neigeuse
for i=1:191
for j=1:280
if MIR(i,j)<481 & SI(i,j)>77 %&
((B0(i,j)-MIR(i,j))./(B0(i,j)-
MIR(i,j)))*1000>=87%&B2(i,j)>=400%& ((B0(i,j)-
B3 (i, j) ) . /B0 (i,j) +B3 (i, j) ) *1000>-773 ;
Surface (i, j ) =Surface (i, j);
else
Surface (i, j ) =0;
end
end end
%CORRECTION DE L'EFFET DE PENTE
surface _corr=(Surface)
./cos((((atan((pente_vgt)/100))*180)/pi)*0.0175);%en
Km2,pi/180=0.0175==>convertion du radian au degrée
% figure;imagesc(pourcent2)
%SURFACE PAR EXPOSITION
Surface _totale _atlas=sum (sum (surface_corr))
fprintf(FID2,'%s \n',' ');
fprintf(FID2, '%s \n',date, 'Surface _totale en Km2=
'); fprintf(FID2, '%6.2f \n' ,Surface _totale _atlas); fprintf(FID2, '%s \t',
'Nord', 'Sud');
%Sur _totale _tranches=(size (find(MNT _vgt>Alt1&MNT
_vgt<Alt2) ) /2) *0.84; Surface _neige _Nord=sum (sum (surface_corr (f ind
(Aspect _vgt==1)))); Surface _neige _Sud=sum(sum(surface _corr(find(Aspect
_vgt==2))));
% EXPORTATION DES RESULTATS DANS UN FICHIER TXT
fprintf(FID2,'%s \n',' ');
fprintf(FID2, '%6.2f \t' ,Surface _neige _Nord);
fprintf(FID2, '%6.2f \t' ,Surface _neige _Sud);
end
end
fclose (FID1); fclose (FID2); %FIN
|