This work presents a comprehensive study of the propagation of femtosecond pulses and the formation and evolution of spatial solitons. The first half (Chapters 2-3) is devoted to the implementation of a novel ultrafast holographic system to capture the nonlinear propagation of laser pulses with femtosecond resolution. Femtosecond pulses are used to record holograms of the ultrafast changes in the material properties. Amplitude and phase changes of the laser beam inside the medium are reconstructed numerically. The strength of the nonlinear material response and the density of free electrons can be recovered from the phase information in the hologram. A single hologram can be captured with fine spatial resolution, or a time-sequence of holograms can be captured in a single shot with reduced spatial resolution. We have observed dramatic differences in the light propagation depending on the material properties. The second part of the thesis (Chapters 4-5) covers the formation and evolution of spatial solitons in a Kerr medium. We have measured the evolution of the beam profile as a function of pulse energy and propagation length. The optical beam breaks up into a pattern of connected lines (constellation) and self-focused spots (solitons). The solitons self-focus to a minimum diameter and release their excess energy through conical emission, which in turn overlaps with the background constellation and seeds the formation of new solitons. The solitons also show a collective self-organizing behavior caused by their mutual interactions. The evolution of 1-D arrays of solitons was captured using Femtosecond Time-resolved Optical Polarigraphy, a technique that measures the transient birefringence induced by the pulses in the medium. When the array was generated in an unstable configuration, the solitons re-arranged themselves into an array with a (larger) more stable period. A transition to a chaotic state is observed when the input power is increased above a threshold level. A time-averaged pulse propagation equation was used to numerically solve for evolution of the beam. There was good agreement between the experimental results and the computer simulation.