Changeset 51

Timestamp:

03/24/08 14:06:03 (10 months ago)

Author:

christophe

Message:

--

Files:

• trunk/mp-geo/doc/doc-mp-geo.tex (modified) (14 diffs)
• trunk/mp-geo/doc/projectionbonne.mp (added)
• trunk/mp-geo/mp-geo.mp (modified) (22 diffs)

trunk/mp-geo/doc/doc-mp-geo.tex

@@ -10 +10 @@
-
-\usepackage{calc,url,subfigure,amsmath}
-
+fourier,
-\input{christ5}
-\usepackage[charter]{mathdesign}
-
-
-\pagestyle{empty}
-
@@ -46 +43 @@
-
-\title{\texttt{mp-geo}}
- 
+
-\date{\today}
-\begin{document}
@@ -59 +56 @@
-\begin{figure}[ht]
-  \centering
-4]{pagepresentation-1
+8]{pagepresentation
-  \caption{Vue de l'Asie ($\theta=90$\degres; $\phi=20$\degres).}
-  \label{Asie}
@@ -117 +114 @@
-\begin{figure}[ht]
-  \centering
-4]{grandslacs-1
+7]{grandslacs
-  \caption{Les grands lacs d'Amérique du Nord ($\theta=-105$\degres;
-    $\phi=45$\degres)}
@@ -155 +152 @@
-  \centering
-  \subfigure[Figure
-4]{payscomplets-1.png}}\hfill\subfigure[Figure \ding{173}]{\includegraphics[scale=0.4]{paysnus-1
+7]{payscomplets.png}}\hfill\subfigure[Figure \ding{173}]{\includegraphics[scale=0.7]{paysnus
-\end{figure}
-Une autre option est celle de colorer uniformément ou pas les
@@ -164 +161 @@
-\begin{figure}[ht]
-  \centering
-4]{mondenoncolore-1
+7]{mondenoncolore
-  \caption{Coloration uniforme ($\theta=105$\degres; $\phi=30$\degres).}
-  \label{unicolore}
@@ -174 +171 @@
-\begin{figure}[ht]
-  \centering
-4]{oceaniemeridien-1
+7]{oceaniemeridien
-  \caption{Méridiens et parallèles ($\theta=120$\degres; $\phi=-20$\degres).}
-  \label{meridien}
@@ -238 +235 @@
-\begin{figure}[ht]
-  \centering
-4
+3
-  \caption{Méthode pour remplir correctement les pays.}
-\end{figure}
@@ -271 +268 @@
-\begin{figure}[ht]
-  \centering
-4]{polesud-1.png
+7]{polesud
-  \caption{L'immensité du pôle sud ($\theta=0$; $\phi=90$\degres).}
-\end{figure}
@@ -280 +277 @@
-\begin{figure}[ht]
-  \centering
-[scale=0.4]{mp-geotrans-1
+{globetransparent
-  \caption{Un peu de transparence.}
-\end{figure}
@@ -342 +339 @@
-\begin{figure}[h]
-  \centering
-[scale=0.4]{mercator-1
+{Mercator
-  \caption{Carte de France selon une projection Mercator.}
-\end{figure}
@@ -379 +376 @@
-\begin{figure}[ht]
-  \centering
-[scale=0.4]{cylindrique-1
+{Cylindrique
-  \caption{Carte de France selon la projection cylindrique.}
-\end{figure}
-\subsection{Projection conique}
-Une projection conique est disponible : il s'agit de la {\em
-  projection conique équidistante tangente}\footnote{D'après
-  \url{http://thierry.hatt.gps.free.fr/projections/images-proj/coniq-equid-lambert.htm}. Pour un point de la sphère terrestre de coordonnées $(\phi;\lambda)$, les formules de projection sont
-  \[\left\{\begin{array}{l}
-      x=\dfrac{\cos\phi\sin(\lambda\sin(45))}{\sin{45}\cos(\phi-45)}\\
-      \\
-      y=\dfrac{\cos\phi\cos(\lambda\sin(45))}{\sin{45}\cos(\phi-45)}\\
-    \end{array}
-  \right.
-  \]
-}; le parallèle de contact étant soit le 45\degres\ Nord soit le
-45\degres\ Sud. On choisit l'hémisphère à représenter par le code
-source
-\begin{verbatim}
-input mp-geo;
-
-figure(-20u,-20u,20u,20u);
-ConiqueH;
-fin;
-
-end
-\end{verbatim}
-pour obtenir
-\begin{figure}[ht]
-  \centering
-  \includegraphics[scale=0.4]{coniqueHN-1.png}
-  \caption{Projection conique - Hémisphère nord.}
-\end{figure}
-ou alors 
-\begin{verbatim}
-input mp-geo;
-
-figure(-20u,-20u,20u,20u);
-ConiqueB;
-fin;
-
-end
-\end{verbatim}
-pour voir l'hémisphère sud.
-\begin{figure}[ht]
-  \centering
-  \includegraphics[scale=0.4]{coniqueBS-1.png}
-  \caption{Projection conique - Hémisphère sud.}
-\end{figure}
-\par\`A noter que cette fois-ci, le maillage n'est pas
-facultatif\footnote{On pourrait peut-être ajouter des graduations}.
-\subsection{Projection de Bonne}
-Cette projection\footnote{Pour un point de la sphère terrestre de
@@ -446 +396 @@
-\begin{figure}[ht]
-  \centering
-[scale=0.4]{bonne-1
+{Bonne
-  \caption{Projection de Bonne}
-\end{figure}
@@ -464 +414 @@
-dans lequel \verb!0,40! représentent $(\lambda_0;\phi_0)$, 1 étant le
-facteur de zoom.
+\subsection{Nouvelle méthode}
+Avec le développement de ces projections, on s'aperçoit qu'elles ont
+toutes le même fonctionnement en commun. Aussi, depuis la version
+\texttt{0.60}, une autre méthode\footnote{Ne vous inquiétez pas :
+  pour des soucis de compatibilité, les autres moyens de projection
+  ont été conservés.} a été implantée.
+\par En effet, il me semble plus simple de procéder de la manière
+suivante :
+\begin{verbatim}
+input mp-geo;
+
+figure(-10u,-10u,10u,10.5u);
+projection:="bonne";
+Projection(0,41,50);
+Lecturevillesp("France");
+fin;
+
+end
+\end{verbatim}
+Ceci afin d'alléger un peu les commandes à connaître. Un autre
+avantage est de choisir simplement le mode de projection choisi.
+\par Pour la création des cartes, il fallait, jusqu'à cette version
+\texttt{0.60}, parcourir tous les points de la planète pour ne garder
+que ceux intéressants.\\Dorénavant, les paramètres \verb!Amnord!,
+\verb!Amsud!, \verb!Amcentrale!, \verb!Caraibes!, \verb!Asie!,
+\verb!Europe!, \verb!Afrique! et \verb!All! font leurs
+apparitions. Par défaut, ils sont tous positionnés à \verb!false! sauf
+\verb!All!. Autrement, le comportement par défaut est celui qui prend
+le plus de temps !
+Aussi, on prendra soin de changer certains de ces positionnements pour
+obtenir des temps de compilation moins grands.
+\subsection{Cartes et échelles}
+Pour pouvoir lire correctement une carte obtenue par
+projection\footnote{Uniquement les projections déjà vues, c'est-à-dire
+  \verb!simple!, \verb!cylindrique!, \verb!mercator! et \verb!bonne!}, il
+m'a semblé nécessaire d'y adjoindre une échelle. On ne peut avoir
+qu'une échelle moyenne; sinon il faudrait déterminer l'échelle de
+chaque parallèle ou méridien.
+\par Aussi, la méthode choisie est la suivante :
+\begin{itemize}
+\item Supposons que la projection soit centrée sur le point
+  $(\lambda_0;\theta_0)$; on détermine alors la longueur {\em réelle}
+  de l'arc de parallèle situé entre les points $(\lambda_0;\theta_0)$
+  et $(\lambda_0;\theta_0+5\mbox{\degres})$;
+\item on détermine ensuite cette même longueur mais {\em sur le dessin};
+\item on termine par une simple relation de proportionnalité.
+\end{itemize}
+On utilisera donc la syntaxe \verb!echelle(2,46.5,50)! pour un calcul
+d'échelle de 1~cm pour 50~km pour une projection centrée sur le point $(46,5;2)$.
+\begin{verbatim}
+input mp-geo;
+
+All:=false;
+Europe:=true;
+
+figure(-10u,-10u,10u,10.5u);
+projection:="bonne";
+Projection(2,46.5,echelle(2,46.5,50));
+fin;
+
+end
+\end{verbatim}
+pour obtenir la figure
+\begin{figure}[ht]
+  \centering
+  \includegraphics{projectionbonne}
+  \caption{Projection de Bonne -- Centrée sur Paris -- \'Echelle $1\,:\,500\,000$.}
+\end{figure}
+\subsection{Projection conique}
+Cette projection est un peu spéciale, je n'ai pas réussi à l'inclure
+avec les autres. La projection conique disponible est la {\em
+  projection conique équidistante tangente}\footnote{D'après
+  \url{http://thierry.hatt.gps.free.fr/projections/images-proj/coniq-equid-lambert.htm}. Pour un point de la sphère terrestre de coordonnées $(\phi;\lambda)$, les formules de projection sont
+  \[\left\{\begin{array}{l}
+      x=\dfrac{\cos\phi\sin(\lambda\sin(45))}{\sin{45}\cos(\phi-45)}\\
+      \\
+      y=\dfrac{\cos\phi\cos(\lambda\sin(45))}{\sin{45}\cos(\phi-45)}\\
+    \end{array}
+  \right.
+  \]
+}; le parallèle de contact étant soit le 45\degres\ Nord soit le
+45\degres\ Sud. On choisit l'hémisphère à représenter par le code
+source
+\begin{verbatim}
+input mp-geo;
+
+figure(-20u,-20u,20u,20u);
+ConiqueH;
+fin;
+
+end
+\end{verbatim}
+pour obtenir
+\begin{figure}[ht]
+  \centering
+  \includegraphics{coniqueHN.png}
+  \caption{Projection conique - Hémisphère nord.}
+\end{figure}
+ou alors 
+\begin{verbatim}
+input mp-geo;
+
+figure(-20u,-20u,20u,20u);
+ConiqueB;
+fin;
+
+end
+\end{verbatim}
+pour voir l'hémisphère sud.
+\begin{figure}[ht]
+  \centering
+  \includegraphics{coniqueBS.png}
+  \caption{Projection conique - Hémisphère sud.}
+\end{figure}
+\par\`A noter que cette fois-ci, le maillage n'est pas
+facultatif\footnote{On pourrait peut-être ajouter des graduations}.
-\newpage
-\section{Historique}
-\begin{description}
+\item[24/03/2008] Version \texttt{0.60} : Regroupement des
+  projections. Ajout d'une échelle sur les cartes obtenues par
+  projection. Ajout des grands fleuves (ou rivières) européens.
-\item[14/03/2008] Version \texttt{0.58} : Ajout de fleuves et lacs
-  supplémentaires.

trunk/mp-geo/mp-geo.mp

@@ -1 +1 @@
-%Package mp-geographie
-%Auteur:Christophe Poulain
-7
-5/03
+5
+8/02
-
-prologues:=2;
@@ -65 +65 @@
-
-def feuille(expr xa,ya,xb,yb) =
+  numeric Xa,Ya;
+  Xa=xa;
+  Ya=ya;
-  feuillet := (xa,ya)--(xa,yb)--(xb,yb)--(xb,ya)--cycle;
-  extra_endfig := "clip currentpicture to feuillet;" & extra_endfig;
@@ -351 +354 @@
-enddef;
-
-
-vardef Lecturelacs=
-  nomfichiermul:=arborescence&"lacs.dat";
@@ -389 +391 @@
-enddef;
-
+vardef Lecturelacssup=
+  nomfichiermul:=arborescence&"lacssup.dat";
+  nblec:=scantokens readfrom nomfichiermul;
+  for p=1 upto nblec:
+    color Coord[],Pays[];
+    numeric ll;
+    ll:=0;
+    nbpts:=scantokens readfrom nomfichiermul;
+    for k=1 upto nbpts:
+      pair latlon;
+      latlon=scantokens readfrom nomfichiermul;
+      Coord[k]=rayon*(cosd(xpart(latlon))*cosd(ypart(latlon)),cosd(xpart(latlon))*sind(ypart(latlon)),sind(xpart(latlon)));
+      if ((xpart(latlon)>phim) and (xpart(latlon)<phip)):
+        if ProduitScalaire(Coord[k]-pte3,Oeil-pte3)>0:
+          ll:=ll+1;
+          Pays[k]=Coord[k];
+        else:
+          Pays[k]=2*Coord[k];
+        fi;
+      else:
+        Pays[k]=2*Coord[k];
+      fi;
+    endfor;
+    path pays;
+    if ll>0:
+      pays=Projgeo(Pays[1])
+      for l=2 upto nbpts:
+        --Projgeo(Pays[l])
+      endfor;
+      fill pays--cycle withcolor cielfonce;
+      draw pays;
+      clip currentpicture to cercles(pte3,pte1,pte3,pte1,pte4);
+    fi;
+  endfor;
+  closefrom nomfichiermul;
+enddef;
+
-vardef Lecturevolcans=
-  drawoptions(withcolor orange);
@@ -427 +466 @@
-      Coord[k]=rayon*(cosd(xpart(latlon/60))*cosd(ypart(latlon/60)),cosd(xpart(latlon/60))*sind(ypart(latlon/60)),sind(xpart(latlon/60)));
-      if ((xpart(latlon/60)>phim) and (xpart(latlon/60)<phip)):
+        if ProduitScalaire(Coord[k]-pte3,Oeil-pte3)>0:
+          ll:=ll+1;
+          Pays[k]=Coord[k];
+        else:
+          Pays[k]=2*Coord[k];
+        fi;
+      else:
+        Pays[k]=2*Coord[k];
+      fi;
+    endfor;
+    path pays;
+    if ll>0:
+      pays=Projgeo(Pays[1])
+      for l=2 upto nbpts:
+        --Projgeo(Pays[l])
+      endfor;
+      draw pays withcolor cielfonce;
+      clip currentpicture to cercles(pte3,pte1,pte3,pte1,pte4);
+    fi;
+  endfor;
+  closefrom nomfichiermul;
+enddef;
+
+vardef Lecturerivieressup=
+  nomfichiermul:=arborescence&"fleuveseurope.dat";
+  nblec:=scantokens readfrom nomfichiermul;
+  for p=1 upto nblec:
+    color Coord[],fond,Pays[];
+    numeric ll;
+    ll:=0;
+    nbpts:=scantokens readfrom nomfichiermul;
+    for k=1 upto nbpts:
+      pair latlon;
+      latlon=scantokens readfrom nomfichiermul;
+      Coord[k]=rayon*(cosd(xpart(latlon))*cosd(ypart(latlon)),cosd(xpart(latlon))*sind(ypart(latlon)),sind(xpart(latlon)));
+      if ((xpart(latlon)>phim) and (xpart(latlon)<phip)):
-        if ProduitScalaire(Coord[k]-pte3,Oeil-pte3)>0:
-          ll:=ll+1;
@@ -594 +669 @@
-  fi;
-  nblec:=scantokens readfrom nomfichiermul;
-o
+n
-    nblec:=nblec-1
-  fi;
@@ -664 +739 @@
-      Coord[p]=mercatorc(ypart(latlon),xpart(latlon));
-    elseif projection="simple":
-xpart(latlon),y
+ypart(latlon),x
-    elseif projection="cylindrique":
-      Coord[p]=cylindriquec(ypart(latlon),xpart(latlon));
@@ -725 +800 @@
-enddef;
-
+vardef Lecturelacspsup=
+  nomfichiermul:=arborescence&"lacssup.dat";
+  nblec:=scantokens readfrom nomfichiermul;
+  for p=1 upto nblec:
+    pair Coord[];
+    nbpts:=scantokens readfrom nomfichiermul;
+    numeric ll;
+    ll=0;
+    for k=1 upto nbpts:
+      pair latlon;
+      latlon=scantokens readfrom nomfichiermul;
+      if projection="mercator":
+        Coord[k]=mercatorc(xpart(latlon),ypart(latlon));
+      elseif projection="simple":
+        Coord[k]=simplec(xpart(latlon),ypart(latlon));
+      elseif projection="cylindrique":
+        Coord[k]=cylindriquec(xpart(latlon),ypart(latlon));
+      elseif projection="bonne":
+        Coord[k]=bonnec(xpart(latlon),ypart(latlon));
+      elseif projection="coniqueh":
+        if (xpart(latlon)>0) or (xpart(latlon)=0):
+          ll:=ll+1;
+          Coord[ll]=coniquec(xpart(latlon),ypart(latlon),45);
+        fi;
+      elseif projection="coniqueb":
+        if (xpart(latlon)<0):
+          ll:=ll+1;
+          Coord[ll]=coniquec(xpart(latlon),ypart(latlon),-45);
+        fi;
+      fi;
+    endfor;
+    path lac;
+    if (projection="mercator") or (projection="cylindrique") or (projection="simple") or (projection="bonne"):
+      ll:=nbpts
+    fi;
+    if ll>0:
+      lac=Coord[1]
+      for l=2 upto ll:
+        --Coord[l]
+      endfor;
+      fill lac--cycle withcolor cielfonce;
+      draw lac;
+    fi;
+  endfor;
+  closefrom nomfichiermul;
+enddef;
+
-vardef Lecturerivieresp=
-  nomfichiermul:=arborescence&"rivieres.dat";
@@ -753 +875 @@
-          ll:=ll+1;
-          Coord[ll]=coniquec(xpart(latlon/60),ypart(latlon/60),-45);
+        fi;
+      fi;
+    endfor;
+    path riv;
+    if (projection="mercator") or (projection="cylindrique") or (projection="simple") or (projection="bonne"):
+      ll:=nbpts;
+    fi;
+    if ll>0:
+      riv=Coord[1]
+      for l=2 upto ll:
+        --Coord[l]
+      endfor;
+      draw riv withcolor cielfonce;
+    fi;
+  endfor;
+  closefrom nomfichiermul;
+enddef;
+
+vardef Lecturerivierespsup=
+  nomfichiermul:=arborescence&"fleuveseurope.dat";
+  nblec:=scantokens readfrom nomfichiermul;
+  for p=1 upto nblec:
+    pair Coord[];
+    nbpts:=scantokens readfrom nomfichiermul;
+    numeric ll;
+    ll:=0;
+    for k=1 upto nbpts:
+      pair latlon;
+      latlon=scantokens readfrom nomfichiermul;
+      if projection="mercator":
+        Coord[k]=mercatorc(xpart(latlon),ypart(latlon));
+      elseif projection="simple":
+        Coord[k]=simplec(xpart(latlon),ypart(latlon));
+      elseif projection="cylindrique":
+        Coord[k]=cylindriquec(xpart(latlon),ypart(latlon));
+      elseif projection="bonne":
+        Coord[k]=bonnec(xpart(latlon),ypart(latlon));
+      elseif projection="coniqueh":
+        if (xpart(latlon)>0) or (xpart(latlon)=0):
+          ll:=ll+1;
+          Coord[ll]=coniquec(xpart(latlon),ypart(latlon),45);
+        fi;
+      elseif projection="coniqueb":
+        if (xpart(latlon)<0):
+          ll:=ll+1;
+          Coord[ll]=coniquec(xpart(latlon),ypart(latlon),-45);
-        fi;
-      fi;
@@ -941 +1109 @@
-enddef;
-
+boolean Echelle;
+Echelle=false;
+
+vardef echelle(expr Th,Ph,long)=%long en km.
+  Echelle:=true;
+  theta:=Th;
+  phi:=Ph;
+  Long:=long;
+  zoom(1);
+  numeric $;
+  if projection="bonne":
+    $=(cm*pi*cosd(Ph)*6340/36)/(long*abs(bonnec(Ph,Th)-bonnec(Ph,Th+5)));
+  elseif projection="mercator":
+    $=(cm*pi*cosd(Ph)*6340/36)/(long*abs(mercatorc(Ph,Th)-mercatorc(Ph,Th+5)));
+  elseif projection="simple":
+    $=(cm*pi*cosd(Ph)*6340/36)/(long*abs(simplec(Ph,Th)-simplec(Ph,Th+5)));
+  elseif projection="cylindrique":
+    $=(cm*pi*cosd(Ph)*6340/36)/(long*abs(cylindriquec(Ph,Th)-cylindriquec(Ph,Th+5)));
+  fi;
+  $
+enddef;
+
+boolean Amsud,Amnord,Amcentrale,Caraibes,Asie,Europe,Afrique,All;
+
+Amsud:=false;
+Amnord:=false;
+Amcentrale:=false;
+Caraibes:=false;
+Asie:=false;
+Europe:=false;
+Afrique:=false;
+All:=true;
+
+vardef Projection(expr TH,PH,Zoom)=
+  theta:=TH;
+  phi:=PH;
+  zoom(Zoom);
+  if All=true:
+    Lecture("Cameriquesud.dat");
+    Lecture("Ccaraibes.dat");
+    Lecture("Cameriquecentrale.dat");
+    Lecture("Cameriquenord.dat");
+    Lecture("CAsie.dat");
+    Lecture("Ceurope.dat");
+    Lecture("Cafrique.dat");
+  else:
+    if Amsud=true:
+      Lecture("Cameriquesud.dat");
+    fi;
+    if Amnord=true:
+      Lecture("Cameriquenord.dat");
+    fi;
+    if Amcentrale=true:
+      Lecture("Cameriquecentrale.dat");
+    fi;
+    if Caraibes=true:
+      Lecture("Ccaraibes.dat");
+    fi;
+    if Asie=true:
+      if projection="bonne":
+        Lecture("Casie.dat");
+      else:
+        Lecture("CAsie.dat");
+      fi;
+    fi;
+    if Europe=true:
+      Lecture("Ceurope.dat");
+    fi;
+    if Afrique=true:
+      Lecture("Cafrique.dat");
+    fi;
+  fi;
+  if lacs=true:
+    Lecturelacsp;
+    Lecturelacspsup;
+  fi;
+  Lectureilesp;
+  if capitales=true:
+    Lecturecapitalesp;
+  fi;
+  if fleuves=true:
+    Lecturerivieresp;
+    Lecturerivierespsup;
+  fi;
+  if maillage=true:
+    drawoptions(withcolor gris);
+    MeridienParallele
+    drawoptions();
+  fi;
+  if Echelle=true:
+    draw ((Xa,Ya)+u*(1,1))--((Xa,Ya)+u*(2,1));
+    label.top(btex 0 etex,(Xa,Ya)+u*(1,1));
+    label.top(TEX(""&decimal(Long)&"~km"),(Xa,Ya)+u*(2,1));
+  fi;
+enddef;
+
-vardef Simple(expr TH,PH,Zoom)=
-  projection:="simple";
@@ -955 +1219 @@
-  if lacs=true:
-    Lecturelacsp;
+    Lecturelacspsup;
-  fi;
-  Lectureilesp;
@@ -962 +1227 @@
-  if fleuves=true:
-    Lecturerivieresp;
-
+ 
-  if maillage=true:
-    MeridienParallele
@@ -982 +1247 @@
-  if lacs=true:
-    Lecturelacsp;
+    Lecturelacspsup;
-  fi;
-  Lectureilesp;
@@ -1009 +1275 @@
-  if lacs=true:
-    Lecturelacsp;
+    Lecturelacspsup;
-  fi;
-  Lectureilesp;
@@ -1037 +1304 @@
-  if lacs=true:
-    Lecturelacsp;
+    Lecturelacspsup;
-  fi;
-  Lectureilesp;
@@ -1064 +1332 @@
-  if lacs=true:
-    Lecturelacsp;
+    Lecturelacspsup;
-  fi;
-  Lectureilesp;
@@ -1071 +1340 @@
-  if fleuves=true:
-    Lecturerivieresp;
+    Lecturerivierespsup;
-  fi;
-  drawoptions(withcolor gris);
@@ -1088 +1358 @@
-  if lacs=true:
-    Lecturelacsp;
+    Lecturelacspsup;
-  fi;
-  Lectureilesp;
@@ -1095 +1366 @@
-  if fleuves=true:
-    Lecturerivieresp;
+    Lecturerivierespsup;
-  fi;
-  drawoptions(withcolor gris);
@@ -1138 +1410 @@
-  fi;
-  fill cercles(pte3,pte1,pte3,pte1,pte4) withcolor ciel;
-
-
-
-
-
+%
+%
+%
+%
+%
-  Lecture("Ceurope.dat");
-
+%
-  if volcans=true:
-    Lecturevolcans;
@@ -1150 +1422 @@
-  if lacs=true:
-    Lecturelacs;
+    Lecturelacssup;
-  fi;
-  Lectureiles;
@@ -1156 +1429 @@
-  fi;
-  if fleuves=true:
-
+
+
-  fi;
-  if maillage=true: